Actual source code: gcreate.c
1: #define PETSCMAT_DLL
3: #include src/mat/matimpl.h
4: #include petscsys.h
8: static PetscErrorCode MatPublish_Base(PetscObject obj)
9: {
11: return(0);
12: }
17: /*@C
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 PETSCMAT_DLLEXPORT 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->m = -1;
81: B->M = -1;
82: B->n = -1;
83: B->N = -1;
84: B->bs = 1;
85: B->preallocated = PETSC_FALSE;
86: B->bops->publish = MatPublish_Base;
87: *A = B;
88: return(0);
89: }
93: /*@
94: MatSetSizes - Sets the local and global sizes, and checks to determine compatibility
96: Collective on Mat
98: Input Parameters:
99: + A - the matrix
100: . m - number of local rows (or PETSC_DECIDE)
101: . n - number of local columns (or PETSC_DECIDE)
102: . M - number of global rows (or PETSC_DETERMINE)
103: - N - number of global columns (or PETSC_DETERMINE)
105: Notes:
106: m (n) and M (N) cannot be both PETSC_DECIDE
107: If one processor calls this with M (N) of PETSC_DECIDE then all processors must, otherwise the program will hang.
109: If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
110: user must ensure that they are chosen to be compatible with the
111: vectors. To do this, one first considers the matrix-vector product
112: 'y = A x'. The 'm' that is used in the above routine must match the
113: local size used in the vector creation routine VecCreateMPI() for 'y'.
114: Likewise, the 'n' used must match that used as the local size in
115: VecCreateMPI() for 'x'.
117: Level: beginner
119: .seealso: MatGetSize(), PetscSplitOwnership()
120: @*/
121: PetscErrorCode PETSCMAT_DLLEXPORT MatSetSizes(Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N)
122: {
125: if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M);
126: if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N);
127: if ((A->m >= 0 || A->M >= 0) && (A->m != m || A->M != 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->m,A->M);
128: if ((A->n >= 0 || A->N >= 0) && (A->n != n || A->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->n,A->N);
129: A->m = m;
130: A->n = n;
131: A->M = M;
132: A->N = N;
133: return(0);
134: }
138: /*@C
139: MatSetFromOptions - Creates a matrix where the type is determined
140: from the options database. Generates a parallel MPI matrix if the
141: communicator has more than one processor. The default matrix type is
142: AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
143: you do not select a type in the options database.
145: Collective on Mat
147: Input Parameter:
148: . A - the matrix
150: Options Database Keys:
151: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
152: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
153: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
154: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
155: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
156: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
157: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
158: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
159: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
161: Even More Options Database Keys:
162: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
163: for additional format-specific options.
165: Level: beginner
167: .keywords: matrix, create
169: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
170: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
171: MatCreateSeqDense(), MatCreateMPIDense(),
172: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
173: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
174: MatConvert()
175: @*/
176: PetscErrorCode PETSCMAT_DLLEXPORT MatSetFromOptions(Mat B)
177: {
179: char mtype[256];
180: PetscTruth flg;
183: PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
184: if (flg) {
185: MatSetType(B,mtype);
186: }
187: if (!B->type_name) {
188: MatSetType(B,MATAIJ);
189: }
190: return(0);
191: }
195: /*@C
196: MatSetUpPreallocation
198: Collective on Mat
200: Input Parameter:
201: . A - the matrix
203: Level: beginner
205: .keywords: matrix, create
207: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
208: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
209: MatCreateSeqDense(), MatCreateMPIDense(),
210: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
211: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
212: MatConvert()
213: @*/
214: PetscErrorCode PETSCMAT_DLLEXPORT MatSetUpPreallocation(Mat B)
215: {
219: if (B->ops->setuppreallocation) {
220: PetscLogInfo((B,"MatSetUpPreallocation: Warning not preallocating matrix storage\n"));
221: (*B->ops->setuppreallocation)(B);
222: B->ops->setuppreallocation = 0;
223: B->preallocated = PETSC_TRUE;
224: }
225: return(0);
226: }
228: /*
229: Copies from Cs header to A
230: */
233: PetscErrorCode MatHeaderCopy(Mat A,Mat C)
234: {
236: PetscInt refct;
237: PetscOps *Abops;
238: MatOps Aops;
239: char *mtype,*mname;
240: void *spptr;
243: /* free all the interior data structures from mat */
244: (*A->ops->destroy)(A);
246: PetscMapDestroy(A->rmap);
247: PetscMapDestroy(A->cmap);
249: /* save the parts of A we need */
250: Abops = A->bops;
251: Aops = A->ops;
252: refct = A->refct;
253: mtype = A->type_name;
254: mname = A->name;
255: spptr = A->spptr;
257: if (C->spptr) {
258: PetscFree(C->spptr);
259: C->spptr = PETSC_NULL;
260: }
262: /* copy C over to A */
263: PetscMemcpy(A,C,sizeof(struct _p_Mat));
265: /* return the parts of A we saved */
266: A->bops = Abops;
267: A->ops = Aops;
268: A->qlist = 0;
269: A->refct = refct;
270: A->type_name = mtype;
271: A->name = mname;
272: A->spptr = spptr;
274: PetscHeaderDestroy(C);
275: return(0);
276: }
277: /*
278: Replace A's header with that of C
279: This is essentially code moved from MatDestroy
280: */
283: PetscErrorCode MatHeaderReplace(Mat A,Mat C)
284: {
288: /* free all the interior data structures from mat */
289: (*A->ops->destroy)(A);
290: PetscHeaderDestroy_Private((PetscObject)A);
291: PetscMapDestroy(A->rmap);
292: PetscMapDestroy(A->cmap);
294: /* copy C over to A */
295: if (C) {
296: PetscMemcpy(A,C,sizeof(struct _p_Mat));
297: PetscLogObjectDestroy((PetscObject)C);
298: PetscFree(C);
299: }
300: return(0);
301: }