1 #include <petsc-private/daimpl.h> 2 3 #undef __FUNCT__ 4 #define __FUNCT__ "DMDACreatePatchIS" 5 6 PetscErrorCode DMDACreatePatchIS(DM da,MatStencil *lower,MatStencil *upper,IS *is) 7 { 8 PetscErrorCode ierr; 9 PetscInt i,j,k,idx; 10 PetscInt ii,jj,kk; 11 Vec v; 12 PetscInt n,pn,bs; 13 PetscMPIInt rank; 14 PetscSF sf,psf; 15 PetscLayout map; 16 MPI_Comm comm; 17 PetscInt *natidx,*globidx,*leafidx; 18 PetscInt *pnatidx,*pleafidx; 19 PetscInt base; 20 PetscInt ox,oy,oz; 21 DM_DA *dd; 22 PetscFunctionBegin; 23 24 comm = ((PetscObject)da)->comm; 25 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 26 dd = (DM_DA *)da->data; 27 28 /* construct the natural mapping */ 29 ierr = DMGetGlobalVector(da,&v);CHKERRQ(ierr); 30 ierr = VecGetLocalSize(v,&n);CHKERRQ(ierr); 31 ierr = VecGetOwnershipRange(v,&base,PETSC_NULL);CHKERRQ(ierr); 32 ierr = VecGetBlockSize(v,&bs);CHKERRQ(ierr); 33 ierr = DMRestoreGlobalVector(da,&v);CHKERRQ(ierr); 34 35 /* construct the layout */ 36 ierr = PetscLayoutCreate(comm,&map);CHKERRQ(ierr); 37 ierr = PetscLayoutSetBlockSize(map,1);CHKERRQ(ierr); 38 ierr = PetscLayoutSetLocalSize(map,n);CHKERRQ(ierr); 39 ierr = PetscLayoutSetUp(map);CHKERRQ(ierr); 40 41 /* construct the list of natural indices on this process when PETSc ordering is considered */ 42 ierr = DMDAGetOffset(da,&ox,&oy,&oz);CHKERRQ(ierr); 43 ierr = PetscMalloc(sizeof(PetscInt)*n,&natidx);CHKERRQ(ierr); 44 ierr = PetscMalloc(sizeof(PetscInt)*n,&globidx);CHKERRQ(ierr); 45 ierr = PetscMalloc(sizeof(PetscInt)*n,&leafidx);CHKERRQ(ierr); 46 idx = 0; 47 for (k=dd->zs;k<dd->ze;k++) { 48 for (j=dd->ys;j<dd->ye;j++) { 49 for (i=dd->xs;i<dd->xe;i++) { 50 natidx[idx] = i + dd->w*(j*dd->M + k*dd->M*dd->N); 51 globidx[idx] = base + idx; 52 leafidx[idx] = 0; 53 idx++; 54 } 55 } 56 } 57 58 if (idx != n) SETERRQ(comm,PETSC_ERR_ARG_WRONGSTATE, "for some reason the count is wrong."); 59 60 /* construct the SF going from the natural indices to the local set of PETSc indices */ 61 ierr = PetscSFCreate(comm,&sf);CHKERRQ(ierr); 62 ierr = PetscSFSetFromOptions(sf);CHKERRQ(ierr); 63 ierr = PetscSFSetGraphLayout(sf,map,n,PETSC_NULL,PETSC_OWN_POINTER,natidx);CHKERRQ(ierr); 64 65 /* broadcast the global indices over to the corresponding natural indices */ 66 ierr = PetscSFGatherBegin(sf,MPIU_INT,globidx,leafidx);CHKERRQ(ierr); 67 ierr = PetscSFGatherEnd(sf,MPIU_INT,globidx,leafidx);CHKERRQ(ierr); 68 ierr = PetscSFDestroy(&sf);CHKERRQ(ierr); 69 70 71 pn = dd->w*(upper->k - lower->k)*(upper->j - lower->j)*(upper->i - lower->i); 72 ierr = PetscMalloc(sizeof(PetscInt)*pn,&pnatidx);CHKERRQ(ierr); 73 ierr = PetscMalloc(sizeof(PetscInt)*pn,&pleafidx);CHKERRQ(ierr); 74 idx = 0; 75 for (k=lower->k-oz;k<upper->k-oz;k++) { 76 for (j=lower->j-oy;j<upper->j-oy;j++) { 77 for (i=dd->w*(lower->i-ox);i<dd->w*(upper->i-ox);i++) { 78 ii = i % (dd->w*dd->M); 79 jj = j % dd->N; 80 kk = k % dd->P; 81 if (ii < 0) ii = dd->w*dd->M + ii; 82 if (jj < 0) jj = dd->N + jj; 83 if (kk < 0) kk = dd->P + kk; 84 pnatidx[idx] = ii + dd->w*(jj*dd->M + kk*dd->M*dd->N); 85 idx++; 86 } 87 } 88 } 89 90 if (idx != pn) SETERRQ(comm,PETSC_ERR_ARG_WRONGSTATE, "for some reason the count is wrong"); 91 92 ierr = PetscSFCreate(comm,&psf);CHKERRQ(ierr); 93 ierr = PetscSFSetFromOptions(psf);CHKERRQ(ierr); 94 ierr = PetscSFSetGraphLayout(psf,map,pn,PETSC_NULL,PETSC_OWN_POINTER,pnatidx);CHKERRQ(ierr); 95 96 /* broadcast the global indices through to the patch */ 97 ierr = PetscSFBcastBegin(psf,MPIU_INT,leafidx,pleafidx);CHKERRQ(ierr); 98 ierr = PetscSFBcastEnd(psf,MPIU_INT,leafidx,pleafidx);CHKERRQ(ierr); 99 100 /* create the IS */ 101 ierr = ISCreateGeneral(comm,pn,pleafidx,PETSC_OWN_POINTER,is);CHKERRQ(ierr); 102 103 ierr = PetscSFDestroy(&psf);CHKERRQ(ierr); 104 105 ierr = PetscLayoutDestroy(&map);CHKERRQ(ierr); 106 107 ierr = PetscFree(globidx);CHKERRQ(ierr); 108 ierr = PetscFree(leafidx);CHKERRQ(ierr); 109 ierr = PetscFree(natidx);CHKERRQ(ierr); 110 ierr = PetscFree(pnatidx);CHKERRQ(ierr); 111 112 PetscFunctionReturn(0); 113 } 114 115 #undef __FUNCT__ 116 #define __FUNCT__ "DMDASubDomainDA_Private" 117 PetscErrorCode DMDASubDomainDA_Private(DM dm, DM *dddm) { 118 DM da; 119 DM_DA *dd; 120 PetscErrorCode ierr; 121 DMDALocalInfo info; 122 PetscReal lmin[3],lmax[3]; 123 PetscInt xsize,ysize,zsize; 124 PetscInt xo,yo,zo; 125 126 PetscFunctionBegin; 127 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 128 ierr = DMDACreate(PETSC_COMM_SELF,&da);CHKERRQ(ierr); 129 ierr = DMSetOptionsPrefix(da,"sub_");CHKERRQ(ierr); 130 131 ierr = DMDASetDim(da, info.dim);CHKERRQ(ierr); 132 ierr = DMDASetDof(da, info.dof);CHKERRQ(ierr); 133 134 ierr = DMDASetStencilType(da,info.st);CHKERRQ(ierr); 135 ierr = DMDASetStencilWidth(da,info.sw);CHKERRQ(ierr); 136 137 dd = (DM_DA *)dm->data; 138 139 /* local with overlap */ 140 xsize = info.xm; 141 ysize = info.ym; 142 zsize = info.zm; 143 xo = info.xs; 144 yo = info.ys; 145 zo = info.zs; 146 if (info.bx == DMDA_BOUNDARY_PERIODIC || (info.xs != 0)) { 147 xsize += dd->overlap; 148 xo -= dd->overlap; 149 } 150 if (info.by == DMDA_BOUNDARY_PERIODIC || (info.ys != 0)) { 151 ysize += dd->overlap; 152 yo -= dd->overlap; 153 } 154 if (info.bz == DMDA_BOUNDARY_PERIODIC || (info.zs != 0)) { 155 zsize += dd->overlap; 156 zo -= dd->overlap; 157 } 158 159 if (info.bx == DMDA_BOUNDARY_PERIODIC || (info.xs+info.xm != info.mx)) { 160 xsize += dd->overlap; 161 } 162 if (info.by == DMDA_BOUNDARY_PERIODIC || (info.ys+info.ym != info.my)) { 163 ysize += dd->overlap; 164 } 165 if (info.bz == DMDA_BOUNDARY_PERIODIC || (info.zs+info.zm != info.mz)) { 166 zsize += dd->overlap; 167 } 168 169 ierr = DMDASetSizes(da, xsize, ysize, zsize);CHKERRQ(ierr); 170 ierr = DMDASetNumProcs(da, 1, 1, 1);CHKERRQ(ierr); 171 ierr = DMDASetBoundaryType(da, DMDA_BOUNDARY_GHOSTED, DMDA_BOUNDARY_GHOSTED, DMDA_BOUNDARY_GHOSTED);CHKERRQ(ierr); 172 173 /* set up as a block instead */ 174 ierr = DMSetUp(da);CHKERRQ(ierr); 175 ierr = DMDASetOffset(da,xo,yo,zo);CHKERRQ(ierr); 176 177 178 /* todo - nonuniform coordinates */ 179 ierr = DMDAGetLocalBoundingBox(dm,lmin,lmax);CHKERRQ(ierr); 180 ierr = DMDASetUniformCoordinates(da,lmin[0],lmax[0],lmin[1],lmax[1],lmin[2],lmax[2]);CHKERRQ(ierr); 181 182 dd = (DM_DA *)da->data; 183 dd->Mo = info.mx; 184 dd->No = info.my; 185 dd->Po = info.mz; 186 187 *dddm = da; 188 189 PetscFunctionReturn(0); 190 } 191 192 #undef __FUNCT__ 193 #define __FUNCT__ "DMCreateDomainDecompositionScatters_DA" 194 /* 195 Fills the local vector problem on the subdomain from the global problem. 196 197 */ 198 PetscErrorCode DMCreateDomainDecompositionScatters_DA(DM dm,PetscInt nsubdms,DM *subdms,VecScatter **iscat,VecScatter **oscat, VecScatter **lscat) { 199 PetscErrorCode ierr; 200 DMDALocalInfo dinfo,sinfo; 201 IS isis,idis,osis,odis,gsis,gdis; 202 PetscInt *ididx,*isidx,*odidx,*osidx,*gdidx,*gsidx,*idx_global,n_global,*idx_sub,n_sub; 203 PetscInt l,i,j,k,d,n_i,n_o,n_g,sl,dl,di,dj,dk,si,sj,sk; 204 Vec dvec,svec,slvec; 205 DM subdm; 206 207 PetscFunctionBegin; 208 209 /* allocate the arrays of scatters */ 210 if (iscat) {ierr = PetscMalloc(sizeof(VecScatter *),iscat);CHKERRQ(ierr);} 211 if (oscat) {ierr = PetscMalloc(sizeof(VecScatter *),oscat);CHKERRQ(ierr);} 212 if (lscat) {ierr = PetscMalloc(sizeof(VecScatter *),lscat);CHKERRQ(ierr);} 213 214 ierr = DMDAGetLocalInfo(dm,&dinfo);CHKERRQ(ierr); 215 ierr = DMDAGetGlobalIndices(dm,&n_global,&idx_global);CHKERRQ(ierr); 216 for (l = 0;l < nsubdms;l++) { 217 n_i = 0; 218 n_o = 0; 219 n_g = 0; 220 subdm = subdms[l]; 221 ierr = DMDAGetLocalInfo(subdm,&sinfo);CHKERRQ(ierr); 222 ierr = DMDAGetGlobalIndices(subdm,&n_sub,&idx_sub);CHKERRQ(ierr); 223 /* count the three region sizes */ 224 for (k=sinfo.gzs;k<sinfo.gzs+sinfo.gzm;k++) { 225 for (j=sinfo.gys;j<sinfo.gys+sinfo.gym;j++) { 226 for (i=sinfo.gxs;i<sinfo.gxs+sinfo.gxm;i++) { 227 for (d=0;d<sinfo.dof;d++) { 228 if (k >= sinfo.zs && j >= sinfo.ys && i >= sinfo.xs && 229 k < sinfo.zs+sinfo.zm && j < sinfo.ys+sinfo.ym && i < sinfo.xs+sinfo.xm) { 230 231 /* interior - subinterior overlap */ 232 if (k >= dinfo.zs && j >= dinfo.ys && i >= dinfo.xs && 233 k < dinfo.zs+dinfo.zm && j < dinfo.ys+dinfo.ym && i < dinfo.xs+dinfo.xm) { 234 n_i++; 235 } 236 /* ghost - subinterior overlap */ 237 if (k >= dinfo.gzs && j >= dinfo.gys && i >= dinfo.gxs && 238 k < dinfo.gzs+dinfo.gzm && j < dinfo.gys+dinfo.gym && i < dinfo.gxs+dinfo.gxm) { 239 n_o++; 240 } 241 } 242 243 /* ghost - subghost overlap */ 244 if (k >= dinfo.gzs && j >= dinfo.gys && i >= dinfo.gxs && 245 k < dinfo.gzs+dinfo.gzm && j < dinfo.gys+dinfo.gym && i < dinfo.gxs+dinfo.gxm) { 246 n_g++; 247 } 248 } 249 } 250 } 251 } 252 253 if (n_g == 0) SETERRQ(((PetscObject)subdm)->comm,PETSC_ERR_ARG_WRONGSTATE,"Processor-local domain and subdomain do not intersect!"); 254 255 /* local and subdomain local index set indices */ 256 ierr = PetscMalloc(sizeof(PetscInt)*n_i,&ididx);CHKERRQ(ierr); 257 ierr = PetscMalloc(sizeof(PetscInt)*n_i,&isidx);CHKERRQ(ierr); 258 259 ierr = PetscMalloc(sizeof(PetscInt)*n_o,&odidx);CHKERRQ(ierr); 260 ierr = PetscMalloc(sizeof(PetscInt)*n_o,&osidx);CHKERRQ(ierr); 261 262 ierr = PetscMalloc(sizeof(PetscInt)*n_g,&gdidx);CHKERRQ(ierr); 263 ierr = PetscMalloc(sizeof(PetscInt)*n_g,&gsidx);CHKERRQ(ierr); 264 265 n_i = 0; n_o = 0;n_g = 0; 266 for (k=sinfo.gzs;k<sinfo.gzs+sinfo.gzm;k++) { 267 for (j=sinfo.gys;j<sinfo.gys+sinfo.gym;j++) { 268 for (i=sinfo.gxs;i<sinfo.gxs+sinfo.gxm;i++) { 269 for (d=0;d<sinfo.dof;d++) { 270 si = i - sinfo.gxs; 271 sj = j - sinfo.gys; 272 sk = k - sinfo.gzs; 273 sl = d + sinfo.dof*(si + sinfo.gxm*(sj + sinfo.gym*sk)); 274 di = i - dinfo.gxs; 275 dj = j - dinfo.gys; 276 dk = k - dinfo.gzs; 277 dl = d + dinfo.dof*(di + dinfo.gxm*(dj + dinfo.gym*dk)); 278 279 if (k >= sinfo.zs && j >= sinfo.ys && i >= sinfo.xs && 280 k < sinfo.zs+sinfo.zm && j < sinfo.ys+sinfo.ym && i < sinfo.xs+sinfo.xm) { 281 282 /* interior - subinterior overlap */ 283 if (k >= dinfo.zs && j >= dinfo.ys && i >= dinfo.xs && 284 k < dinfo.zs+dinfo.zm && j < dinfo.ys+dinfo.ym && i < dinfo.xs+dinfo.xm) { 285 ididx[n_i] = idx_global[dl]; 286 isidx[n_i] = idx_sub[sl]; 287 n_i++; 288 } 289 /* ghost - subinterior overlap */ 290 if (k >= dinfo.gzs && j >= dinfo.gys && i >= dinfo.gxs && 291 k < dinfo.gzs+dinfo.gzm && j < dinfo.gys+dinfo.gym && i < dinfo.gxs+dinfo.gxm) { 292 odidx[n_o] = idx_global[dl]; 293 osidx[n_o] = idx_sub[sl]; 294 n_o++; 295 } 296 } 297 298 /* ghost - subghost overlap */ 299 if (k >= dinfo.gzs && j >= dinfo.gys && i >= dinfo.gxs && 300 k < dinfo.gzs+dinfo.gzm && j < dinfo.gys+dinfo.gym && i < dinfo.gxs+dinfo.gxm) { 301 gdidx[n_g] = idx_global[dl]; 302 gsidx[n_g] = sl; 303 n_g++; 304 } 305 } 306 } 307 } 308 } 309 310 ierr = ISCreateGeneral(PETSC_COMM_SELF,n_i,ididx,PETSC_OWN_POINTER,&idis);CHKERRQ(ierr); 311 ierr = ISCreateGeneral(PETSC_COMM_SELF,n_i,isidx,PETSC_OWN_POINTER,&isis);CHKERRQ(ierr); 312 313 ierr = ISCreateGeneral(PETSC_COMM_SELF,n_o,odidx,PETSC_OWN_POINTER,&odis);CHKERRQ(ierr); 314 ierr = ISCreateGeneral(PETSC_COMM_SELF,n_o,osidx,PETSC_OWN_POINTER,&osis);CHKERRQ(ierr); 315 316 ierr = ISCreateGeneral(PETSC_COMM_SELF,n_g,gdidx,PETSC_OWN_POINTER,&gdis);CHKERRQ(ierr); 317 ierr = ISCreateGeneral(PETSC_COMM_SELF,n_g,gsidx,PETSC_OWN_POINTER,&gsis);CHKERRQ(ierr); 318 319 /* form the scatter */ 320 ierr = DMGetGlobalVector(dm,&dvec);CHKERRQ(ierr); 321 ierr = DMGetGlobalVector(subdm,&svec);CHKERRQ(ierr); 322 ierr = DMGetLocalVector(subdm,&slvec);CHKERRQ(ierr); 323 324 if (iscat) {ierr = VecScatterCreate(dvec,idis,svec,isis,&(*iscat)[l]);CHKERRQ(ierr);} 325 if (oscat) {ierr = VecScatterCreate(dvec,odis,svec,osis,&(*oscat)[l]);CHKERRQ(ierr);} 326 if (lscat) {ierr = VecScatterCreate(dvec,gdis,slvec,gsis,&(*lscat)[l]);CHKERRQ(ierr);} 327 328 ierr = DMRestoreGlobalVector(dm,&dvec);CHKERRQ(ierr); 329 ierr = DMRestoreGlobalVector(subdm,&svec);CHKERRQ(ierr); 330 ierr = DMRestoreLocalVector(subdm,&slvec);CHKERRQ(ierr); 331 332 ierr = ISDestroy(&idis);CHKERRQ(ierr); 333 ierr = ISDestroy(&isis);CHKERRQ(ierr); 334 335 ierr = ISDestroy(&odis);CHKERRQ(ierr); 336 ierr = ISDestroy(&osis);CHKERRQ(ierr); 337 338 ierr = ISDestroy(&gdis);CHKERRQ(ierr); 339 ierr = ISDestroy(&gsis);CHKERRQ(ierr); 340 } 341 PetscFunctionReturn(0); 342 343 } 344 345 #undef __FUNCT__ 346 #define __FUNCT__ "DMDASubDomainIS_Private" 347 /* We have that the interior regions are going to be the same, but the ghost regions might not match up 348 349 ---------- 350 ---------- 351 --++++++o= 352 --++++++o= 353 --++++++o= 354 --++++++o= 355 --ooooooo= 356 --======== 357 358 Therefore, for each point in the overall, we must check if it's: 359 360 1. +: In the interior of the global dm; it lines up 361 2. o: In the overlap region -- for now the same as 1; no overlap 362 3. =: In the shared ghost region -- handled by DMCreateDomainDecompositionLocalScatter() 363 4. -: In the nonshared ghost region 364 */ 365 366 PetscErrorCode DMDASubDomainIS_Private(DM dm,DM subdm,IS *iis,IS *ois) { 367 PetscErrorCode ierr; 368 DMDALocalInfo info,subinfo; 369 PetscInt *iiidx,*oiidx,*gidx,gindx; 370 PetscInt i,j,k,d,n,nsub,nover,llindx,lindx,li,lj,lk,gi,gj,gk; 371 372 PetscFunctionBegin; 373 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 374 ierr = DMDAGetLocalInfo(subdm,&subinfo);CHKERRQ(ierr); 375 ierr = DMDAGetGlobalIndices(dm,&n,&gidx);CHKERRQ(ierr); 376 /* todo -- overlap */ 377 nsub = info.xm*info.ym*info.zm*info.dof; 378 nover = subinfo.xm*subinfo.ym*subinfo.zm*subinfo.dof; 379 /* iis is going to have size of the local problem's global part but have a lot of fill-in */ 380 ierr = PetscMalloc(sizeof(PetscInt)*nsub,&iiidx);CHKERRQ(ierr); 381 /* ois is going to have size of the local problem's global part */ 382 ierr = PetscMalloc(sizeof(PetscInt)*nover,&oiidx);CHKERRQ(ierr); 383 /* loop over the ghost region of the subdm and fill in the indices */ 384 for (k=subinfo.gzs;k<subinfo.gzs+subinfo.gzm;k++) { 385 for (j=subinfo.gys;j<subinfo.gys+subinfo.gym;j++) { 386 for (i=subinfo.gxs;i<subinfo.gxs+subinfo.gxm;i++) { 387 for (d=0;d<subinfo.dof;d++) { 388 li = i - subinfo.xs; 389 lj = j - subinfo.ys; 390 lk = k - subinfo.zs; 391 lindx = d + subinfo.dof*(li + subinfo.xm*(lj + subinfo.ym*lk)); 392 li = i - info.xs; 393 lj = j - info.ys; 394 lk = k - info.zs; 395 llindx = d + info.dof*(li + info.xm*(lj + info.ym*lk)); 396 gi = i - info.gxs; 397 gj = j - info.gys; 398 gk = k - info.gzs; 399 gindx = d + info.dof*(gi + info.gxm*(gj + info.gym*gk)); 400 401 /* check if the current point is inside the interior region */ 402 if (k >= info.zs && j >= info.ys && i >= info.xs && 403 k < info.zs+info.zm && j < info.ys+info.ym && i < info.xs+info.xm) { 404 iiidx[llindx] = gidx[gindx]; 405 oiidx[lindx] = gidx[gindx]; 406 /* overlap region */ 407 } else if (k >= subinfo.zs && j >= subinfo.ys && i >= subinfo.xs && 408 k < subinfo.zs+subinfo.zm && j < subinfo.ys+subinfo.ym && i < subinfo.xs+subinfo.xm) { 409 oiidx[lindx] = gidx[gindx]; 410 } 411 } 412 } 413 } 414 } 415 416 /* create the index sets */ 417 ierr = ISCreateGeneral(PETSC_COMM_SELF,nsub,iiidx,PETSC_OWN_POINTER,iis);CHKERRQ(ierr); 418 ierr = ISCreateGeneral(PETSC_COMM_SELF,nover,oiidx,PETSC_OWN_POINTER,ois);CHKERRQ(ierr); 419 PetscFunctionReturn(0); 420 } 421 422 #undef __FUNCT__ 423 #define __FUNCT__ "DMCreateDomainDecomposition_DA" 424 PetscErrorCode DMCreateDomainDecomposition_DA(DM dm,PetscInt *len,char ***names,IS **iis,IS **ois,DM **subdm) { 425 PetscErrorCode ierr; 426 IS iis0,ois0; 427 DM subdm0; 428 PetscFunctionBegin; 429 if (len)*len = 1; 430 431 if (iis) {ierr = PetscMalloc(sizeof(IS *),iis);CHKERRQ(ierr);} 432 if (ois) {ierr = PetscMalloc(sizeof(IS *),ois);CHKERRQ(ierr);} 433 if (subdm) {ierr = PetscMalloc(sizeof(DM *),subdm);CHKERRQ(ierr);} 434 if (names) {ierr = PetscMalloc(sizeof(char *),names);CHKERRQ(ierr);} 435 ierr = DMDASubDomainDA_Private(dm,&subdm0);CHKERRQ(ierr); 436 ierr = DMDASubDomainIS_Private(dm,subdm0,&iis0,&ois0);CHKERRQ(ierr); 437 if (iis) { 438 (*iis)[0] = iis0; 439 } else { 440 ierr = ISDestroy(&iis0);CHKERRQ(ierr); 441 } 442 if (ois) { 443 (*ois)[0] = ois0; 444 } else { 445 ierr = ISDestroy(&ois0);CHKERRQ(ierr); 446 } 447 if (subdm) (*subdm)[0] = subdm0; 448 if (names) (*names)[0] = 0; 449 PetscFunctionReturn(0); 450 } 451