xref: /petsc/src/mat/utils/gcreate.c (revision 300a7f5b1a2fd87c9550c4ae2ba6f40f37aef7e5)
1 #define PETSCMAT_DLL
2 
3 #include "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   PetscFunctionReturn(0);
138 }
139 
140 #undef __FUNCT__
141 #define __FUNCT__ "MatSetFromOptions"
142 /*@
143    MatSetFromOptions - Creates a matrix where the type is determined
144    from the options database. Generates a parallel MPI matrix if the
145    communicator has more than one processor.  The default matrix type is
146    AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
147    you do not select a type in the options database.
148 
149    Collective on Mat
150 
151    Input Parameter:
152 .  A - the matrix
153 
154    Options Database Keys:
155 +    -mat_type seqaij   - AIJ type, uses MatCreateSeqAIJ()
156 .    -mat_type mpiaij   - AIJ type, uses MatCreateMPIAIJ()
157 .    -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
158 .    -mat_type seqdense - dense type, uses MatCreateSeqDense()
159 .    -mat_type mpidense - dense type, uses MatCreateMPIDense()
160 .    -mat_type seqbaij  - block AIJ type, uses MatCreateSeqBAIJ()
161 -    -mat_type mpibaij  - block AIJ type, uses MatCreateMPIBAIJ()
162 
163    Even More Options Database Keys:
164    See the manpages for particular formats (e.g., MatCreateSeqAIJ())
165    for additional format-specific options.
166 
167    Level: beginner
168 
169 .keywords: matrix, create
170 
171 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
172           MatCreateSeqDense(), MatCreateMPIDense(),
173           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
174           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
175           MatConvert()
176 @*/
177 PetscErrorCode PETSCMAT_DLLEXPORT MatSetFromOptions(Mat B)
178 {
179   PetscErrorCode ierr;
180   const char     *deft = MATAIJ;
181   char           type[256];
182   PetscTruth     flg;
183 
184   PetscFunctionBegin;
185   PetscValidHeaderSpecific(B,MAT_COOKIE,1);
186 
187   ierr = PetscOptionsBegin(((PetscObject)B)->comm,((PetscObject)B)->prefix,"Matrix options","Mat");CHKERRQ(ierr);
188     ierr = PetscOptionsList("-mat_type","Matrix type","MatSetType",MatList,deft,type,256,&flg);CHKERRQ(ierr);
189     if (flg) {
190       ierr = MatSetType(B,type);CHKERRQ(ierr);
191     } else if (!((PetscObject)B)->type_name) {
192       ierr = MatSetType(B,deft);CHKERRQ(ierr);
193     }
194 
195     if (B->ops->setfromoptions) {
196       ierr = (*B->ops->setfromoptions)(B);CHKERRQ(ierr);
197     }
198 
199   ierr = PetscOptionsEnd();CHKERRQ(ierr);
200 
201   PetscFunctionReturn(0);
202 }
203 
204 #undef __FUNCT__
205 #define __FUNCT__ "MatSetUpPreallocation"
206 /*@
207    MatSetUpPreallocation
208 
209    Collective on Mat
210 
211    Input Parameter:
212 .  A - the matrix
213 
214    Level: beginner
215 
216 .keywords: matrix, create
217 
218 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
219           MatCreateSeqDense(), MatCreateMPIDense(),
220           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
221           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
222           MatConvert()
223 @*/
224 PetscErrorCode PETSCMAT_DLLEXPORT MatSetUpPreallocation(Mat B)
225 {
226   PetscErrorCode ierr;
227 
228   PetscFunctionBegin;
229   if (!B->preallocated && B->ops->setuppreallocation) {
230     ierr = PetscInfo(B,"Warning not preallocating matrix storage\n");CHKERRQ(ierr);
231     ierr = (*B->ops->setuppreallocation)(B);CHKERRQ(ierr);
232   }
233   B->preallocated = PETSC_TRUE;
234   PetscFunctionReturn(0);
235 }
236 
237 /*
238         Copies from Cs header to A
239 
240         This is somewhat different from MatHeaderReplace() it would be nice to merge the code
241 */
242 #undef __FUNCT__
243 #define __FUNCT__ "MatHeaderCopy"
244 PetscErrorCode MatHeaderCopy(Mat A,Mat C)
245 {
246   PetscErrorCode ierr;
247   PetscInt       refct;
248   PetscOps       *Abops;
249   MatOps         Aops;
250   char           *mtype,*mname;
251   void           *spptr;
252 
253   PetscFunctionBegin;
254   /* save the parts of A we need */
255   Abops = ((PetscObject)A)->bops;
256   Aops  = A->ops;
257   refct = ((PetscObject)A)->refct;
258   mtype = ((PetscObject)A)->type_name;
259   mname = ((PetscObject)A)->name;
260   spptr = A->spptr;
261 
262   /* zero these so the destroy below does not free them */
263   ((PetscObject)A)->type_name = 0;
264   ((PetscObject)A)->name      = 0;
265 
266   /* free all the interior data structures from mat */
267   ierr = (*A->ops->destroy)(A);CHKERRQ(ierr);
268 
269   ierr = PetscFree(C->spptr);CHKERRQ(ierr);
270 
271   ierr = PetscMapDestroy(A->rmap);CHKERRQ(ierr);
272   ierr = PetscMapDestroy(A->cmap);CHKERRQ(ierr);
273   ierr = PetscFListDestroy(&((PetscObject)A)->qlist);CHKERRQ(ierr);
274   ierr = PetscOListDestroy(((PetscObject)A)->olist);CHKERRQ(ierr);
275 
276   /* copy C over to A */
277   ierr  = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr);
278 
279   /* return the parts of A we saved */
280   ((PetscObject)A)->bops      = Abops;
281   A->ops                      = Aops;
282   ((PetscObject)A)->refct     = refct;
283   ((PetscObject)A)->type_name = mtype;
284   ((PetscObject)A)->name      = mname;
285   A->spptr                    = spptr;
286 
287   /* since these two are copied into A we do not want them destroyed in C */
288   ((PetscObject)C)->qlist = 0;
289   ((PetscObject)C)->olist = 0;
290   ierr = PetscHeaderDestroy(C);CHKERRQ(ierr);
291   PetscFunctionReturn(0);
292 }
293 /*
294         Replace A's header with that of C
295         This is essentially code moved from MatDestroy
296 
297         This is somewhat different from MatHeaderCopy() it would be nice to merge the code
298 */
299 #undef __FUNCT__
300 #define __FUNCT__ "MatHeaderReplace"
301 PetscErrorCode MatHeaderReplace(Mat A,Mat C)
302 {
303   PetscErrorCode ierr;
304 
305   PetscFunctionBegin;
306   /* free all the interior data structures from mat */
307   ierr = (*A->ops->destroy)(A);CHKERRQ(ierr);
308   ierr = PetscHeaderDestroy_Private((PetscObject)A);CHKERRQ(ierr);
309   ierr = PetscFree(A->ops);CHKERRQ(ierr);
310   ierr = PetscMapDestroy(A->rmap);CHKERRQ(ierr);
311   ierr = PetscMapDestroy(A->cmap);CHKERRQ(ierr);
312   ierr = PetscFree(A->spptr);CHKERRQ(ierr);
313 
314   /* copy C over to A */
315   if (C) {
316     ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr);
317     ierr = PetscLogObjectDestroy((PetscObject)C);CHKERRQ(ierr);
318     ierr = PetscFree(C);CHKERRQ(ierr);
319   }
320   PetscFunctionReturn(0);
321 }
322