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(PetscObjectComm((PetscObject)dm),&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(PetscObjectComm((PetscObject)dm),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(PetscObjectComm((PetscObject)dm),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(PetscObjectComm((PetscObject)dm),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(PetscObjectComm((PetscObject)dm),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 ierr = VecSetDM(*gvec,dm);CHKERRQ(ierr); 106 PetscFunctionReturn(0); 107 } 108 109 #undef __FUNCT__ 110 #define __FUNCT__ "DMShellSetMatrix" 111 /*@ 112 DMShellSetMatrix - sets a template matrix associated with the DMShell 113 114 Collective 115 116 Input Arguments: 117 + dm - shell DM 118 - J - template matrix 119 120 Level: advanced 121 122 .seealso: DMCreateMatrix(), DMShellSetCreateMatrix() 123 @*/ 124 PetscErrorCode DMShellSetMatrix(DM dm,Mat J) 125 { 126 DM_Shell *shell = (DM_Shell*)dm->data; 127 PetscErrorCode ierr; 128 PetscBool isshell; 129 130 PetscFunctionBegin; 131 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 132 PetscValidHeaderSpecific(J,MAT_CLASSID,2); 133 ierr = PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);CHKERRQ(ierr); 134 if (!isshell) PetscFunctionReturn(0); 135 ierr = PetscObjectReference((PetscObject)J);CHKERRQ(ierr); 136 ierr = MatDestroy(&shell->A);CHKERRQ(ierr); 137 shell->A = J; 138 PetscFunctionReturn(0); 139 } 140 141 #undef __FUNCT__ 142 #define __FUNCT__ "DMShellSetCreateMatrix" 143 /*@C 144 DMShellSetCreateMatrix - sets the routine to create a matrix associated with the shell DM 145 146 Logically Collective on DM 147 148 Input Arguments: 149 + dm - the shell DM 150 - func - the function to create a matrix 151 152 Level: advanced 153 154 .seealso: DMCreateMatrix(), DMShellSetMatrix() 155 @*/ 156 PetscErrorCode DMShellSetCreateMatrix(DM dm,PetscErrorCode (*func)(DM,MatType,Mat*)) 157 { 158 159 PetscFunctionBegin; 160 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 161 dm->ops->creatematrix = func; 162 PetscFunctionReturn(0); 163 } 164 165 #undef __FUNCT__ 166 #define __FUNCT__ "DMShellSetGlobalVector" 167 /*@ 168 DMShellSetGlobalVector - sets a template global vector associated with the DMShell 169 170 Logically Collective on DM 171 172 Input Arguments: 173 + dm - shell DM 174 - X - template vector 175 176 Level: advanced 177 178 .seealso: DMCreateGlobalVector(), DMShellSetMatrix(), DMShellSetCreateGlobalVector() 179 @*/ 180 PetscErrorCode DMShellSetGlobalVector(DM dm,Vec X) 181 { 182 DM_Shell *shell = (DM_Shell*)dm->data; 183 PetscErrorCode ierr; 184 PetscBool isshell; 185 186 PetscFunctionBegin; 187 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 188 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 189 ierr = PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);CHKERRQ(ierr); 190 if (!isshell) PetscFunctionReturn(0); 191 ierr = PetscObjectReference((PetscObject)X);CHKERRQ(ierr); 192 ierr = VecDestroy(&shell->Xglobal);CHKERRQ(ierr); 193 shell->Xglobal = X; 194 PetscFunctionReturn(0); 195 } 196 197 #undef __FUNCT__ 198 #define __FUNCT__ "DMShellSetCreateGlobalVector" 199 /*@C 200 DMShellSetCreateGlobalVector - sets the routine to create a global vector associated with the shell DM 201 202 Logically Collective 203 204 Input Arguments: 205 + dm - the shell DM 206 - func - the creation routine 207 208 Level: advanced 209 210 .seealso: DMShellSetGlobalVector(), DMShellSetCreateMatrix() 211 @*/ 212 PetscErrorCode DMShellSetCreateGlobalVector(DM dm,PetscErrorCode (*func)(DM,Vec*)) 213 { 214 215 PetscFunctionBegin; 216 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 217 dm->ops->createglobalvector = func; 218 PetscFunctionReturn(0); 219 } 220 221 #undef __FUNCT__ 222 #define __FUNCT__ "DMShellSetLocalVector" 223 /*@ 224 DMShellSetLocalVector - sets a template local vector associated with the DMShell 225 226 Logically Collective on DM 227 228 Input Arguments: 229 + dm - shell DM 230 - X - template vector 231 232 Level: advanced 233 234 .seealso: DMCreateLocalVector(), DMShellSetMatrix(), DMShellSetCreateLocalVector() 235 @*/ 236 PetscErrorCode DMShellSetLocalVector(DM dm,Vec X) 237 { 238 DM_Shell *shell = (DM_Shell*)dm->data; 239 PetscErrorCode ierr; 240 PetscBool isshell; 241 242 PetscFunctionBegin; 243 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 244 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 245 ierr = PetscObjectTypeCompare((PetscObject)dm,DMSHELL,&isshell);CHKERRQ(ierr); 246 if (!isshell) PetscFunctionReturn(0); 247 ierr = PetscObjectReference((PetscObject)X);CHKERRQ(ierr); 248 ierr = VecDestroy(&shell->Xlocal);CHKERRQ(ierr); 249 shell->Xlocal = X; 250 PetscFunctionReturn(0); 251 } 252 253 #undef __FUNCT__ 254 #define __FUNCT__ "DMShellSetCreateLocalVector" 255 /*@C 256 DMShellSetCreateLocalVector - sets the routine to create a local vector associated with the shell DM 257 258 Logically Collective 259 260 Input Arguments: 261 + dm - the shell DM 262 - func - the creation routine 263 264 Level: advanced 265 266 .seealso: DMShellSetLocalVector(), DMShellSetCreateMatrix() 267 @*/ 268 PetscErrorCode DMShellSetCreateLocalVector(DM dm,PetscErrorCode (*func)(DM,Vec*)) 269 { 270 271 PetscFunctionBegin; 272 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 273 dm->ops->createlocalvector = func; 274 PetscFunctionReturn(0); 275 } 276 277 #undef __FUNCT__ 278 #define __FUNCT__ "DMDestroy_Shell" 279 static PetscErrorCode DMDestroy_Shell(DM dm) 280 { 281 PetscErrorCode ierr; 282 DM_Shell *shell = (DM_Shell*)dm->data; 283 284 PetscFunctionBegin; 285 ierr = MatDestroy(&shell->A);CHKERRQ(ierr); 286 ierr = VecDestroy(&shell->Xglobal);CHKERRQ(ierr); 287 ierr = VecDestroy(&shell->Xlocal);CHKERRQ(ierr); 288 /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */ 289 ierr = PetscFree(shell);CHKERRQ(ierr); 290 PetscFunctionReturn(0); 291 } 292 293 #undef __FUNCT__ 294 #define __FUNCT__ "DMView_Shell" 295 static PetscErrorCode DMView_Shell(DM dm,PetscViewer v) 296 { 297 PetscErrorCode ierr; 298 DM_Shell *shell = (DM_Shell*)dm->data; 299 300 PetscFunctionBegin; 301 ierr = VecView(shell->Xglobal,v);CHKERRQ(ierr); 302 PetscFunctionReturn(0); 303 } 304 305 #undef __FUNCT__ 306 #define __FUNCT__ "DMLoad_Shell" 307 static PetscErrorCode DMLoad_Shell(DM dm,PetscViewer v) 308 { 309 PetscErrorCode ierr; 310 DM_Shell *shell = (DM_Shell*)dm->data; 311 312 PetscFunctionBegin; 313 ierr = VecCreate(PetscObjectComm((PetscObject)dm),&shell->Xglobal);CHKERRQ(ierr); 314 ierr = VecLoad(shell->Xglobal,v);CHKERRQ(ierr); 315 PetscFunctionReturn(0); 316 } 317 318 #undef __FUNCT__ 319 #define __FUNCT__ "DMCreate_Shell" 320 PETSC_EXTERN_C PetscErrorCode DMCreate_Shell(DM dm) 321 { 322 PetscErrorCode ierr; 323 DM_Shell *shell; 324 325 PetscFunctionBegin; 326 ierr = PetscNewLog(dm,DM_Shell,&shell);CHKERRQ(ierr); 327 dm->data = shell; 328 329 ierr = PetscObjectChangeTypeName((PetscObject)dm,DMSHELL);CHKERRQ(ierr); 330 331 dm->ops->destroy = DMDestroy_Shell; 332 dm->ops->createglobalvector = DMCreateGlobalVector_Shell; 333 dm->ops->createlocalvector = DMCreateLocalVector_Shell; 334 dm->ops->creatematrix = DMCreateMatrix_Shell; 335 dm->ops->view = DMView_Shell; 336 dm->ops->load = DMLoad_Shell; 337 PetscFunctionReturn(0); 338 } 339 340 #undef __FUNCT__ 341 #define __FUNCT__ "DMShellCreate" 342 /*@ 343 DMShellCreate - Creates a shell DM object, used to manage user-defined problem data 344 345 Collective on MPI_Comm 346 347 Input Parameter: 348 . comm - the processors that will share the global vector 349 350 Output Parameters: 351 . shell - the shell DM 352 353 Level: advanced 354 355 .seealso DMDestroy(), DMCreateGlobalVector(), DMCreateLocalVector() 356 @*/ 357 PetscErrorCode DMShellCreate(MPI_Comm comm,DM *dm) 358 { 359 PetscErrorCode ierr; 360 361 PetscFunctionBegin; 362 PetscValidPointer(dm,2); 363 ierr = DMCreate(comm,dm);CHKERRQ(ierr); 364 ierr = DMSetType(*dm,DMSHELL);CHKERRQ(ierr); 365 PetscFunctionReturn(0); 366 } 367