1 2 #include "src/mat/matimpl.h" /*I "petscmat.h" I*/ 3 #include "petscsys.h" 4 5 #undef __FUNCT__ 6 #define __FUNCT__ "MatPublish_Base" 7 static PetscErrorCode MatPublish_Base(PetscObject obj) 8 { 9 PetscFunctionBegin; 10 PetscFunctionReturn(0); 11 } 12 13 14 #undef __FUNCT__ 15 #define __FUNCT__ "MatCreate" 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. 24 25 Collective on MPI_Comm 26 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 36 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() 47 48 Even More Options Database Keys: 49 See the manpages for particular formats (e.g., MatCreateSeqAIJ()) 50 for additional format-specific options. 51 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'. 60 61 Level: beginner 62 63 User manual sections: 64 + sec_matcreate 65 - chapter_matrices 66 67 .keywords: matrix, create 68 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 PetscErrorCode ierr; 80 81 PetscFunctionBegin; 82 PetscValidPointer(A,6); 83 if (M > 0 && m > M) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local column size %D cannot be larger than global column size %D",m,M); 84 if (N > 0 && n > N) SETERRQ2(PETSC_ERR_ARG_INCOMP,"Local row size %D cannot be larger than global row size %D",n,N); 85 86 *A = PETSC_NULL; 87 #ifndef PETSC_USE_DYNAMIC_LIBRARIES 88 ierr = MatInitializePackage(PETSC_NULL);CHKERRQ(ierr); 89 #endif 90 91 ierr = PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);CHKERRQ(ierr); 92 93 B->m = m; 94 B->n = n; 95 B->M = M; 96 B->N = N; 97 B->bs = 1; 98 B->preallocated = PETSC_FALSE; 99 B->bops->publish = MatPublish_Base; 100 *A = B; 101 PetscFunctionReturn(0); 102 } 103 104 #undef __FUNCT__ 105 #define __FUNCT__ "MatSetFromOptions" 106 /*@C 107 MatSetFromOptions - Creates a matrix where the type is determined 108 from the options database. Generates a parallel MPI matrix if the 109 communicator has more than one processor. The default matrix type is 110 AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ() if 111 you do not select a type in the options database. 112 113 Collective on Mat 114 115 Input Parameter: 116 . A - the matrix 117 118 Options Database Keys: 119 + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ() 120 . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ() 121 . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag() 122 . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag() 123 . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs() 124 . -mat_type seqdense - dense type, uses MatCreateSeqDense() 125 . -mat_type mpidense - dense type, uses MatCreateMPIDense() 126 . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ() 127 - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ() 128 129 Even More Options Database Keys: 130 See the manpages for particular formats (e.g., MatCreateSeqAIJ()) 131 for additional format-specific options. 132 133 Level: beginner 134 135 .keywords: matrix, create 136 137 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 138 MatCreateSeqBDiag(),MatCreateMPIBDiag(), 139 MatCreateSeqDense(), MatCreateMPIDense(), 140 MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), 141 MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(), 142 MatConvert() 143 @*/ 144 PetscErrorCode MatSetFromOptions(Mat B) 145 { 146 PetscErrorCode ierr; 147 char mtype[256]; 148 PetscTruth flg; 149 150 PetscFunctionBegin; 151 ierr = PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);CHKERRQ(ierr); 152 if (flg) { 153 ierr = MatSetType(B,mtype);CHKERRQ(ierr); 154 } 155 if (!B->type_name) { 156 ierr = MatSetType(B,MATAIJ);CHKERRQ(ierr); 157 } 158 PetscFunctionReturn(0); 159 } 160 161 #undef __FUNCT__ 162 #define __FUNCT__ "MatSetUpPreallocation" 163 /*@C 164 MatSetUpPreallocation 165 166 Collective on Mat 167 168 Input Parameter: 169 . A - the matrix 170 171 Level: beginner 172 173 .keywords: matrix, create 174 175 .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(), 176 MatCreateSeqBDiag(),MatCreateMPIBDiag(), 177 MatCreateSeqDense(), MatCreateMPIDense(), 178 MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(), 179 MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(), 180 MatConvert() 181 @*/ 182 PetscErrorCode MatSetUpPreallocation(Mat B) 183 { 184 PetscErrorCode ierr; 185 186 PetscFunctionBegin; 187 if (B->ops->setuppreallocation) { 188 PetscLogInfo(B,"MatSetUpPreallocation: Warning not preallocating matrix storage"); 189 ierr = (*B->ops->setuppreallocation)(B);CHKERRQ(ierr); 190 B->ops->setuppreallocation = 0; 191 B->preallocated = PETSC_TRUE; 192 } 193 PetscFunctionReturn(0); 194 } 195 196 /* 197 Copies from Cs header to A 198 */ 199 #undef __FUNCT__ 200 #define __FUNCT__ "MatHeaderCopy" 201 PetscErrorCode MatHeaderCopy(Mat A,Mat C) 202 { 203 PetscErrorCode ierr; 204 PetscInt refct; 205 PetscOps *Abops; 206 MatOps Aops; 207 char *mtype,*mname; 208 void *spptr; 209 210 PetscFunctionBegin; 211 /* free all the interior data structures from mat */ 212 ierr = (*A->ops->destroy)(A);CHKERRQ(ierr); 213 214 ierr = PetscMapDestroy(A->rmap);CHKERRQ(ierr); 215 ierr = PetscMapDestroy(A->cmap);CHKERRQ(ierr); 216 217 /* save the parts of A we need */ 218 Abops = A->bops; 219 Aops = A->ops; 220 refct = A->refct; 221 mtype = A->type_name; 222 mname = A->name; 223 spptr = A->spptr; 224 225 if (C->spptr) { 226 ierr = PetscFree(C->spptr);CHKERRQ(ierr); 227 C->spptr = PETSC_NULL; 228 } 229 230 /* copy C over to A */ 231 ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr); 232 233 /* return the parts of A we saved */ 234 A->bops = Abops; 235 A->ops = Aops; 236 A->qlist = 0; 237 A->refct = refct; 238 A->type_name = mtype; 239 A->name = mname; 240 A->spptr = spptr; 241 242 ierr = PetscHeaderDestroy(C);CHKERRQ(ierr); 243 PetscFunctionReturn(0); 244 } 245 /* 246 Replace A's header with that of C 247 This is essentially code moved from MatDestroy 248 */ 249 #undef __FUNCT__ 250 #define __FUNCT__ "MatHeaderReplace" 251 PetscErrorCode MatHeaderReplace(Mat A,Mat C) 252 { 253 PetscErrorCode ierr; 254 255 PetscFunctionBegin; 256 /* free all the interior data structures from mat */ 257 ierr = (*A->ops->destroy)(A);CHKERRQ(ierr); 258 ierr = PetscHeaderDestroy_Private((PetscObject)A);CHKERRQ(ierr); 259 260 /* copy C over to A */ 261 if (C) { 262 ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));CHKERRQ(ierr); 263 ierr = PetscLogObjectDestroy((PetscObject)C);CHKERRQ(ierr); 264 ierr = PetscFree(C);CHKERRQ(ierr); 265 } 266 PetscFunctionReturn(0); 267 } 268