Actual source code: gcreate.c

 2:  #include src/mat/matimpl.h
 3:  #include petscsys.h

  7: static PetscErrorCode MatPublish_Base(PetscObject obj)
  8: {
 10:   return(0);
 11: }


 16: /*@C
 17:    MatCreate - Creates a matrix where the type is determined
 18:    from either a call to MatSetType() or from the options database
 19:    with a call to MatSetFromOptions(). The default matrix type is
 20:    AIJ, using the routines MatCreateSeqAIJ() or MatCreateMPIAIJ()
 21:    if you do not set a type in the options database. If you never
 22:    call MatSetType() or MatSetFromOptions() it will generate an 
 23:    error when you try to use the matrix.

 25:    Collective on MPI_Comm

 27:    Input Parameters:
 28: +  m - number of local rows (or PETSC_DECIDE)
 29: .  n - number of local columns (or PETSC_DECIDE)
 30: .  M - number of global rows (or PETSC_DETERMINE)
 31: .  N - number of global columns (or PETSC_DETERMINE)
 32: -  comm - MPI communicator
 33:  
 34:    Output Parameter:
 35: .  A - the matrix

 37:    Options Database Keys:
 38: +    -mat_type seqaij   - AIJ type, uses MatCreateSeqAIJ()
 39: .    -mat_type mpiaij   - AIJ type, uses MatCreateMPIAIJ()
 40: .    -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
 41: .    -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
 42: .    -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
 43: .    -mat_type seqdense - dense type, uses MatCreateSeqDense()
 44: .    -mat_type mpidense - dense type, uses MatCreateMPIDense()
 45: .    -mat_type seqbaij  - block AIJ type, uses MatCreateSeqBAIJ()
 46: -    -mat_type mpibaij  - block AIJ type, uses MatCreateMPIBAIJ()

 48:    Even More Options Database Keys:
 49:    See the manpages for particular formats (e.g., MatCreateSeqAIJ())
 50:    for additional format-specific options.

 52:    Notes:
 53:    If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
 54:    user must ensure that they are chosen to be compatible with the
 55:    vectors. To do this, one first considers the matrix-vector product 
 56:    'y = A x'. The 'm' that is used in the above routine must match the 
 57:    local size used in the vector creation routine VecCreateMPI() for 'y'.
 58:    Likewise, the 'n' used must match that used as the local size in
 59:    VecCreateMPI() for 'x'.

 61:    Level: beginner

 63:    User manual sections:
 64: +   Section 3.1 Creating and Assembling Matrices
 65: -   Chapter 3 Matrices

 67: .keywords: matrix, create

 69: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 
 70:           MatCreateSeqBDiag(),MatCreateMPIBDiag(),
 71:           MatCreateSeqDense(), MatCreateMPIDense(), 
 72:           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
 73:           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
 74:           MatConvert()
 75: @*/
 76: PetscErrorCode MatCreate(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt M,PetscInt N,Mat *A)
 77: {
 78:   Mat            B;
 79: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
 81: #endif

 85:   if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M);
 86:   if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N);

 88:   *A = PETSC_NULL;
 89: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
 90:   MatInitializePackage(PETSC_NULL);
 91: #endif

 93:   PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);
 94:   PetscLogObjectCreate(B);

 96:   B->m             = m;
 97:   B->n             = n;
 98:   B->M             = M;
 99:   B->N             = N;
100:   B->bs            = 1;
101:   B->preallocated  = PETSC_FALSE;
102:   B->bops->publish = MatPublish_Base;
103:   *A               = B;
104:   return(0);
105: }

109: /*@C
110:    MatSetFromOptions - Creates a matrix where the type is determined
111:    from the options database. Generates a parallel MPI matrix if the
112:    communicator has more than one processor.  The default matrix type is
113:    AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if
114:    you do not select a type in the options database.

116:    Collective on Mat

118:    Input Parameter:
119: .  A - the matrix

121:    Options Database Keys:
122: +    -mat_type seqaij   - AIJ type, uses MatCreateSeqAIJ()
123: .    -mat_type mpiaij   - AIJ type, uses MatCreateMPIAIJ()
124: .    -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
125: .    -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
126: .    -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
127: .    -mat_type seqdense - dense type, uses MatCreateSeqDense()
128: .    -mat_type mpidense - dense type, uses MatCreateMPIDense()
129: .    -mat_type seqbaij  - block AIJ type, uses MatCreateSeqBAIJ()
130: -    -mat_type mpibaij  - block AIJ type, uses MatCreateMPIBAIJ()

132:    Even More Options Database Keys:
133:    See the manpages for particular formats (e.g., MatCreateSeqAIJ())
134:    for additional format-specific options.

136:    Level: beginner

138: .keywords: matrix, create

140: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 
141:           MatCreateSeqBDiag(),MatCreateMPIBDiag(),
142:           MatCreateSeqDense(), MatCreateMPIDense(), 
143:           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
144:           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
145:           MatConvert()
146: @*/
147: PetscErrorCode MatSetFromOptions(Mat B)
148: {
150:   PetscMPIInt    size;
151:   char           mtype[256];
152:   PetscTruth     flg;

155:   PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
156:   if (flg) {
157:     MatSetType(B,mtype);
158:   }
159:   if (!B->type_name) {
160:     MPI_Comm_size(B->comm,&size);
161:     MatSetType(B,MATAIJ);
162:   }
163:   return(0);
164: }

168: /*@C
169:    MatSetUpPreallocation

171:    Collective on Mat

173:    Input Parameter:
174: .  A - the matrix

176:    Level: beginner

178: .keywords: matrix, create

180: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 
181:           MatCreateSeqBDiag(),MatCreateMPIBDiag(),
182:           MatCreateSeqDense(), MatCreateMPIDense(), 
183:           MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
184:           MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
185:           MatConvert()
186: @*/
187: PetscErrorCode MatSetUpPreallocation(Mat B)
188: {

192:   if (B->ops->setuppreallocation) {
193:     PetscLogInfo(B,"MatSetUpPreallocation: Warning not preallocating matrix storage");
194:     (*B->ops->setuppreallocation)(B);
195:     B->ops->setuppreallocation = 0;
196:     B->preallocated            = PETSC_TRUE;
197:   }
198:   return(0);
199: }

201: /*
202:         Copies from Cs header to A
203: */
206: PetscErrorCode MatHeaderCopy(Mat A,Mat C)
207: {
209:   PetscInt       refct;
210:   PetscOps       *Abops;
211:   MatOps         Aops;
212:   char           *mtype,*mname;

215:   /* free all the interior data structures from mat */
216:   (*A->ops->destroy)(A);

218:   PetscMapDestroy(A->rmap);
219:   PetscMapDestroy(A->cmap);

221:   /* save the parts of A we need */
222:   Abops = A->bops;
223:   Aops  = A->ops;
224:   refct = A->refct;
225:   mtype = A->type_name;
226:   mname = A->name;

228:   /* copy C over to A */
229:   PetscMemcpy(A,C,sizeof(struct _p_Mat));

231:   /* return the parts of A we saved */
232:   A->bops      = Abops;
233:   A->ops       = Aops;
234:   A->qlist     = 0;
235:   A->refct     = refct;
236:   A->type_name = mtype;
237:   A->name      = mname;

239:   PetscLogObjectDestroy(C);
240:   PetscHeaderDestroy(C);
241:   return(0);
242: }