Actual source code: gcreate.c
1: /*$Id: gcreate.c,v 1.131 2001/07/20 21:22:13 bsmith Exp $*/
3: #include src/mat/matimpl.h
4: #include petscsys.h
6: #undef __FUNCT__
8: static int MatPublish_Base(PetscObject obj)
9: {
10: #if defined(PETSC_HAVE_AMS)
11: Mat mat = (Mat)obj;
13: #endif
16: #if defined(PETSC_HAVE_AMS)
17: /* if it is already published then return */
18: if (mat->amem >=0) return(0);
20: PetscObjectPublishBaseBegin(obj);
21: AMS_Memory_add_field((AMS_Memory)mat->amem,"globalsizes",&mat->M,2,AMS_INT,AMS_READ,
22: AMS_COMMON,AMS_REDUCT_UNDEF);
23: AMS_Memory_add_field((AMS_Memory)mat->amem,"localsizes",&mat->m,2,AMS_INT,AMS_READ,
24: AMS_DISTRIBUTED,AMS_REDUCT_UNDEF);
25: PetscObjectPublishBaseEnd(obj);
26: #endif
28: return(0);
29: }
32: #undef __FUNCT__
34: /*@C
35: MatCreate - Creates a matrix where the type is determined
36: from either a call to MatSetType() or from the options database
37: with a call to MatSetFromOptions(). The default matrix type is
38: AIJ, using the routines MatCreateSeqAIJ() or MatCreateMPIAIJ()
39: if you do not set a type in the options database. If you never
40: call MatSetType() or MatSetFromOptions() it will generate an
41: error when you try to use the matrix.
43: Collective on MPI_Comm
45: Input Parameters:
46: + m - number of local rows (or PETSC_DECIDE)
47: . n - number of local columns (or PETSC_DECIDE)
48: . M - number of global rows (or PETSC_DETERMINE)
49: . N - number of global columns (or PETSC_DETERMINE)
50: - comm - MPI communicator
51:
52: Output Parameter:
53: . A - the matrix
55: Options Database Keys:
56: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
57: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
58: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
59: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
60: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
61: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
62: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
63: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
64: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
66: Even More Options Database Keys:
67: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
68: for additional format-specific options.
70: Notes:
71: If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
72: user must ensure that they are chosen to be compatible with the
73: vectors. To do this, one first considers the matrix-vector product
74: 'y = A x'. The 'm' that is used in the above routine must match the
75: local size used in the vector creation routine VecCreateMPI() for 'y'.
76: Likewise, the 'n' used must match that used as the local size in
77: VecCreateMPI() for 'x'.
79: Level: beginner
81: User manual sections:
82: + Section 3.1 Creating and Assembling Matrices
83: - Chapter 3 Matrices
85: .keywords: matrix, create
87: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
88: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
89: MatCreateSeqDense(), MatCreateMPIDense(),
90: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
91: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
92: MatConvert()
93: @*/
94: int MatCreate(MPI_Comm comm,int m,int n,int M,int N,Mat *A)
95: {
96: Mat B;
97: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
99: #endif
103: *A = PETSC_NULL;
104: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
105: MatInitializePackage(PETSC_NULL);
106: #endif
108: PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);
109: PetscLogObjectCreate(B);
111: B->m = m;
112: B->n = n;
113: B->M = M;
114: B->N = N;
116: B->preallocated = PETSC_FALSE;
117: B->bops->publish = MatPublish_Base;
118: *A = B;
119: return(0);
120: }
122: #undef __FUNCT__
124: /*@C
125: MatSetFromOptions - Creates a matrix where the type is determined
126: from the options database. Generates a parallel MPI matrix if the
127: communicator has more than one processor. The default matrix type is
128: AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
129: you do not select a type in the options database.
131: Collective on Mat
133: Input Parameter:
134: . A - the matrix
136: Options Database Keys:
137: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
138: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
139: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
140: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
141: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
142: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
143: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
144: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
145: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
147: Even More Options Database Keys:
148: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
149: for additional format-specific options.
151: Level: beginner
153: .keywords: matrix, create
155: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
156: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
157: MatCreateSeqDense(), MatCreateMPIDense(),
158: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
159: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
160: MatConvert()
161: @*/
162: int MatSetFromOptions(Mat B)
163: {
164: int ierr,size;
165: char mtype[256];
166: PetscTruth flg;
169: PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
170: if (flg) {
171: MatSetType(B,mtype);
172: }
173: if (!B->type_name) {
174: MPI_Comm_size(B->comm,&size);
175: if (size == 1) {
176: MatSetType(B,MATSEQAIJ);
177: } else {
178: MatSetType(B,MATMPIAIJ);
179: }
180: }
181: #if defined(__cplusplus) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE) && defined(PETSC_HAVE_CXX_NAMESPACE)
182: MatESISetFromOptions(B);
183: #endif
184: return(0);
185: }
187: #undef __FUNCT__
189: /*@C
190: MatSetUpPreallocation
192: Collective on Mat
194: Input Parameter:
195: . A - the matrix
197: Level: beginner
199: .keywords: matrix, create
201: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
202: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
203: MatCreateSeqDense(), MatCreateMPIDense(),
204: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
205: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
206: MatConvert()
207: @*/
208: int MatSetUpPreallocation(Mat B)
209: {
210: int ierr;
213: if (B->ops->setuppreallocation) {
214: PetscLogInfo(B,"MatSetUpPreallocation: Warning not preallocating matrix storage");
215: (*B->ops->setuppreallocation)(B);
216: B->ops->setuppreallocation = 0;
217: B->preallocated = PETSC_TRUE;
218: }
219: return(0);
220: }
222: /*
223: Copies from Cs header to A
224: */
225: #undef __FUNCT__
227: int MatHeaderCopy(Mat A,Mat C)
228: {
229: int ierr,refct;
230: PetscOps *Abops;
231: MatOps Aops;
232: char *mtype,*mname;
235: /* free all the interior data structures from mat */
236: (*A->ops->destroy)(A);
238: PetscMapDestroy(A->rmap);
239: PetscMapDestroy(A->cmap);
241: /* save the parts of A we need */
242: Abops = A->bops;
243: Aops = A->ops;
244: refct = A->refct;
245: mtype = A->type_name;
246: mname = A->name;
248: /* copy C over to A */
249: ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));
251: /* return the parts of A we saved */
252: A->bops = Abops;
253: A->ops = Aops;
254: A->qlist = 0;
255: A->refct = refct;
256: A->type_name = mtype;
257: A->name = mname;
259: PetscLogObjectDestroy(C);
260: PetscHeaderDestroy(C);
261: return(0);
262: }