1 #include <petscdmshell.h> /*I "petscdmshell.h" I*/ 2 #include <petscmat.h> 3 #include <petsc-private/dmimpl.h> 4 5 typedef struct { 6 Vec Xglobal; 7 Mat A; 8 } DM_Shell; 9 10 #undef __FUNCT__ 11 #define __FUNCT__ "DMCreateMatrix_Shell" 12 static PetscErrorCode DMCreateMatrix_Shell(DM dm,MatType mtype,Mat *J) 13 { 14 PetscErrorCode ierr; 15 DM_Shell *shell = (DM_Shell*)dm->data; 16 Mat A; 17 18 PetscFunctionBegin; 19 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 20 PetscValidPointer(J,3); 21 if (!shell->A) { 22 if (shell->Xglobal) { 23 PetscInt m,M; 24 ierr = PetscInfo(dm,"Naively creating matrix using global vector distribution without preallocation");CHKERRQ(ierr); 25 ierr = VecGetSize(shell->Xglobal,&M);CHKERRQ(ierr); 26 ierr = VecGetLocalSize(shell->Xglobal,&m);CHKERRQ(ierr); 27 ierr = MatCreate(((PetscObject)dm)->comm,&shell->A);CHKERRQ(ierr); 28 ierr = MatSetSizes(shell->A,m,m,M,M);CHKERRQ(ierr); 29 if (mtype) {ierr = MatSetType(shell->A,mtype);CHKERRQ(ierr);} 30 ierr = MatSetUp(shell->A);CHKERRQ(ierr); 31 } else SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_USER,"Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector"); 32 } 33 A = shell->A; 34 /* the check below is tacky and incomplete */ 35 if (mtype) { 36 PetscBool flg,aij,seqaij,mpiaij; 37 ierr = PetscObjectTypeCompare((PetscObject)A,mtype,&flg);CHKERRQ(ierr); 38 ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&seqaij);CHKERRQ(ierr); 39 ierr = PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&mpiaij);CHKERRQ(ierr); 40 ierr = PetscStrcmp(mtype,MATAIJ,&aij);CHKERRQ(ierr); 41 if (!flg) { 42 if (!(aij & (seqaij || mpiaij))) SETERRQ2(((PetscObject)dm)->comm,PETSC_ERR_ARG_NOTSAMETYPE,"Requested matrix of type %s, but only %s available",mtype,((PetscObject)A)->type_name); 43 } 44 } 45 if (((PetscObject)A)->refct < 2) { /* We have an exclusive reference so we can give it out */ 46 ierr = PetscObjectReference((PetscObject)A);CHKERRQ(ierr); 47 ierr = MatZeroEntries(A);CHKERRQ(ierr); 48 *J = A; 49 } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */ 50 ierr = MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,J);CHKERRQ(ierr); 51 ierr = MatZeroEntries(*J);CHKERRQ(ierr); 52 } 53 PetscFunctionReturn(0); 54 } 55 56 #undef __FUNCT__ 57 #define __FUNCT__ "DMCreateGlobalVector_Shell" 58 PetscErrorCode DMCreateGlobalVector_Shell(DM dm,Vec *gvec) 59 { 60 PetscErrorCode ierr; 61 DM_Shell *shell = (DM_Shell*)dm->data; 62 Vec X; 63 64 PetscFunctionBegin; 65 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 66 PetscValidPointer(gvec,2); 67 *gvec = 0; 68 X = shell->Xglobal; 69 if (!X) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_USER,"Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()"); 70 if (((PetscObject)X)->refct < 2) { /* We have an exclusive reference so we can give it out */ 71 ierr = PetscObjectReference((PetscObject)X);CHKERRQ(ierr); 72 ierr = VecZeroEntries(X);CHKERRQ(ierr); 73 *gvec = X; 74 } else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */ 75 ierr = VecDuplicate(X,gvec);CHKERRQ(ierr); 76 ierr = VecZeroEntries(*gvec);CHKERRQ(ierr); 77 } 78 ierr = VecSetDM(*gvec,dm);CHKERRQ(ierr); 79 PetscFunctionReturn(0); 80 } 81 82 #undef __FUNCT__ 83 #define __FUNCT__ "DMShellSetMatrix" 84 /*@ 85 DMShellSetMatrix - sets a template matrix associated with the DMShell 86 87 Collective 88 89 Input Arguments: 90 + dm - shell DM 91 - J - template matrix 92 93 Level: advanced 94 95 .seealso: DMCreateMatrix(), DMShellSetCreateMatrix() 96 @*/ 97 PetscErrorCode DMShellSetMatrix(DM dm,Mat J) 98 { 99 DM_Shell *shell = (DM_Shell*)dm->data; 100 PetscErrorCode ierr; 101 PetscBool isshell; 102 103 PetscFunctionBegin; 104 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 105 PetscValidHeaderSpecific(J,MAT_CLASSID,2); 106 ierr = PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);CHKERRQ(ierr); 107 if (!isshell) PetscFunctionReturn(0); 108 ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr); 109 ierr = MatDestroy(&shell->A);CHKERRQ(ierr); 110 shell->A = J; 111 PetscFunctionReturn(0); 112 } 113 114 #undef __FUNCT__ 115 #define __FUNCT__ "DMShellSetCreateMatrix" 116 /*@C 117 DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM 118 119 Logically Collective on DM 120 121 Input Arguments: 122 + dm - the shell DM 123 - func - the function to create a matrix 124 125 Level: advanced 126 127 .seealso: DMCreateMatrix(), DMShellSetMatrix() 128 @*/ 129 PetscErrorCode DMShellSetCreateMatrix(DM dm,PetscErrorCode (*func)(DM,MatType,Mat*)) 130 { 131 132 PetscFunctionBegin; 133 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 134 dm->ops->creatematrix = func; 135 PetscFunctionReturn(0); 136 } 137 138 #undef __FUNCT__ 139 #define __FUNCT__ "DMShellSetGlobalVector" 140 /*@ 141 DMShellSetGlobalVector - sets a template global vector associated with the DMShell 142 143 Logically Collective on DM 144 145 Input Arguments: 146 + dm - shell DM 147 - X - template vector 148 149 Level: advanced 150 151 .seealso: DMCreateGlobalVector(), DMShellSetMatrix(), DMShellSetCreateGlobalVector() 152 @*/ 153 PetscErrorCode DMShellSetGlobalVector(DM dm,Vec X) 154 { 155 DM_Shell *shell = (DM_Shell*)dm->data; 156 PetscErrorCode ierr; 157 PetscBool isshell; 158 159 PetscFunctionBegin; 160 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 161 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 162 ierr = PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);CHKERRQ(ierr); 163 if (!isshell) PetscFunctionReturn(0); 164 ierr = PetscObjectReference((PetscObject)X);CHKERRQ(ierr); 165 ierr = VecDestroy(&shell->Xglobal);CHKERRQ(ierr); 166 shell->Xglobal = X; 167 PetscFunctionReturn(0); 168 } 169 170 #undef __FUNCT__ 171 #define __FUNCT__ "DMShellSetCreateGlobalVector" 172 /*@C 173 DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM 174 175 Logically Collective 176 177 Input Arguments: 178 + dm - the shell DM 179 - func - the creation routine 180 181 Level: advanced 182 183 .seealso: DMShellSetGlobalVector(), DMShellSetCreateMatrix() 184 @*/ 185 PetscErrorCode DMShellSetCreateGlobalVector(DM dm,PetscErrorCode (*func)(DM,Vec*)) 186 { 187 188 PetscFunctionBegin; 189 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 190 dm->ops->createglobalvector = func; 191 PetscFunctionReturn(0); 192 } 193 194 #undef __FUNCT__ 195 #define __FUNCT__ "DMDestroy_Shell" 196 static PetscErrorCode DMDestroy_Shell(DM dm) 197 { 198 PetscErrorCode ierr; 199 DM_Shell *shell = (DM_Shell*)dm->data; 200 201 PetscFunctionBegin; 202 ierr = MatDestroy(&shell->A);CHKERRQ(ierr); 203 ierr = VecDestroy(&shell->Xglobal);CHKERRQ(ierr); 204 /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */ 205 ierr = PetscFree(shell);CHKERRQ(ierr); 206 PetscFunctionReturn(0); 207 } 208 209 #undef __FUNCT__ 210 #define __FUNCT__ "DMView_Shell" 211 static PetscErrorCode DMView_Shell(DM dm,PetscViewer v) 212 { 213 PetscErrorCode ierr; 214 DM_Shell *shell = (DM_Shell*)dm->data; 215 216 PetscFunctionBegin; 217 ierr = VecView(shell->Xglobal,v);CHKERRQ(ierr); 218 PetscFunctionReturn(0); 219 } 220 221 #undef __FUNCT__ 222 #define __FUNCT__ "DMLoad_Shell" 223 static PetscErrorCode DMLoad_Shell(DM dm,PetscViewer v) 224 { 225 PetscErrorCode ierr; 226 DM_Shell *shell = (DM_Shell*)dm->data; 227 228 PetscFunctionBegin; 229 ierr = VecCreate(((PetscObject)dm)->comm,&shell->Xglobal);CHKERRQ(ierr); 230 ierr = VecLoad(shell->Xglobal,v);CHKERRQ(ierr); 231 PetscFunctionReturn(0); 232 } 233 234 #undef __FUNCT__ 235 #define __FUNCT__ "DMCreate_Shell" 236 PETSC_EXTERN_C PetscErrorCode DMCreate_Shell(DM dm) 237 { 238 PetscErrorCode ierr; 239 DM_Shell *shell; 240 241 PetscFunctionBegin; 242 ierr = PetscNewLog(dm,DM_Shell,&shell);CHKERRQ(ierr); 243 dm->data = shell; 244 245 ierr = PetscObjectChangeTypeName((PetscObject)dm,DMSHELL);CHKERRQ(ierr); 246 dm->ops->destroy = DMDestroy_Shell; 247 dm->ops->createglobalvector = DMCreateGlobalVector_Shell; 248 dm->ops->creatematrix = DMCreateMatrix_Shell; 249 dm->ops->view = DMView_Shell; 250 dm->ops->load = DMLoad_Shell; 251 PetscFunctionReturn(0); 252 } 253 254 #undef __FUNCT__ 255 #define __FUNCT__ "DMShellCreate" 256 /*@ 257 DMShellCreate - Creates a shell DM object, used to manage user-defined problem data 258 259 Collective on MPI_Comm 260 261 Input Parameter: 262 . comm - the processors that will share the global vector 263 264 Output Parameters: 265 . shell - the shell DM 266 267 Level: advanced 268 269 .seealso DMDestroy(), DMCreateGlobalVector() 270 @*/ 271 PetscErrorCode DMShellCreate(MPI_Comm comm,DM *dm) 272 { 273 PetscErrorCode ierr; 274 275 PetscFunctionBegin; 276 PetscValidPointer(dm,2); 277 ierr = DMCreate(comm,dm);CHKERRQ(ierr); 278 ierr = DMSetType(*dm,DMSHELL);CHKERRQ(ierr); 279 PetscFunctionReturn(0); 280 } 281