Actual source code: gcreate.c
1: #define PETSCMAT_DLL
3: #include include/private/matimpl.h
4: #include petscsys.h
8: static PetscErrorCode MatPublish_Base(PetscObject obj)
9: {
11: return(0);
12: }
17: /*@
18: MatCreate - Creates a matrix where the type is determined
19: from either a call to MatSetType() or from the options database
20: with a call to MatSetFromOptions(). The default matrix type is
21: AIJ, using the routines MatCreateSeqAIJ() or MatCreateMPIAIJ()
22: if you do not set a type in the options database. If you never
23: call MatSetType() or MatSetFromOptions() it will generate an
24: error when you try to use the matrix.
26: Collective on MPI_Comm
28: Input Parameter:
29: . comm - MPI communicator
30:
31: Output Parameter:
32: . A - the matrix
34: Options Database Keys:
35: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
36: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
37: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
38: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
39: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
40: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
41: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
42: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
43: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
45: Even More Options Database Keys:
46: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
47: for additional format-specific options.
49: Notes:
51: Level: beginner
53: User manual sections:
54: + Section 3.1 Creating and Assembling Matrices
55: - Chapter 3 Matrices
57: .keywords: matrix, create
59: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
60: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
61: MatCreateSeqDense(), MatCreateMPIDense(),
62: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
63: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
64: MatConvert()
65: @*/
66: PetscErrorCode MatCreate(MPI_Comm comm,Mat *A)
67: {
68: Mat B;
74: *A = PETSC_NULL;
75: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
76: MatInitializePackage(PETSC_NULL);
77: #endif
79: PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);
80: B->rmap.n = -1;
81: B->rmap.N = -1;
82: B->cmap.n = -1;
83: B->cmap.N = -1;
84: B->rmap.bs = 1;
85: B->cmap.bs = 1;
86: B->preallocated = PETSC_FALSE;
87: B->bops->publish = MatPublish_Base;
88: *A = B;
89: return(0);
90: }
94: /*@
95: MatSetSizes - Sets the local and global sizes, and checks to determine compatibility
97: Collective on Mat
99: Input Parameters:
100: + A - the matrix
101: . m - number of local rows (or PETSC_DECIDE)
102: . n - number of local columns (or PETSC_DECIDE)
103: . M - number of global rows (or PETSC_DETERMINE)
104: - N - number of global columns (or PETSC_DETERMINE)
106: Notes:
107: m (n) and M (N) cannot be both PETSC_DECIDE
108: If one processor calls this with M (N) of PETSC_DECIDE then all processors must, otherwise the program will hang.
110: If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
111: user must ensure that they are chosen to be compatible with the
112: vectors. To do this, one first considers the matrix-vector product
113: 'y = A x'. The 'm' that is used in the above routine must match the
114: local size used in the vector creation routine VecCreateMPI() for 'y'.
115: Likewise, the 'n' used must match that used as the local size in
116: VecCreateMPI() for 'x'.
118: Level: beginner
120: .seealso: MatGetSize(), PetscSplitOwnership()
121: @*/
122: PetscErrorCode MatSetSizes(Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N)
123: {
128: if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M);
129: if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N);
130: if (A->ops->setsizes) {
131: /* Since this will not be set until the type has been set, this will NOT be called on the initial
132: call of MatSetSizes() (which must be called BEFORE MatSetType() */
133: (*A->ops->setsizes)(A,m,n,M,N);
134: } else {
135: 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);
136: 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);
137: }
138: A->rmap.n = m;
139: A->cmap.n = n;
140: A->rmap.N = M;
141: A->cmap.N = N;
142: return(0);
143: }
147: /*@
148: MatSetFromOptions - Creates a matrix where the type is determined
149: from the options database. Generates a parallel MPI matrix if the
150: communicator has more than one processor. The default matrix type is
151: AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
152: you do not select a type in the options database.
154: Collective on Mat
156: Input Parameter:
157: . A - the matrix
159: Options Database Keys:
160: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
161: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
162: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
163: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
164: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
165: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
166: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
167: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
168: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
170: Even More Options Database Keys:
171: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
172: for additional format-specific options.
174: Level: beginner
176: .keywords: matrix, create
178: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
179: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
180: MatCreateSeqDense(), MatCreateMPIDense(),
181: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
182: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
183: MatConvert()
184: @*/
185: PetscErrorCode MatSetFromOptions(Mat B)
186: {
188: char mtype[256];
189: PetscTruth flg;
192: PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
193: if (flg) {
194: MatSetType(B,mtype);
195: }
196: if (!B->type_name) {
197: MatSetType(B,MATAIJ);
198: }
199: return(0);
200: }
204: /*@C
205: MatSetUpPreallocation
207: Collective on Mat
209: Input Parameter:
210: . A - the matrix
212: Level: beginner
214: .keywords: matrix, create
216: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
217: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
218: MatCreateSeqDense(), MatCreateMPIDense(),
219: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
220: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
221: MatConvert()
222: @*/
223: PetscErrorCode MatSetUpPreallocation(Mat B)
224: {
228: if (!B->preallocated && B->ops->setuppreallocation) {
229: PetscInfo(B,"Warning not preallocating matrix storage\n");
230: (*B->ops->setuppreallocation)(B);
231: }
232: B->preallocated = PETSC_TRUE;
233: return(0);
234: }
236: /*
237: Copies from Cs header to A
238: */
241: PetscErrorCode MatHeaderCopy(Mat A,Mat C)
242: {
244: PetscInt refct;
245: PetscOps *Abops;
246: MatOps Aops;
247: char *mtype,*mname;
248: void *spptr;
251: /* save the parts of A we need */
252: Abops = A->bops;
253: Aops = A->ops;
254: refct = A->refct;
255: mtype = A->type_name; A->type_name = 0;
256: mname = A->name; A->name = 0;
257: spptr = A->spptr;
259: /* free all the interior data structures from mat */
260: (*A->ops->destroy)(A);
262: PetscFree(C->spptr);
264: PetscFree(A->rmap.range);
265: PetscFree(A->cmap.range);
267: /* copy C over to A */
268: PetscMemcpy(A,C,sizeof(struct _p_Mat));
270: /* return the parts of A we saved */
271: A->bops = Abops;
272: A->ops = Aops;
273: A->qlist = 0;
274: A->refct = refct;
275: A->type_name = mtype;
276: A->name = mname;
277: A->spptr = spptr;
279: PetscHeaderDestroy(C);
280: return(0);
281: }
282: /*
283: Replace A's header with that of C
284: This is essentially code moved from MatDestroy
285: */
288: PetscErrorCode MatHeaderReplace(Mat A,Mat C)
289: {
293: /* free all the interior data structures from mat */
294: (*A->ops->destroy)(A);
295: PetscHeaderDestroy_Private((PetscObject)A);
296: PetscFree(A->rmap.range);
297: PetscFree(A->cmap.range);
298: PetscFree(A->spptr);
299:
300: /* copy C over to A */
301: if (C) {
302: PetscMemcpy(A,C,sizeof(struct _p_Mat));
303: PetscLogObjectDestroy((PetscObject)C);
304: PetscFree(C);
305: }
306: return(0);
307: }