1 2 #include <petsc/private/dmdaimpl.h> /*I "petscdmda.h" I*/ 3 4 static PetscErrorCode DMDAGetElements_1D(DM dm,PetscInt *nel,PetscInt *nen,const PetscInt *e[]) 5 { 6 PetscErrorCode ierr; 7 DM_DA *da = (DM_DA*)dm->data; 8 PetscInt i,xs,xe,Xs,Xe; 9 PetscInt cnt=0; 10 11 PetscFunctionBegin; 12 if (!da->e) { 13 PetscInt corners[2]; 14 15 if (!da->s) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Cannot get elements for DMDA with zero stencil width"); 16 ierr = DMDAGetCorners(dm,&xs,0,0,&xe,0,0);CHKERRQ(ierr); 17 ierr = DMDAGetGhostCorners(dm,&Xs,0,0,&Xe,0,0);CHKERRQ(ierr); 18 xe += xs; Xe += Xs; if (xs != Xs) xs -= 1; 19 da->ne = 1*(xe - xs - 1); 20 ierr = PetscMalloc1(1 + 2*da->ne,&da->e);CHKERRQ(ierr); 21 for (i=xs; i<xe-1; i++) { 22 da->e[cnt++] = (i-Xs); 23 da->e[cnt++] = (i-Xs+1); 24 } 25 da->nen = 2; 26 27 corners[0] = (xs -Xs); 28 corners[1] = (xe-1-Xs); 29 ierr = ISCreateGeneral(PETSC_COMM_SELF,2,corners,PETSC_COPY_VALUES,&da->ecorners);CHKERRQ(ierr); 30 } 31 *nel = da->ne; 32 *nen = da->nen; 33 *e = da->e; 34 PetscFunctionReturn(0); 35 } 36 37 static PetscErrorCode DMDAGetElements_2D(DM dm,PetscInt *nel,PetscInt *nen,const PetscInt *e[]) 38 { 39 PetscErrorCode ierr; 40 DM_DA *da = (DM_DA*)dm->data; 41 PetscInt i,xs,xe,Xs,Xe; 42 PetscInt j,ys,ye,Ys,Ye; 43 PetscInt cnt=0, cell[4], ns=2; 44 PetscInt c, split[] = {0,1,3, 45 2,3,1}; 46 47 PetscFunctionBegin; 48 if (!da->e) { 49 PetscInt corners[4],nn = 0; 50 51 if (!da->s) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Cannot get elements for DMDA with zero stencil width"); 52 53 switch (da->elementtype) { 54 case DMDA_ELEMENT_Q1: 55 da->nen = 4; 56 break; 57 case DMDA_ELEMENT_P1: 58 da->nen = 3; 59 break; 60 default: 61 SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Unknown element type %d",da->elementtype); 62 break; 63 } 64 nn = da->nen; 65 66 if (da->elementtype == DMDA_ELEMENT_P1) {ns=2;} 67 if (da->elementtype == DMDA_ELEMENT_Q1) {ns=1;} 68 ierr = DMDAGetCorners(dm,&xs,&ys,0,&xe,&ye,0);CHKERRQ(ierr); 69 ierr = DMDAGetGhostCorners(dm,&Xs,&Ys,0,&Xe,&Ye,0);CHKERRQ(ierr); 70 xe += xs; Xe += Xs; if (xs != Xs) xs -= 1; 71 ye += ys; Ye += Ys; if (ys != Ys) ys -= 1; 72 da->ne = ns*(xe - xs - 1)*(ye - ys - 1); 73 ierr = PetscMalloc1(1 + nn*da->ne,&da->e);CHKERRQ(ierr); 74 for (j=ys; j<ye-1; j++) { 75 for (i=xs; i<xe-1; i++) { 76 cell[0] = (i-Xs ) + (j-Ys )*(Xe-Xs); 77 cell[1] = (i-Xs+1) + (j-Ys )*(Xe-Xs); 78 cell[2] = (i-Xs+1) + (j-Ys+1)*(Xe-Xs); 79 cell[3] = (i-Xs ) + (j-Ys+1)*(Xe-Xs); 80 if (da->elementtype == DMDA_ELEMENT_P1) { 81 for (c=0; c<ns*nn; c++) da->e[cnt++] = cell[split[c]]; 82 } 83 if (da->elementtype == DMDA_ELEMENT_Q1) { 84 for (c=0; c<ns*nn; c++) da->e[cnt++] = cell[c]; 85 } 86 } 87 } 88 89 corners[0] = (xs -Xs) + (ys -Ys)*(Xe-Xs); 90 corners[1] = (xe-1-Xs) + (ys -Ys)*(Xe-Xs); 91 corners[2] = (xs -Xs) + (ye-1-Ys)*(Xe-Xs); 92 corners[3] = (xe-1-Xs) + (ye-1-Ys)*(Xe-Xs); 93 ierr = ISCreateGeneral(PETSC_COMM_SELF,4,corners,PETSC_COPY_VALUES,&da->ecorners);CHKERRQ(ierr); 94 } 95 *nel = da->ne; 96 *nen = da->nen; 97 *e = da->e; 98 PetscFunctionReturn(0); 99 } 100 101 static PetscErrorCode DMDAGetElements_3D(DM dm,PetscInt *nel,PetscInt *nen,const PetscInt *e[]) 102 { 103 PetscErrorCode ierr; 104 DM_DA *da = (DM_DA*)dm->data; 105 PetscInt i,xs,xe,Xs,Xe; 106 PetscInt j,ys,ye,Ys,Ye; 107 PetscInt k,zs,ze,Zs,Ze; 108 PetscInt cnt=0, cell[8], ns=6; 109 PetscInt c, split[] = {0,1,3,7, 110 0,1,7,4, 111 1,2,3,7, 112 1,2,7,6, 113 1,4,5,7, 114 1,5,6,7}; 115 116 PetscFunctionBegin; 117 if (!da->e) { 118 PetscInt corners[8],nn = 0; 119 120 if (!da->s) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Cannot get elements for DMDA with zero stencil width"); 121 122 switch (da->elementtype) { 123 case DMDA_ELEMENT_Q1: 124 da->nen = 8; 125 break; 126 case DMDA_ELEMENT_P1: 127 da->nen = 4; 128 break; 129 default: 130 SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Unknown element type %d",da->elementtype); 131 break; 132 } 133 nn = da->nen; 134 135 if (da->elementtype == DMDA_ELEMENT_P1) {ns=6;} 136 if (da->elementtype == DMDA_ELEMENT_Q1) {ns=1;} 137 ierr = DMDAGetCorners(dm,&xs,&ys,&zs,&xe,&ye,&ze);CHKERRQ(ierr); 138 ierr = DMDAGetGhostCorners(dm,&Xs,&Ys,&Zs,&Xe,&Ye,&Ze);CHKERRQ(ierr); 139 xe += xs; Xe += Xs; if (xs != Xs) xs -= 1; 140 ye += ys; Ye += Ys; if (ys != Ys) ys -= 1; 141 ze += zs; Ze += Zs; if (zs != Zs) zs -= 1; 142 da->ne = ns*(xe - xs - 1)*(ye - ys - 1)*(ze - zs - 1); 143 ierr = PetscMalloc1(1 + nn*da->ne,&da->e);CHKERRQ(ierr); 144 for (k=zs; k<ze-1; k++) { 145 for (j=ys; j<ye-1; j++) { 146 for (i=xs; i<xe-1; i++) { 147 cell[0] = (i-Xs ) + (j-Ys )*(Xe-Xs) + (k-Zs )*(Xe-Xs)*(Ye-Ys); 148 cell[1] = (i-Xs+1) + (j-Ys )*(Xe-Xs) + (k-Zs )*(Xe-Xs)*(Ye-Ys); 149 cell[2] = (i-Xs+1) + (j-Ys+1)*(Xe-Xs) + (k-Zs )*(Xe-Xs)*(Ye-Ys); 150 cell[3] = (i-Xs ) + (j-Ys+1)*(Xe-Xs) + (k-Zs )*(Xe-Xs)*(Ye-Ys); 151 cell[4] = (i-Xs ) + (j-Ys )*(Xe-Xs) + (k-Zs+1)*(Xe-Xs)*(Ye-Ys); 152 cell[5] = (i-Xs+1) + (j-Ys )*(Xe-Xs) + (k-Zs+1)*(Xe-Xs)*(Ye-Ys); 153 cell[6] = (i-Xs+1) + (j-Ys+1)*(Xe-Xs) + (k-Zs+1)*(Xe-Xs)*(Ye-Ys); 154 cell[7] = (i-Xs ) + (j-Ys+1)*(Xe-Xs) + (k-Zs+1)*(Xe-Xs)*(Ye-Ys); 155 if (da->elementtype == DMDA_ELEMENT_P1) { 156 for (c=0; c<ns*nn; c++) da->e[cnt++] = cell[split[c]]; 157 } 158 if (da->elementtype == DMDA_ELEMENT_Q1) { 159 for (c=0; c<ns*nn; c++) da->e[cnt++] = cell[c]; 160 } 161 } 162 } 163 } 164 165 corners[0] = (xs -Xs) + (ys -Ys)*(Xe-Xs) + (zs-Zs )*(Xe-Xs)*(Ye-Ys); 166 corners[1] = (xe-1-Xs) + (ys -Ys)*(Xe-Xs) + (zs-Zs )*(Xe-Xs)*(Ye-Ys); 167 corners[2] = (xs -Xs) + (ye-1-Ys)*(Xe-Xs) + (zs-Zs )*(Xe-Xs)*(Ye-Ys); 168 corners[3] = (xe-1-Xs) + (ye-1-Ys)*(Xe-Xs) + (zs-Zs )*(Xe-Xs)*(Ye-Ys); 169 corners[4] = (xs -Xs) + (ys -Ys)*(Xe-Xs) + (ze-1-Zs)*(Xe-Xs)*(Ye-Ys); 170 corners[5] = (xe-1-Xs) + (ys -Ys)*(Xe-Xs) + (ze-1-Zs)*(Xe-Xs)*(Ye-Ys); 171 corners[6] = (xs -Xs) + (ye-1-Ys)*(Xe-Xs) + (ze-1-Zs)*(Xe-Xs)*(Ye-Ys); 172 corners[7] = (xe-1-Xs) + (ye-1-Ys)*(Xe-Xs) + (ze-1-Zs)*(Xe-Xs)*(Ye-Ys); 173 ierr = ISCreateGeneral(PETSC_COMM_SELF,8,corners,PETSC_COPY_VALUES,&da->ecorners);CHKERRQ(ierr); 174 } 175 *nel = da->ne; 176 *nen = da->nen; 177 *e = da->e; 178 PetscFunctionReturn(0); 179 } 180 181 /*@ 182 DMDAGetElementsCorners - Returns the global (x,y,z) indices of the lower left 183 corner of the non-overlapping decomposition identified by DMDAGetElements() 184 185 Not Collective 186 187 Input Parameter: 188 . da - the DM object 189 190 Output Parameters: 191 + gx - the x index 192 . gy - the y index 193 - gz - the z index 194 195 Level: intermediate 196 197 Notes: 198 199 .seealso: DMDAElementType, DMDASetElementType(), DMDAGetElements() 200 @*/ 201 PetscErrorCode DMDAGetElementsCorners(DM da, PetscInt *gx, PetscInt *gy, PetscInt *gz) 202 { 203 PetscInt xs,Xs; 204 PetscInt ys,Ys; 205 PetscInt zs,Zs; 206 PetscBool isda; 207 PetscErrorCode ierr; 208 209 PetscFunctionBegin; 210 PetscValidHeaderSpecificType(da,DM_CLASSID,1,DMDA); 211 if (gx) PetscValidIntPointer(gx,2); 212 if (gy) PetscValidIntPointer(gy,3); 213 if (gz) PetscValidIntPointer(gz,4); 214 ierr = PetscObjectTypeCompare((PetscObject)da,DMDA,&isda);CHKERRQ(ierr); 215 if (!isda) SETERRQ1(PetscObjectComm((PetscObject)da),PETSC_ERR_USER,"Not for DM type %s",((PetscObject)da)->type_name); 216 ierr = DMDAGetCorners(da,&xs,&ys,&zs,NULL,NULL,NULL);CHKERRQ(ierr); 217 ierr = DMDAGetGhostCorners(da,&Xs,&Ys,&Zs,NULL,NULL,NULL);CHKERRQ(ierr); 218 if (xs != Xs) xs -= 1; 219 if (ys != Ys) ys -= 1; 220 if (zs != Zs) zs -= 1; 221 if (gx) *gx = xs; 222 if (gy) *gy = ys; 223 if (gz) *gz = zs; 224 PetscFunctionReturn(0); 225 } 226 227 /*@ 228 DMDAGetElementsSizes - Gets the local number of elements per direction for the non-overlapping decomposition identified by DMDAGetElements() 229 230 Not Collective 231 232 Input Parameter: 233 . da - the DM object 234 235 Output Parameters: 236 + mx - number of local elements in x-direction 237 . my - number of local elements in y-direction 238 - mz - number of local elements in z-direction 239 240 Level: intermediate 241 242 Notes: 243 It returns the same number of elements, irrespective of the DMDAElementType 244 245 .seealso: DMDAElementType, DMDASetElementType(), DMDAGetElements 246 @*/ 247 PetscErrorCode DMDAGetElementsSizes(DM da, PetscInt *mx, PetscInt *my, PetscInt *mz) 248 { 249 PetscInt xs,xe,Xs; 250 PetscInt ys,ye,Ys; 251 PetscInt zs,ze,Zs; 252 PetscInt dim; 253 PetscBool isda; 254 PetscErrorCode ierr; 255 256 PetscFunctionBegin; 257 PetscValidHeaderSpecificType(da,DM_CLASSID,1,DMDA); 258 if (mx) PetscValidIntPointer(mx,2); 259 if (my) PetscValidIntPointer(my,3); 260 if (mz) PetscValidIntPointer(mz,4); 261 ierr = PetscObjectTypeCompare((PetscObject)da,DMDA,&isda);CHKERRQ(ierr); 262 if (!isda) SETERRQ1(PetscObjectComm((PetscObject)da),PETSC_ERR_USER,"Not for DM type %s",((PetscObject)da)->type_name); 263 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xe,&ye,&ze);CHKERRQ(ierr); 264 ierr = DMDAGetGhostCorners(da,&Xs,&Ys,&Zs,NULL,NULL,NULL);CHKERRQ(ierr); 265 xe += xs; if (xs != Xs) xs -= 1; 266 ye += ys; if (ys != Ys) ys -= 1; 267 ze += zs; if (zs != Zs) zs -= 1; 268 if (mx) *mx = 0; 269 if (my) *my = 0; 270 if (mz) *mz = 0; 271 ierr = DMGetDimension(da,&dim);CHKERRQ(ierr); 272 switch (dim) { 273 case 3: 274 if (mz) *mz = ze - zs - 1; 275 case 2: 276 if (my) *my = ye - ys - 1; 277 case 1: 278 if (mx) *mx = xe - xs - 1; 279 break; 280 } 281 PetscFunctionReturn(0); 282 } 283 284 /*@ 285 DMDASetElementType - Sets the element type to be returned by DMDAGetElements() 286 287 Not Collective 288 289 Input Parameter: 290 . da - the DMDA object 291 292 Output Parameters: 293 . etype - the element type, currently either DMDA_ELEMENT_P1 or DMDA_ELEMENT_Q1 294 295 Level: intermediate 296 297 .seealso: DMDAElementType, DMDAGetElementType(), DMDAGetElements(), DMDARestoreElements() 298 @*/ 299 PetscErrorCode DMDASetElementType(DM da, DMDAElementType etype) 300 { 301 DM_DA *dd = (DM_DA*)da->data; 302 PetscErrorCode ierr; 303 PetscBool isda; 304 305 PetscFunctionBegin; 306 PetscValidHeaderSpecificType(da,DM_CLASSID,1,DMDA); 307 PetscValidLogicalCollectiveEnum(da,etype,2); 308 ierr = PetscObjectTypeCompare((PetscObject)da,DMDA,&isda);CHKERRQ(ierr); 309 if (!isda) PetscFunctionReturn(0); 310 if (dd->elementtype != etype) { 311 ierr = PetscFree(dd->e);CHKERRQ(ierr); 312 ierr = ISDestroy(&dd->ecorners);CHKERRQ(ierr); 313 314 dd->elementtype = etype; 315 dd->ne = 0; 316 dd->nen = 0; 317 dd->e = NULL; 318 } 319 PetscFunctionReturn(0); 320 } 321 322 /*@ 323 DMDAGetElementType - Gets the element type to be returned by DMDAGetElements() 324 325 Not Collective 326 327 Input Parameter: 328 . da - the DMDA object 329 330 Output Parameters: 331 . etype - the element type, currently either DMDA_ELEMENT_P1 or DMDA_ELEMENT_Q1 332 333 Level: intermediate 334 335 .seealso: DMDAElementType, DMDASetElementType(), DMDAGetElements(), DMDARestoreElements() 336 @*/ 337 PetscErrorCode DMDAGetElementType(DM da, DMDAElementType *etype) 338 { 339 DM_DA *dd = (DM_DA*)da->data; 340 PetscErrorCode ierr; 341 PetscBool isda; 342 343 PetscFunctionBegin; 344 PetscValidHeaderSpecificType(da,DM_CLASSID,1,DMDA); 345 PetscValidPointer(etype,2); 346 ierr = PetscObjectTypeCompare((PetscObject)da,DMDA,&isda);CHKERRQ(ierr); 347 if (!isda) SETERRQ1(PetscObjectComm((PetscObject)da),PETSC_ERR_USER,"Not for DM type %s",((PetscObject)da)->type_name); 348 *etype = dd->elementtype; 349 PetscFunctionReturn(0); 350 } 351 352 /*@C 353 DMDAGetElements - Gets an array containing the indices (in local coordinates) 354 of all the local elements 355 356 Not Collective 357 358 Input Parameter: 359 . dm - the DM object 360 361 Output Parameters: 362 + nel - number of local elements 363 . nen - number of element nodes 364 - e - the local indices of the elements' vertices 365 366 Level: intermediate 367 368 Notes: 369 Call DMDARestoreElements() once you have finished accessing the elements. 370 371 Each process uniquely owns a subset of the elements. That is no element is owned by two or more processes. 372 373 If on each process you integrate over its owned elements and use ADD_VALUES in Vec/MatSetValuesLocal() then you'll obtain the correct result. 374 375 Not supported in Fortran 376 377 .seealso: DMDAElementType, DMDASetElementType(), VecSetValuesLocal(), MatSetValuesLocal(), DMGlobalToLocalBegin(), DMLocalToGlobalBegin() 378 @*/ 379 PetscErrorCode DMDAGetElements(DM dm,PetscInt *nel,PetscInt *nen,const PetscInt *e[]) 380 { 381 PetscInt dim; 382 PetscErrorCode ierr; 383 DM_DA *dd = (DM_DA*)dm->data; 384 PetscBool isda; 385 386 PetscFunctionBegin; 387 PetscValidHeaderSpecificType(dm,DM_CLASSID,1,DMDA); 388 PetscValidIntPointer(nel,2); 389 PetscValidIntPointer(nen,3); 390 PetscValidPointer(e,4); 391 ierr = PetscObjectTypeCompare((PetscObject)dm,DMDA,&isda);CHKERRQ(ierr); 392 if (!isda) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Not for DM type %s",((PetscObject)dm)->type_name); 393 if (dd->stencil_type == DMDA_STENCIL_STAR) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"DMDAGetElements() requires you use a stencil type of DMDA_STENCIL_BOX"); 394 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 395 if (dd->e) { 396 *nel = dd->ne; 397 *nen = dd->nen; 398 *e = dd->e; 399 PetscFunctionReturn(0); 400 } 401 if (dim==-1) { 402 *nel = 0; *nen = 0; *e = NULL; 403 } else if (dim==1) { 404 ierr = DMDAGetElements_1D(dm,nel,nen,e);CHKERRQ(ierr); 405 } else if (dim==2) { 406 ierr = DMDAGetElements_2D(dm,nel,nen,e);CHKERRQ(ierr); 407 } else if (dim==3) { 408 ierr = DMDAGetElements_3D(dm,nel,nen,e);CHKERRQ(ierr); 409 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 410 PetscFunctionReturn(0); 411 } 412 413 /*@ 414 DMDAGetSubdomainCornersIS - Gets an index set containing the corner indices (in local coordinates) 415 of the non-overlapping decomposition identified by DMDAGetElements 416 417 Not Collective 418 419 Input Parameter: 420 . dm - the DM object 421 422 Output Parameters: 423 . is - the index set 424 425 Level: intermediate 426 427 Notes: 428 Call DMDARestoreSubdomainCornersIS() once you have finished accessing the index set. 429 430 .seealso: DMDAElementType, DMDASetElementType(), DMDAGetElements(), DMDARestoreElementsCornersIS() 431 @*/ 432 PetscErrorCode DMDAGetSubdomainCornersIS(DM dm,IS *is) 433 { 434 PetscErrorCode ierr; 435 DM_DA *dd = (DM_DA*)dm->data; 436 PetscBool isda; 437 438 PetscFunctionBegin; 439 PetscValidHeaderSpecificType(dm,DM_CLASSID,1,DMDA); 440 PetscValidPointer(is,2); 441 ierr = PetscObjectTypeCompare((PetscObject)dm,DMDA,&isda);CHKERRQ(ierr); 442 if (!isda) SETERRQ1(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Not for DM type %s",((PetscObject)dm)->type_name); 443 if (dd->stencil_type == DMDA_STENCIL_STAR) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"DMDAGetElement() requires you use a stencil type of DMDA_STENCIL_BOX"); 444 if (!dd->ecorners) { /* compute elements if not yet done */ 445 const PetscInt *e; 446 PetscInt nel,nen; 447 448 ierr = DMDAGetElements(dm,&nel,&nen,&e);CHKERRQ(ierr); 449 ierr = DMDARestoreElements(dm,&nel,&nen,&e);CHKERRQ(ierr); 450 } 451 *is = dd->ecorners; 452 PetscFunctionReturn(0); 453 } 454 455 /*@C 456 DMDARestoreElements - Restores the array obtained with DMDAGetElements() 457 458 Not Collective 459 460 Input Parameter: 461 + dm - the DM object 462 . nel - number of local elements 463 . nen - number of element nodes 464 - e - the local indices of the elements' vertices 465 466 Level: intermediate 467 468 Note: You should not access these values after you have called this routine. 469 470 This restore signals the DMDA object that you no longer need access to the array information. 471 472 Not supported in Fortran 473 474 .seealso: DMDAElementType, DMDASetElementType(), DMDAGetElements() 475 @*/ 476 PetscErrorCode DMDARestoreElements(DM dm,PetscInt *nel,PetscInt *nen,const PetscInt *e[]) 477 { 478 PetscFunctionBegin; 479 PetscValidHeaderSpecificType(dm,DM_CLASSID,1,DMDA); 480 PetscValidIntPointer(nel,2); 481 PetscValidIntPointer(nen,3); 482 PetscValidPointer(e,4); 483 *nel = 0; 484 *nen = -1; 485 *e = NULL; 486 PetscFunctionReturn(0); 487 } 488 489 /*@ 490 DMDARestoreSubdomainCornersIS - Restores the IS obtained with DMDAGetSubdomainCornersIS() 491 492 Not Collective 493 494 Input Parameter: 495 + dm - the DM object 496 - is - the index set 497 498 Level: intermediate 499 500 Note: 501 502 .seealso: DMDAElementType, DMDASetElementType(), DMDAGetSubdomainCornersIS() 503 @*/ 504 PetscErrorCode DMDARestoreSubdomainCornersIS(DM dm,IS *is) 505 { 506 PetscFunctionBegin; 507 PetscValidHeaderSpecificType(dm,DM_CLASSID,1,DMDA); 508 PetscValidHeaderSpecific(*is,IS_CLASSID,2); 509 *is = NULL; 510 PetscFunctionReturn(0); 511 } 512