Actual source code: gcreate.c
1: /*$Id: gcreate.c,v 1.131 2001/07/20 21:22:13 bsmith Exp $*/
3: #include src/mat/matimpl.h
4: #include petscsys.h
6: static int MatPublish_Base(PetscObject obj)
7: {
8: #if defined(PETSC_HAVE_AMS)
9: Mat mat = (Mat)obj;
11: #endif
14: #if defined(PETSC_HAVE_AMS)
15: /* if it is already published then return */
16: if (mat->amem >=0) return(0);
18: PetscObjectPublishBaseBegin(obj);
19: AMS_Memory_add_field((AMS_Memory)mat->amem,"globalsizes",&mat->M,2,AMS_INT,AMS_READ,
20: AMS_COMMON,AMS_REDUCT_UNDEF);
21: AMS_Memory_add_field((AMS_Memory)mat->amem,"localsizes",&mat->m,2,AMS_INT,AMS_READ,
22: AMS_DISTRIBUTED,AMS_REDUCT_UNDEF);
23: PetscObjectPublishBaseEnd(obj);
24: #endif
26: return(0);
27: }
30: /*@C
31: MatCreate - Creates a matrix where the type is determined
32: from either a call to MatSetType() or from the options database
33: with a call to MatSetFromOptions(). The default matrix type is
34: AIJ, using the routines MatCreateSeqAIJ() or MatCreateMPIAIJ()
35: if you do not set a type in the options database. If you never
36: call MatSetType() or MatSetFromOptions() it will generate an
37: error when you try to use the matrix.
39: Collective on MPI_Comm
41: Input Parameters:
42: + m - number of local rows (or PETSC_DECIDE)
43: . n - number of local columns (or PETSC_DECIDE)
44: . M - number of global rows (or PETSC_DETERMINE)
45: . N - number of global columns (or PETSC_DETERMINE)
46: - comm - MPI communicator
47:
48: Output Parameter:
49: . A - the matrix
51: Options Database Keys:
52: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
53: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
54: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
55: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
56: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
57: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
58: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
59: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
60: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
62: Even More Options Database Keys:
63: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
64: for additional format-specific options.
66: Notes:
67: If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
68: user must ensure that they are chosen to be compatible with the
69: vectors. To do this, one first considers the matrix-vector product
70: 'y = A x'. The 'm' that is used in the above routine must match the
71: local size used in the vector creation routine VecCreateMPI() for 'y'.
72: Likewise, the 'n' used must match that used as the local size in
73: VecCreateMPI() for 'x'.
75: Level: beginner
77: User manual sections:
78: + Section 3.1 Creating and Assembling Matrices
79: - Chapter 3 Matrices
81: .keywords: matrix, create
83: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
84: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
85: MatCreateSeqDense(), MatCreateMPIDense(),
86: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
87: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
88: MatConvert()
89: @*/
90: int MatCreate(MPI_Comm comm,int m,int n,int M,int N,Mat *A)
91: {
92: Mat B;
93: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
95: #endif
99: *A = PETSC_NULL;
100: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
101: MatInitializePackage(PETSC_NULL);
102: #endif
104: PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);
105: PetscLogObjectCreate(B);
107: B->m = m;
108: B->n = n;
109: B->M = M;
110: B->N = N;
112: B->preallocated = PETSC_FALSE;
113: B->bops->publish = MatPublish_Base;
114: *A = B;
115: return(0);
116: }
118: /*@C
119: MatSetFromOptions - Creates a matrix where the type is determined
120: from the options database. Generates a parallel MPI matrix if the
121: communicator has more than one processor. The default matrix type is
122: AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
123: you do not select a type in the options database.
125: Collective on Mat
127: Input Parameter:
128: . A - the matrix
130: Options Database Keys:
131: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
132: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
133: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
134: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
135: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
136: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
137: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
138: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
139: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
141: Even More Options Database Keys:
142: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
143: for additional format-specific options.
145: Level: beginner
147: .keywords: matrix, create
149: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
150: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
151: MatCreateSeqDense(), MatCreateMPIDense(),
152: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
153: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
154: MatConvert()
155: @*/
156: int MatSetFromOptions(Mat B)
157: {
158: int ierr,size;
159: char mtype[256];
160: PetscTruth flg;
163: PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
164: if (flg) {
165: MatSetType(B,mtype);
166: }
167: if (!B->type_name) {
168: MPI_Comm_size(B->comm,&size);
169: if (size == 1) {
170: MatSetType(B,MATSEQAIJ);
171: } else {
172: MatSetType(B,MATMPIAIJ);
173: }
174: }
175: #if defined(__cplusplus) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_MATSINGLE) && defined(PETSC_HAVE_CXX_NAMESPACE)
176: MatESISetFromOptions(B);
177: #endif
178: return(0);
179: }
181: /*@C
182: MatSetUpPreallocation
184: Collective on Mat
186: Input Parameter:
187: . A - the matrix
189: Level: beginner
191: .keywords: matrix, create
193: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
194: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
195: MatCreateSeqDense(), MatCreateMPIDense(),
196: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
197: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
198: MatConvert()
199: @*/
200: int MatSetUpPreallocation(Mat B)
201: {
202: int ierr;
205: if (B->ops->setuppreallocation) {
206: PetscLogInfo(B,"MatSetUpPreallocation: Warning not preallocating matrix storage");
207: (*B->ops->setuppreallocation)(B);
208: B->ops->setuppreallocation = 0;
209: B->preallocated = PETSC_TRUE;
210: }
211: return(0);
212: }
214: /*
215: Copies from Cs header to A
216: */
217: int MatHeaderCopy(Mat A,Mat C)
218: {
219: int ierr,refct;
220: PetscOps *Abops;
221: MatOps Aops;
222: char *mtype,*mname;
225: /* free all the interior data structures from mat */
226: (*A->ops->destroy)(A);
228: PetscMapDestroy(A->rmap);
229: PetscMapDestroy(A->cmap);
231: /* save the parts of A we need */
232: Abops = A->bops;
233: Aops = A->ops;
234: refct = A->refct;
235: mtype = A->type_name;
236: mname = A->name;
238: /* copy C over to A */
239: ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));
241: /* return the parts of A we saved */
242: A->bops = Abops;
243: A->ops = Aops;
244: A->qlist = 0;
245: A->refct = refct;
246: A->type_name = mtype;
247: A->name = mname;
249: PetscLogObjectDestroy(C);
250: PetscHeaderDestroy(C);
251: return(0);
252: }