Actual source code: gcreate.c
1: #define PETSCMAT_DLL
3: #include include/private/matimpl.h
4: #include petscsys.h
6: #if 0
9: static PetscErrorCode MatPublish_Base(PetscObject obj)
10: {
12: return(0);
13: }
14: #endif
18: /*@
19: MatCreate - Creates a matrix where the type is determined
20: from either a call to MatSetType() or from the options database
21: with a call to MatSetFromOptions(). The default matrix type is
22: AIJ, using the routines MatCreateSeqAIJ() or MatCreateMPIAIJ()
23: if you do not set a type in the options database. If you never
24: call MatSetType() or MatSetFromOptions() it will generate an
25: error when you try to use the matrix.
27: Collective on MPI_Comm
29: Input Parameter:
30: . comm - MPI communicator
31:
32: Output Parameter:
33: . A - the matrix
35: Options Database Keys:
36: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
37: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
38: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
39: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
40: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
41: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
42: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
43: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
44: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
46: Even More Options Database Keys:
47: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
48: for additional format-specific options.
50: Notes:
52: Level: beginner
54: User manual sections:
55: + Section 3.1 Creating and Assembling Matrices
56: - Chapter 3 Matrices
58: .keywords: matrix, create
60: .seealso: MatCreateSeqAIJ(), MatCreateMPIAIJ(),
61: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
62: MatCreateSeqDense(), MatCreateMPIDense(),
63: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
64: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
65: MatConvert()
66: @*/
67: PetscErrorCode MatCreate(MPI_Comm comm,Mat *A)
68: {
69: Mat B;
75: *A = PETSC_NULL;
76: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
77: MatInitializePackage(PETSC_NULL);
78: #endif
80: PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);
81: PetscMapInitialize(comm,&B->rmap);
82: PetscMapInitialize(comm,&B->cmap);
83: B->preallocated = PETSC_FALSE;
84: *A = B;
85: return(0);
86: }
90: /*@
91: MatSetSizes - Sets the local and global sizes, and checks to determine compatibility
93: Collective on Mat
95: Input Parameters:
96: + A - the matrix
97: . m - number of local rows (or PETSC_DECIDE)
98: . n - number of local columns (or PETSC_DECIDE)
99: . M - number of global rows (or PETSC_DETERMINE)
100: - N - number of global columns (or PETSC_DETERMINE)
102: Notes:
103: m (n) and M (N) cannot be both PETSC_DECIDE
104: If one processor calls this with M (N) of PETSC_DECIDE then all processors must, otherwise the program will hang.
106: If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
107: user must ensure that they are chosen to be compatible with the
108: vectors. To do this, one first considers the matrix-vector product
109: 'y = A x'. The 'm' that is used in the above routine must match the
110: local size used in the vector creation routine VecCreateMPI() for 'y'.
111: Likewise, the 'n' used must match that used as the local size in
112: VecCreateMPI() for 'x'.
114: Level: beginner
116: .seealso: MatGetSize(), PetscSplitOwnership()
117: @*/
118: PetscErrorCode MatSetSizes(Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N)
119: {
124: if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M);
125: if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N);
126: if (A->ops->setsizes) {
127: /* Since this will not be set until the type has been set, this will NOT be called on the initial
128: call of MatSetSizes() (which must be called BEFORE MatSetType() */
129: (*A->ops->setsizes)(A,m,n,M,N);
130: } else {
131: if ((A->rmap.n >= 0 || A->rmap.N >= 0) && (A->rmap.n != m || A->rmap.N != M)) SETERRQ4(PETSC_ERR_SUP,"Cannot change/reset row sizes to %D local %D global after previously setting them to %D local %D global",m,M,A->rmap.n,A->rmap.N);
132: if ((A->cmap.n >= 0 || A->cmap.N >= 0) && (A->cmap.n != n || A->cmap.N != N)) SETERRQ4(PETSC_ERR_SUP,"Cannot change/reset column sizes to %D local %D global after previously setting them to %D local %D global",n,N,A->cmap.n,A->cmap.N);
133: }
134: A->rmap.n = m;
135: A->cmap.n = n;
136: A->rmap.N = M;
137: A->cmap.N = N;
138: return(0);
139: }
143: /*@
144: MatSetFromOptions - Creates a matrix where the type is determined
145: from the options database. Generates a parallel MPI matrix if the
146: communicator has more than one processor. The default matrix type is
147: AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
148: you do not select a type in the options database.
150: Collective on Mat
152: Input Parameter:
153: . A - the matrix
155: Options Database Keys:
156: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
157: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
158: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
159: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
160: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
161: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
162: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
163: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
164: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
166: Even More Options Database Keys:
167: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
168: for additional format-specific options.
170: Level: beginner
172: .keywords: matrix, create
174: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
175: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
176: MatCreateSeqDense(), MatCreateMPIDense(),
177: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
178: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
179: MatConvert()
180: @*/
181: PetscErrorCode MatSetFromOptions(Mat B)
182: {
184: const char *deft = MATAIJ;
185: char type[256];
186: PetscTruth flg;
191: PetscOptionsBegin(((PetscObject)B)->comm,((PetscObject)B)->prefix,"Matrix options","Mat");
193: if (((PetscObject)B)->type_name) { deft = ((PetscObject)B)->type_name; }
194: PetscOptionsList("-mat_type","Matrix type","MatSetType",MatList,deft,type,256,&flg);
195: if (flg) {
196: MatSetType(B,type);
197: } else if (!((PetscObject)B)->type_name) {
198: MatSetType(B,deft);
199: }
201: if (B->ops->setfromoptions) {
202: (*B->ops->setfromoptions)(B);
203: }
205: PetscOptionsEnd();
207: return(0);
208: }
212: /*@
213: MatSetUpPreallocation
215: Collective on Mat
217: Input Parameter:
218: . A - the matrix
220: Level: beginner
222: .keywords: matrix, create
224: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
225: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
226: MatCreateSeqDense(), MatCreateMPIDense(),
227: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
228: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
229: MatConvert()
230: @*/
231: PetscErrorCode MatSetUpPreallocation(Mat B)
232: {
236: if (!B->preallocated && B->ops->setuppreallocation) {
237: PetscInfo(B,"Warning not preallocating matrix storage\n");
238: (*B->ops->setuppreallocation)(B);
239: }
240: B->preallocated = PETSC_TRUE;
241: return(0);
242: }
244: /*
245: Copies from Cs header to A
246: */
249: PetscErrorCode MatHeaderCopy(Mat A,Mat C)
250: {
252: PetscInt refct;
253: PetscOps *Abops;
254: MatOps Aops;
255: char *mtype,*mname;
256: void *spptr;
259: /* save the parts of A we need */
260: Abops = ((PetscObject)A)->bops;
261: Aops = A->ops;
262: refct = ((PetscObject)A)->refct;
263: mtype = ((PetscObject)A)->type_name; ((PetscObject)A)->type_name = 0;
264: mname = ((PetscObject)A)->name; ((PetscObject)A)->name = 0;
265: spptr = A->spptr;
267: /* free all the interior data structures from mat */
268: (*A->ops->destroy)(A);
270: PetscFree(C->spptr);
272: PetscFree(A->rmap.range);
273: PetscFree(A->cmap.range);
275: /* copy C over to A */
276: PetscMemcpy(A,C,sizeof(struct _p_Mat));
278: /* return the parts of A we saved */
279: ((PetscObject)A)->bops = Abops;
280: A->ops = Aops;
281: ((PetscObject)A)->qlist = 0;
282: ((PetscObject)A)->refct = refct;
283: ((PetscObject)A)->type_name = mtype;
284: ((PetscObject)A)->name = mname;
285: A->spptr = spptr;
287: PetscHeaderDestroy(C);
288: return(0);
289: }
290: /*
291: Replace A's header with that of C
292: This is essentially code moved from MatDestroy
293: */
296: PetscErrorCode MatHeaderReplace(Mat A,Mat C)
297: {
301: /* free all the interior data structures from mat */
302: (*A->ops->destroy)(A);
303: PetscHeaderDestroy_Private((PetscObject)A);
304: PetscFree(A->ops);
305: PetscFree(A->rmap.range);
306: PetscFree(A->cmap.range);
307: PetscFree(A->spptr);
308:
309: /* copy C over to A */
310: if (C) {
311: PetscMemcpy(A,C,sizeof(struct _p_Mat));
312: PetscLogObjectDestroy((PetscObject)C);
313: PetscFree(C);
314: }
315: return(0);
316: }