xref: /petsc/src/mat/utils/gcreate.c (revision 7c4f633dc6bb6149cca88d301ead35a99e103cbb)
1 #define PETSCMAT_DLL
2 
3 #include "private/matimpl.h"       /*I "petscmat.h"  I*/
4 #include "petscsys.h"
5 
6 #if 0
7 #undef __FUNCT__
8 #define __FUNCT__ "MatPublish_Base"
9 static PetscErrorCode MatPublish_Base(PetscObject obj)
10 {
11   PetscFunctionBegin;
12   PetscFunctionReturn(0);
13 }
14 #endif
15 
16 #undef __FUNCT__
17 #define __FUNCT__ "MatCreate"
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.
26 
27    Collective on MPI_Comm
28 
29    Input Parameter:
30 .  comm - MPI communicator
31 
32    Output Parameter:
33 .  A - the matrix
34 
35    Options Database Keys:
36 +    -mat_type seqaij   - AIJ type, uses MatCreateSeqAIJ()
37 .    -mat_type mpiaij   - AIJ type, uses MatCreateMPIAIJ()
38 .    -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
39 .    -mat_type seqdense - dense type, uses MatCreateSeqDense()
40 .    -mat_type mpidense - dense type, uses MatCreateMPIDense()
41 .    -mat_type seqbaij  - block AIJ type, uses MatCreateSeqBAIJ()
42 -    -mat_type mpibaij  - block AIJ type, uses MatCreateMPIBAIJ()
43 
44    Even More Options Database Keys:
45    See the manpages for particular formats (e.g., MatCreateSeqAIJ())
46    for additional format-specific options.
47 
48    Notes:
49 
50    Level: beginner
51 
52    User manual sections:
53 +   sec_matcreate
54 -   chapter_matrices
55 
56 .keywords: matrix, create
57 
58 .seealso: MatCreateSeqAIJ(), MatCreateMPIAIJ(),
59           MatCreateSeqDense(), MatCreateMPIDense(),
60           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
61           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
62           MatConvert()
63 @*/
64 PetscErrorCode PETSCMAT_DLLEXPORT MatCreate(MPI_Comm comm,Mat *A)
65 {
66   Mat            B;
67   PetscErrorCode ierr;
68 
69   PetscFunctionBegin;
70   PetscValidPointer(A,2);
71 
72   *A = PETSC_NULL;
73 #ifndef PETSC_USE_DYNAMIC_LIBRARIES
74   ierr = MatInitializePackage(PETSC_NULL);CHKERRQ(ierr);
75 #endif
76 
77   ierr = PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);CHKERRQ(ierr);
78   ierr = PetscNew(PetscMap,&B->rmap);CHKERRQ(ierr);
79   ierr = PetscNew(PetscMap,&B->cmap);CHKERRQ(ierr);
80   ierr = PetscMapInitialize(comm,B->rmap);CHKERRQ(ierr);
81   ierr = PetscMapInitialize(comm,B->cmap);CHKERRQ(ierr);
82   B->preallocated  = PETSC_FALSE;
83   *A               = B;
84   PetscFunctionReturn(0);
85 }
86 
87 #undef __FUNCT__
88 #define __FUNCT__ "MatSetSizes"
89 /*@
90   MatSetSizes - Sets the local and global sizes, and checks to determine compatibility
91 
92   Collective on Mat
93 
94   Input Parameters:
95 +  A - the matrix
96 .  m - number of local rows (or PETSC_DECIDE)
97 .  n - number of local columns (or PETSC_DECIDE)
98 .  M - number of global rows (or PETSC_DETERMINE)
99 -  N - number of global columns (or PETSC_DETERMINE)
100 
101    Notes:
102    m (n) and M (N) cannot be both PETSC_DECIDE
103    If one processor calls this with M (N) of PETSC_DECIDE then all processors must, otherwise the program will hang.
104 
105    If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
106    user must ensure that they are chosen to be compatible with the
107    vectors. To do this, one first considers the matrix-vector product
108    'y = A x'. The 'm' that is used in the above routine must match the
109    local size used in the vector creation routine VecCreateMPI() for 'y'.
110    Likewise, the 'n' used must match that used as the local size in
111    VecCreateMPI() for 'x'.
112 
113   Level: beginner
114 
115 .seealso: MatGetSize(), PetscSplitOwnership()
116 @*/
117 PetscErrorCode PETSCMAT_DLLEXPORT MatSetSizes(Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N)
118 {
119   PetscErrorCode ierr;
120 
121   PetscFunctionBegin;
122   PetscValidHeaderSpecific(A,MAT_COOKIE,1);
123   if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M);
124   if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N);
125   if (A->ops->setsizes) {
126     /* Since this will not be set until the type has been set, this will NOT be called on the initial
127        call of MatSetSizes() (which must be called BEFORE MatSetType() */
128     ierr = (*A->ops->setsizes)(A,m,n,M,N);CHKERRQ(ierr);
129   } else {
130     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);
131     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);
132   }
133   A->rmap->n = m;
134   A->cmap->n = n;
135   A->rmap->N = M;
136   A->cmap->N = N;
137   if (A->ops->create) {
138     ierr = (*A->ops->create)(A);CHKERRQ(ierr);
139     A->ops->create = 0;
140   }
141 
142   PetscFunctionReturn(0);
143 }
144 
145 #undef __FUNCT__
146 #define __FUNCT__ "MatSetFromOptions"
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.
153 
154    Collective on Mat
155 
156    Input Parameter:
157 .  A - the matrix
158 
159    Options Database Keys:
160 +    -mat_type seqaij   - AIJ type, uses MatCreateSeqAIJ()
161 .    -mat_type mpiaij   - AIJ type, uses MatCreateMPIAIJ()
162 .    -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
163 .    -mat_type seqdense - dense type, uses MatCreateSeqDense()
164 .    -mat_type mpidense - dense type, uses MatCreateMPIDense()
165 .    -mat_type seqbaij  - block AIJ type, uses MatCreateSeqBAIJ()
166 -    -mat_type mpibaij  - block AIJ type, uses MatCreateMPIBAIJ()
167 
168    Even More Options Database Keys:
169    See the manpages for particular formats (e.g., MatCreateSeqAIJ())
170    for additional format-specific options.
171 
172    Level: beginner
173 
174 .keywords: matrix, create
175 
176 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
177           MatCreateSeqDense(), MatCreateMPIDense(),
178           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
179           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
180           MatConvert()
181 @*/
182 PetscErrorCode PETSCMAT_DLLEXPORT MatSetFromOptions(Mat B)
183 {
184   PetscErrorCode ierr;
185   const char     *deft = MATAIJ;
186   char           type[256];
187   PetscTruth     flg;
188 
189   PetscFunctionBegin;
190   PetscValidHeaderSpecific(B,MAT_COOKIE,1);
191 
192   ierr = PetscOptionsBegin(((PetscObject)B)->comm,((PetscObject)B)->prefix,"Matrix options","Mat");CHKERRQ(ierr);
193     ierr = PetscOptionsList("-mat_type","Matrix type","MatSetType",MatList,deft,type,256,&flg);CHKERRQ(ierr);
194     if (flg) {
195       ierr = MatSetType(B,type);CHKERRQ(ierr);
196     } else if (!((PetscObject)B)->type_name) {
197       ierr = MatSetType(B,deft);CHKERRQ(ierr);
198     }
199 
200     if (B->ops->setfromoptions) {
201       ierr = (*B->ops->setfromoptions)(B);CHKERRQ(ierr);
202     }
203 
204   ierr = PetscOptionsEnd();CHKERRQ(ierr);
205 
206   PetscFunctionReturn(0);
207 }
208 
209 #undef __FUNCT__
210 #define __FUNCT__ "MatSetUpPreallocation"
211 /*@
212    MatSetUpPreallocation
213 
214    Collective on Mat
215 
216    Input Parameter:
217 .  A - the matrix
218 
219    Level: beginner
220 
221 .keywords: matrix, create
222 
223 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
224           MatCreateSeqDense(), MatCreateMPIDense(),
225           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
226           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
227           MatConvert()
228 @*/
229 PetscErrorCode PETSCMAT_DLLEXPORT MatSetUpPreallocation(Mat B)
230 {
231   PetscErrorCode ierr;
232 
233   PetscFunctionBegin;
234   if (!B->preallocated && B->ops->setuppreallocation) {
235     ierr = PetscInfo(B,"Warning not preallocating matrix storage\n");CHKERRQ(ierr);
236     ierr = (*B->ops->setuppreallocation)(B);CHKERRQ(ierr);
237   }
238   B->preallocated = PETSC_TRUE;
239   PetscFunctionReturn(0);
240 }
241 
242 /*
243         Copies from Cs header to A
244 
245         This is somewhat different from MatHeaderReplace() it would be nice to merge the code
246 */
247 #undef __FUNCT__
248 #define __FUNCT__ "MatHeaderCopy"
249 PetscErrorCode MatHeaderCopy(Mat A,Mat C)
250 {
251   PetscErrorCode ierr;
252   PetscInt       refct;
253   PetscOps       *Abops;
254   MatOps         Aops;
255   char           *mtype,*mname;
256   void           *spptr;
257 
258   PetscFunctionBegin;
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;
264   mname = ((PetscObject)A)->name;
265   spptr = A->spptr;
266 
267   /* zero these so the destroy below does not free them */
268   ((PetscObject)A)->type_name = 0;
269   ((PetscObject)A)->name      = 0;
270 
271   /* free all the interior data structures from mat */
272   ierr = (*A->ops->destroy)(A);CHKERRQ(ierr);
273 
274   ierr = PetscFree(C->spptr);CHKERRQ(ierr);
275 
276   ierr = PetscMapDestroy(A->rmap);CHKERRQ(ierr);
277   ierr = PetscMapDestroy(A->cmap);CHKERRQ(ierr);
278   ierr = PetscFListDestroy(&((PetscObject)A)->qlist);CHKERRQ(ierr);
279   ierr = PetscOListDestroy(((PetscObject)A)->olist);CHKERRQ(ierr);
280 
281   /* copy C over to A */
282   ierr  = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr);
283 
284   /* return the parts of A we saved */
285   ((PetscObject)A)->bops      = Abops;
286   A->ops                      = Aops;
287   ((PetscObject)A)->refct     = refct;
288   ((PetscObject)A)->type_name = mtype;
289   ((PetscObject)A)->name      = mname;
290   A->spptr                    = spptr;
291 
292   /* since these two are copied into A we do not want them destroyed in C */
293   ((PetscObject)C)->qlist = 0;
294   ((PetscObject)C)->olist = 0;
295   ierr = PetscHeaderDestroy(C);CHKERRQ(ierr);
296   PetscFunctionReturn(0);
297 }
298 /*
299         Replace A's header with that of C
300         This is essentially code moved from MatDestroy
301 
302         This is somewhat different from MatHeaderCopy() it would be nice to merge the code
303 */
304 #undef __FUNCT__
305 #define __FUNCT__ "MatHeaderReplace"
306 PetscErrorCode MatHeaderReplace(Mat A,Mat C)
307 {
308   PetscErrorCode ierr;
309 
310   PetscFunctionBegin;
311   if (A == C) PetscFunctionReturn(0);
312 
313   /* free all the interior data structures from mat */
314   ierr = (*A->ops->destroy)(A);CHKERRQ(ierr);
315   ierr = PetscHeaderDestroy_Private((PetscObject)A);CHKERRQ(ierr);
316   ierr = PetscFree(A->ops);CHKERRQ(ierr);
317   ierr = PetscMapDestroy(A->rmap);CHKERRQ(ierr);
318   ierr = PetscMapDestroy(A->cmap);CHKERRQ(ierr);
319   ierr = PetscFree(A->spptr);CHKERRQ(ierr);
320 
321   /* copy C over to A */
322   if (C) {
323     ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr);
324     ierr = PetscLogObjectDestroy((PetscObject)C);CHKERRQ(ierr);
325     ierr = PetscFree(C);CHKERRQ(ierr);
326   }
327   PetscFunctionReturn(0);
328 }
329