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