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