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 = PetscLayoutCreate(comm,&B->rmap);CHKERRQ(ierr); 79 ierr = PetscLayoutCreate(comm,&B->cmap);CHKERRQ(ierr); 80 B->preallocated = PETSC_FALSE; 81 *A = B; 82 PetscFunctionReturn(0); 83 } 84 85 #undef __FUNCT__ 86 #define __FUNCT__ "MatSetSizes" 87 /*@ 88 MatSetSizes - Sets the local and global sizes, and checks to determine compatibility 89 90 Collective on Mat 91 92 Input Parameters: 93 + A - the matrix 94 . m - number of local rows (or PETSC_DECIDE) 95 . n - number of local columns (or PETSC_DECIDE) 96 . M - number of global rows (or PETSC_DETERMINE) 97 - N - number of global columns (or PETSC_DETERMINE) 98 99 Notes: 100 m (n) and M (N) cannot be both PETSC_DECIDE 101 If one processor calls this with M (N) of PETSC_DECIDE then all processors must, otherwise the program will hang. 102 103 If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the 104 user must ensure that they are chosen to be compatible with the 105 vectors. To do this, one first considers the matrix-vector product 106 'y = A x'. The 'm' that is used in the above routine must match the 107 local size used in the vector creation routine VecCreateMPI() for 'y'. 108 Likewise, the 'n' used must match that used as the local size in 109 VecCreateMPI() for 'x'. 110 111 Level: beginner 112 113 .seealso: MatGetSize(), PetscSplitOwnership() 114 @*/ 115 PetscErrorCode PETSCMAT_DLLEXPORT MatSetSizes(Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N) 116 { 117 PetscErrorCode ierr; 118 119 PetscFunctionBegin; 120 PetscValidHeaderSpecific(A,MAT_COOKIE,1); 121 if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M); 122 if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N); 123 if (A->ops->setsizes) { 124 /* Since this will not be set until the type has been set, this will NOT be called on the initial 125 call of MatSetSizes() (which must be called BEFORE MatSetType() */ 126 ierr = (*A->ops->setsizes)(A,m,n,M,N);CHKERRQ(ierr); 127 } else { 128 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); 129 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); 130 } 131 A->rmap->n = m; 132 A->cmap->n = n; 133 A->rmap->N = M; 134 A->cmap->N = N; 135 if (A->ops->create) { 136 ierr = (*A->ops->create)(A);CHKERRQ(ierr); 137 A->ops->create = 0; 138 } 139 140 PetscFunctionReturn(0); 141 } 142 143 #undef __FUNCT__ 144 #define __FUNCT__ "MatSetFromOptions" 145 /*@ 146 MatSetFromOptions - Creates a matrix where the type is determined 147 from the options database. Generates a parallel MPI matrix if the 148 communicator has more than one processor. The default matrix type is 149 AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if 150 you do not select a type in the options database. 151 152 Collective on Mat 153 154 Input Parameter: 155 . A - the matrix 156 157 Options Database Keys: 158 + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ() 159 . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ() 160 . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs() 161 . -mat_type seqdense - dense type, uses MatCreateSeqDense() 162 . -mat_type mpidense - dense type, uses MatCreateMPIDense() 163 . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ() 164 - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ() 165 166 Even More Options Database Keys: 167 See the manpages for particular formats (e.g., MatCreateSeqAIJ()) 168 for additional format-specific options. 169 170 Level: beginner 171 172 .keywords: matrix, create 173 174 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 175 MatCreateSeqDense(), MatCreateMPIDense(), 176 MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), 177 MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(), 178 MatConvert() 179 @*/ 180 PetscErrorCode PETSCMAT_DLLEXPORT MatSetFromOptions(Mat B) 181 { 182 PetscErrorCode ierr; 183 const char *deft = MATAIJ; 184 char type[256]; 185 PetscTruth flg; 186 187 PetscFunctionBegin; 188 PetscValidHeaderSpecific(B,MAT_COOKIE,1); 189 190 ierr = PetscOptionsBegin(((PetscObject)B)->comm,((PetscObject)B)->prefix,"Matrix options","Mat");CHKERRQ(ierr); 191 ierr = PetscOptionsList("-mat_type","Matrix type","MatSetType",MatList,deft,type,256,&flg);CHKERRQ(ierr); 192 if (flg) { 193 ierr = MatSetType(B,type);CHKERRQ(ierr); 194 } else if (!((PetscObject)B)->type_name) { 195 ierr = MatSetType(B,deft);CHKERRQ(ierr); 196 } 197 198 if (B->ops->setfromoptions) { 199 ierr = (*B->ops->setfromoptions)(B);CHKERRQ(ierr); 200 } 201 202 ierr = PetscOptionsEnd();CHKERRQ(ierr); 203 204 PetscFunctionReturn(0); 205 } 206 207 #undef __FUNCT__ 208 #define __FUNCT__ "MatSetUpPreallocation" 209 /*@ 210 MatSetUpPreallocation 211 212 Collective on Mat 213 214 Input Parameter: 215 . A - the matrix 216 217 Level: beginner 218 219 .keywords: matrix, create 220 221 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 222 MatCreateSeqDense(), MatCreateMPIDense(), 223 MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), 224 MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(), 225 MatConvert() 226 @*/ 227 PetscErrorCode PETSCMAT_DLLEXPORT MatSetUpPreallocation(Mat B) 228 { 229 PetscErrorCode ierr; 230 231 PetscFunctionBegin; 232 if (!B->preallocated && B->ops->setuppreallocation) { 233 ierr = PetscInfo(B,"Warning not preallocating matrix storage\n");CHKERRQ(ierr); 234 ierr = (*B->ops->setuppreallocation)(B);CHKERRQ(ierr); 235 } 236 B->preallocated = PETSC_TRUE; 237 PetscFunctionReturn(0); 238 } 239 240 /* 241 Copies from Cs header to A 242 243 This is somewhat different from MatHeaderReplace() it would be nice to merge the code 244 */ 245 #undef __FUNCT__ 246 #define __FUNCT__ "MatHeaderCopy" 247 PetscErrorCode MatHeaderCopy(Mat A,Mat C) 248 { 249 PetscErrorCode ierr; 250 PetscInt refct; 251 PetscOps *Abops; 252 MatOps Aops; 253 char *mtype,*mname; 254 void *spptr; 255 256 PetscFunctionBegin; 257 /* save the parts of A we need */ 258 Abops = ((PetscObject)A)->bops; 259 Aops = A->ops; 260 refct = ((PetscObject)A)->refct; 261 mtype = ((PetscObject)A)->type_name; 262 mname = ((PetscObject)A)->name; 263 spptr = A->spptr; 264 265 /* zero these so the destroy below does not free them */ 266 ((PetscObject)A)->type_name = 0; 267 ((PetscObject)A)->name = 0; 268 269 /* free all the interior data structures from mat */ 270 ierr = (*A->ops->destroy)(A);CHKERRQ(ierr); 271 272 ierr = PetscFree(C->spptr);CHKERRQ(ierr); 273 274 ierr = PetscLayoutDestroy(A->rmap);CHKERRQ(ierr); 275 ierr = PetscLayoutDestroy(A->cmap);CHKERRQ(ierr); 276 ierr = PetscFListDestroy(&((PetscObject)A)->qlist);CHKERRQ(ierr); 277 ierr = PetscOListDestroy(((PetscObject)A)->olist);CHKERRQ(ierr); 278 279 /* copy C over to A */ 280 ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr); 281 282 /* return the parts of A we saved */ 283 ((PetscObject)A)->bops = Abops; 284 A->ops = Aops; 285 ((PetscObject)A)->refct = refct; 286 ((PetscObject)A)->type_name = mtype; 287 ((PetscObject)A)->name = mname; 288 A->spptr = spptr; 289 290 /* since these two are copied into A we do not want them destroyed in C */ 291 ((PetscObject)C)->qlist = 0; 292 ((PetscObject)C)->olist = 0; 293 ierr = PetscHeaderDestroy(C);CHKERRQ(ierr); 294 PetscFunctionReturn(0); 295 } 296 /* 297 Replace A's header with that of C 298 This is essentially code moved from MatDestroy 299 300 This is somewhat different from MatHeaderCopy() it would be nice to merge the code 301 */ 302 #undef __FUNCT__ 303 #define __FUNCT__ "MatHeaderReplace" 304 PetscErrorCode MatHeaderReplace(Mat A,Mat C) 305 { 306 PetscErrorCode ierr; 307 308 PetscFunctionBegin; 309 if (A == C) PetscFunctionReturn(0); 310 311 /* free all the interior data structures from mat */ 312 ierr = (*A->ops->destroy)(A);CHKERRQ(ierr); 313 ierr = PetscHeaderDestroy_Private((PetscObject)A);CHKERRQ(ierr); 314 ierr = PetscFree(A->ops);CHKERRQ(ierr); 315 ierr = PetscLayoutDestroy(A->rmap);CHKERRQ(ierr); 316 ierr = PetscLayoutDestroy(A->cmap);CHKERRQ(ierr); 317 ierr = PetscFree(A->spptr);CHKERRQ(ierr); 318 319 /* copy C over to A */ 320 if (C) { 321 ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr); 322 ierr = PetscLogObjectDestroy((PetscObject)C);CHKERRQ(ierr); 323 ierr = PetscFree(C);CHKERRQ(ierr); 324 } 325 PetscFunctionReturn(0); 326 } 327