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