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 = PetscObjectCompose((PetscObject)*gvec,"DM",(PetscObject)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 204 #undef __FUNCT__ 205 #define __FUNCT__ "DMCreate_Shell" 206 PETSC_EXTERN_C PetscErrorCode DMCreate_Shell(DM dm) 207 { 208 PetscErrorCode ierr; 209 DM_Shell *shell; 210 211 PetscFunctionBegin; 212 ierr = PetscNewLog(dm,DM_Shell,&shell);CHKERRQ(ierr); 213 dm->data = shell; 214 215 ierr = PetscObjectChangeTypeName((PetscObject)dm,DMSHELL);CHKERRQ(ierr); 216 dm->ops->destroy = DMDestroy_Shell; 217 dm->ops->createglobalvector = DMCreateGlobalVector_Shell; 218 dm->ops->creatematrix = DMCreateMatrix_Shell; 219 PetscFunctionReturn(0); 220 } 221 222 #undef __FUNCT__ 223 #define __FUNCT__ "DMShellCreate" 224 /*@ 225 DMShellCreate - Creates a shell DM object, used to manage user-defined problem data 226 227 Collective on MPI_Comm 228 229 Input Parameter: 230 . comm - the processors that will share the global vector 231 232 Output Parameters: 233 . shell - the shell DM 234 235 Level: advanced 236 237 .seealso DMDestroy(), DMCreateGlobalVector() 238 @*/ 239 PetscErrorCode DMShellCreate(MPI_Comm comm,DM *dm) 240 { 241 PetscErrorCode ierr; 242 243 PetscFunctionBegin; 244 PetscValidPointer(dm,2); 245 ierr = DMCreate(comm,dm);CHKERRQ(ierr); 246 ierr = DMSetType(*dm,DMSHELL);CHKERRQ(ierr); 247 PetscFunctionReturn(0); 248 } 249