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