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: PetscMapInitialize(comm,&B->rmap);
81: PetscMapInitialize(comm,&B->cmap);
82: B->preallocated = PETSC_FALSE;
83: B->bops->publish = MatPublish_Base;
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: char mtype[256];
185: PetscTruth flg;
188: PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
189: if (flg) {
190: MatSetType(B,mtype);
191: }
192: if (!B->type_name) {
193: MatSetType(B,MATAIJ);
194: }
195: return(0);
196: }
200: /*@C
201: MatSetUpPreallocation
203: Collective on Mat
205: Input Parameter:
206: . A - the matrix
208: Level: beginner
210: .keywords: matrix, create
212: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
213: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
214: MatCreateSeqDense(), MatCreateMPIDense(),
215: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
216: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
217: MatConvert()
218: @*/
219: PetscErrorCode MatSetUpPreallocation(Mat B)
220: {
224: if (!B->preallocated && B->ops->setuppreallocation) {
225: PetscInfo(B,"Warning not preallocating matrix storage\n");
226: (*B->ops->setuppreallocation)(B);
227: }
228: B->preallocated = PETSC_TRUE;
229: return(0);
230: }
232: /*
233: Copies from Cs header to A
234: */
237: PetscErrorCode MatHeaderCopy(Mat A,Mat C)
238: {
240: PetscInt refct;
241: PetscOps *Abops;
242: MatOps Aops;
243: char *mtype,*mname;
244: void *spptr;
247: /* save the parts of A we need */
248: Abops = A->bops;
249: Aops = A->ops;
250: refct = A->refct;
251: mtype = A->type_name; A->type_name = 0;
252: mname = A->name; A->name = 0;
253: spptr = A->spptr;
255: /* free all the interior data structures from mat */
256: (*A->ops->destroy)(A);
258: PetscFree(C->spptr);
260: PetscFree(A->rmap.range);
261: PetscFree(A->cmap.range);
263: /* copy C over to A */
264: PetscMemcpy(A,C,sizeof(struct _p_Mat));
266: /* return the parts of A we saved */
267: A->bops = Abops;
268: A->ops = Aops;
269: A->qlist = 0;
270: A->refct = refct;
271: A->type_name = mtype;
272: A->name = mname;
273: A->spptr = spptr;
275: PetscHeaderDestroy(C);
276: return(0);
277: }
278: /*
279: Replace A's header with that of C
280: This is essentially code moved from MatDestroy
281: */
284: PetscErrorCode MatHeaderReplace(Mat A,Mat C)
285: {
289: /* free all the interior data structures from mat */
290: (*A->ops->destroy)(A);
291: PetscHeaderDestroy_Private((PetscObject)A);
292: PetscFree(A->rmap.range);
293: PetscFree(A->cmap.range);
294: PetscFree(A->spptr);
295:
296: /* copy C over to A */
297: if (C) {
298: PetscMemcpy(A,C,sizeof(struct _p_Mat));
299: PetscLogObjectDestroy((PetscObject)C);
300: PetscFree(C);
301: }
302: return(0);
303: }