C++ ISDestroy函数代码示例

PetscErrorCode  PCISDestroy(PC pc)
{
  PC_IS          *pcis = (PC_IS*)(pc->data);
  PetscErrorCode ierr;

  PetscFunctionBegin;
  ierr = ISDestroy(&pcis->is_B_local);CHKERRQ(ierr);
  ierr = ISDestroy(&pcis->is_I_local);CHKERRQ(ierr);
  ierr = ISDestroy(&pcis->is_B_global);CHKERRQ(ierr);
  ierr = ISDestroy(&pcis->is_I_global);CHKERRQ(ierr);
  ierr = MatDestroy(&pcis->A_II);CHKERRQ(ierr);
  ierr = MatDestroy(&pcis->A_IB);CHKERRQ(ierr);
  ierr = MatDestroy(&pcis->A_BI);CHKERRQ(ierr);
  ierr = MatDestroy(&pcis->A_BB);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->D);CHKERRQ(ierr);
  ierr = KSPDestroy(&pcis->ksp_N);CHKERRQ(ierr);
  ierr = KSPDestroy(&pcis->ksp_D);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec1_N);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec2_N);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec1_D);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec2_D);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec3_D);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec1_B);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec2_B);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec3_B);CHKERRQ(ierr);
  ierr = VecDestroy(&pcis->vec1_global);CHKERRQ(ierr);
  ierr = VecScatterDestroy(&pcis->global_to_D);CHKERRQ(ierr);
  ierr = VecScatterDestroy(&pcis->N_to_B);CHKERRQ(ierr);
  ierr = VecScatterDestroy(&pcis->global_to_B);CHKERRQ(ierr);
  ierr = PetscFree(pcis->work_N);CHKERRQ(ierr);
  if (pcis->ISLocalToGlobalMappingGetInfoWasCalled) {
    ierr = ISLocalToGlobalMappingRestoreInfo((ISLocalToGlobalMapping)0,&(pcis->n_neigh),&(pcis->neigh),&(pcis->n_shared),&(pcis->shared));CHKERRQ(ierr);
  }
  ierr = PetscObjectComposeFunction((PetscObject)pc,"PCISSetUseStiffnessScaling_C",NULL);CHKERRQ(ierr);
  ierr = PetscObjectComposeFunction((PetscObject)pc,"PCISSetSubdomainScalingFactor_C",NULL);CHKERRQ(ierr);
  ierr = PetscObjectComposeFunction((PetscObject)pc,"PCISSetSubdomainDiagonalScaling_C",NULL);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

int main(int argc,char **argv)
{
    PetscErrorCode ierr;
    PetscInt       i,n,first,step;
    IS             set;
    const PetscInt *indices;

    ierr = PetscInitialize(&argc,&argv,(char*)0,help);
    CHKERRQ(ierr);

    n     = 10;
    first = 3;
    step  = 2;

    /*
      Create stride index set, starting at 3 with a stride of 2
      Note each processor is generating its own index set
      (in this case they are all identical)
    */
    ierr = ISCreateStride(PETSC_COMM_SELF,n,first,step,&set);
    CHKERRQ(ierr);
    ierr = ISView(set,PETSC_VIEWER_STDOUT_SELF);
    CHKERRQ(ierr);

    /*
      Extract indices from set.
    */
    ierr = ISGetIndices(set,&indices);
    CHKERRQ(ierr);
    ierr = PetscPrintf(PETSC_COMM_WORLD,"Printing indices directly\n");
    CHKERRQ(ierr);
    for (i=0; i<n; i++) {
        ierr = PetscPrintf(PETSC_COMM_WORLD,"%D\n",indices[i]);
        CHKERRQ(ierr);
    }

    ierr = ISRestoreIndices(set,&indices);
    CHKERRQ(ierr);

    /*
        Determine information on stride
    */
    ierr = ISStrideGetInfo(set,&first,&step);
    CHKERRQ(ierr);
    if (first != 3 || step != 2) SETERRQ(PETSC_COMM_SELF,1,"Stride info not correct!\n");
    ierr = ISDestroy(&set);
    CHKERRQ(ierr);
    ierr = PetscFinalize();
    return 0;
}

int main(int argc,char **argv)
{
  PetscInt       n = 5;
  PetscErrorCode ierr;
  PetscMPIInt    rank,size;
  IS             ispetsc,isapp;
  AO             ao;

  ierr = PetscInitialize(&argc,&argv,(char*)0,help);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,"-n",&n,NULL);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);

  /* create the index sets */
  ierr = ISCreateStride(PETSC_COMM_WORLD,n,rank,size,&ispetsc);CHKERRQ(ierr);
  ierr = ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&isapp);CHKERRQ(ierr);

  /* create the application ordering */
  ierr = AOCreateBasicIS(isapp,ispetsc,&ao);CHKERRQ(ierr);


  ierr = AOView(ao,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);

  ierr = ISView(ispetsc,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = ISView(isapp,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = AOPetscToApplicationIS(ao,ispetsc);CHKERRQ(ierr);
  ierr = ISView(isapp,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = ISView(ispetsc,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);


  ierr = ISDestroy(&ispetsc);CHKERRQ(ierr);
  ierr = ISDestroy(&isapp);CHKERRQ(ierr);

  ierr = AODestroy(&ao);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}

int main(int argc,char **args)
{
  Mat            BAIJ,SBAIJ,*subBAIJ,*subSBAIJ;
  PetscViewer    viewer;
  char           file[PETSC_MAX_PATH_LEN];
  PetscBool      flg;
  PetscErrorCode ierr;
  PetscInt       n = 2,issize;
  PetscMPIInt    rank;
  IS             is,iss[2];

  PetscInitialize(&argc,&args,(char*)0,help);
  ierr = PetscOptionsGetString(NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr);
  ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer);CHKERRQ(ierr);
  ierr = MatCreate(PETSC_COMM_WORLD,&BAIJ);CHKERRQ(ierr);
  ierr = MatSetType(BAIJ,MATMPIBAIJ);CHKERRQ(ierr);
  ierr = MatLoad(BAIJ,viewer);CHKERRQ(ierr);
  ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
  ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer);CHKERRQ(ierr);
  ierr = MatCreate(PETSC_COMM_WORLD,&SBAIJ);CHKERRQ(ierr);
  ierr = MatSetType(SBAIJ,MATMPISBAIJ);CHKERRQ(ierr);
  ierr = MatLoad(SBAIJ,viewer);CHKERRQ(ierr);
  ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);

  ierr   = MatGetSize(BAIJ,&issize,0);CHKERRQ(ierr);
  issize = 9;
  ierr   = ISCreateStride(PETSC_COMM_SELF,issize,0,1,&is);CHKERRQ(ierr);
  iss[0] = is;iss[1] = is;
  ierr   = MatGetSubMatrices(BAIJ,n,iss,iss,MAT_INITIAL_MATRIX,&subBAIJ);CHKERRQ(ierr);
  ierr   = MatGetSubMatrices(SBAIJ,n,iss,iss,MAT_INITIAL_MATRIX,&subSBAIJ);CHKERRQ(ierr);

  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
#if defined(PETSC_USE_SOCKET_VIEWER)
  if (!rank) {
    ierr = MatView(subBAIJ[0],PETSC_VIEWER_SOCKET_SELF);CHKERRQ(ierr);
    ierr = MatView(subSBAIJ[0],PETSC_VIEWER_SOCKET_SELF);CHKERRQ(ierr);
  }
#endif

  /* Free data structures */
  ierr = ISDestroy(&is);CHKERRQ(ierr);
  ierr = MatDestroyMatrices(n,&subBAIJ);CHKERRQ(ierr);
  ierr = MatDestroyMatrices(n,&subSBAIJ);CHKERRQ(ierr);
  ierr = MatDestroy(&BAIJ);CHKERRQ(ierr);
  ierr = MatDestroy(&SBAIJ);CHKERRQ(ierr);

  ierr = PetscFinalize();
  return 0;
}

int main(int argc,char **argv)
{
  PetscErrorCode ierr;
  PetscInt       i,n,*indices;
  PetscInt       rank,size;
  IS             is,newis;

  ierr = PetscInitialize(&argc,&argv,(char*)0,help);CHKERRQ(ierr); 
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);

  /*
     Create IS
  */
  n = 4 + rank;
  ierr = PetscMalloc(n*sizeof(PetscInt),&indices);CHKERRQ(ierr);
  for (i=0; i<n; i++) {
    indices[i] = rank + i;
  }
  ierr = ISCreateGeneral(PETSC_COMM_WORLD,n,indices,PETSC_COPY_VALUES,&is);CHKERRQ(ierr);
  ierr = PetscFree(indices);CHKERRQ(ierr);

  /*
      Stick them together from all processors 
  */
  ierr = ISAllGather(is,&newis);CHKERRQ(ierr);

  if (!rank) {
    ierr = ISView(newis,PETSC_VIEWER_STDOUT_SELF);CHKERRQ(ierr);
  }

  ierr = ISDestroy(&newis);CHKERRQ(ierr);
  ierr = ISDestroy(&is);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}

static PetscErrorCode TaoDestroy_ASILS(Tao tao)
{
  TAO_SSLS       *ssls = (TAO_SSLS *)tao->data;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  ierr = VecDestroy(&ssls->ff);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->dpsi);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->da);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->db);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->w);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->t1);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->t2);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->r1);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->r2);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->r3);CHKERRQ(ierr);
  ierr = VecDestroy(&ssls->dxfree);CHKERRQ(ierr);
  ierr = MatDestroy(&ssls->J_sub);CHKERRQ(ierr);
  ierr = MatDestroy(&ssls->Jpre_sub);CHKERRQ(ierr);
  ierr = ISDestroy(&ssls->fixed);CHKERRQ(ierr);
  ierr = ISDestroy(&ssls->free);CHKERRQ(ierr);
  ierr = PetscFree(tao->data);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

PetscErrorCode MatPartitioningHierarchical_AssembleSubdomain(Mat adj,IS destination,Mat *sadj, ISLocalToGlobalMapping *mapping)
{
  IS              irows,icols;
  PetscInt        irows_ln;
  PetscMPIInt     rank;
  const PetscInt *irows_indices;
  MPI_Comm        comm;
  PetscErrorCode  ierr;

  PetscFunctionBegin;
  ierr = PetscObjectGetComm((PetscObject)adj,&comm);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
  /* figure out where data comes from  */
  ierr = ISBuildTwoSided(destination,NULL,&irows);CHKERRQ(ierr);
  ierr = ISDuplicate(irows,&icols);CHKERRQ(ierr);
  ierr = ISGetLocalSize(irows,&irows_ln);CHKERRQ(ierr);
  ierr = ISGetIndices(irows,&irows_indices);CHKERRQ(ierr);
  ierr = ISLocalToGlobalMappingCreate(comm,1,irows_ln,irows_indices,PETSC_COPY_VALUES,mapping);CHKERRQ(ierr);
  ierr = ISRestoreIndices(irows,&irows_indices);CHKERRQ(ierr);
  ierr = MatGetSubMatrices(adj,1,&irows,&icols,MAT_INITIAL_MATRIX,&sadj);CHKERRQ(ierr);
  ierr = ISDestroy(&irows);CHKERRQ(ierr);
  ierr = ISDestroy(&icols);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

static PetscErrorCode ISView_Block(IS is, PetscViewer viewer)
{
  IS_Block       *sub = (IS_Block*)is->data;
  PetscErrorCode ierr;
  PetscInt       i,bs,n,*idx = sub->idx;
  PetscBool      iascii;

  PetscFunctionBegin;
  ierr = PetscLayoutGetBlockSize(is->map, &bs);CHKERRQ(ierr);
  ierr = PetscLayoutGetLocalSize(is->map, &n);CHKERRQ(ierr);
  n   /= bs;
  ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
  if (iascii) {
    PetscViewerFormat fmt;

    ierr = PetscViewerGetFormat(viewer,&fmt);CHKERRQ(ierr);
    if (fmt == PETSC_VIEWER_ASCII_MATLAB) {
      IS             ist;
      const char     *name;
      const PetscInt *idx;
      PetscInt       n;

      ierr = PetscObjectGetName((PetscObject)is,&name);CHKERRQ(ierr);
      ierr = ISGetLocalSize(is,&n);CHKERRQ(ierr);
      ierr = ISGetIndices(is,&idx);CHKERRQ(ierr);
      ierr = ISCreateGeneral(PetscObjectComm((PetscObject)is),n,idx,PETSC_USE_POINTER,&ist);CHKERRQ(ierr);
      ierr = PetscObjectSetName((PetscObject)ist,name);CHKERRQ(ierr);
      ierr = ISView(ist,viewer);CHKERRQ(ierr);
      ierr = ISDestroy(&ist);CHKERRQ(ierr);
      ierr = ISRestoreIndices(is,&idx);CHKERRQ(ierr);
    } else {
      ierr = PetscViewerASCIIPushSynchronized(viewer);CHKERRQ(ierr);
      if (is->isperm) {
        ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Block Index set is permutation\n");CHKERRQ(ierr);
      }
      ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Block size %D\n",bs);CHKERRQ(ierr);
      ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Number of block indices in set %D\n",n);CHKERRQ(ierr);
      ierr = PetscViewerASCIISynchronizedPrintf(viewer,"The first indices of each block are\n");CHKERRQ(ierr);
      for (i=0; i<n; i++) {
        ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Block %D Index %D\n",i,idx[i]);CHKERRQ(ierr);
      }
      ierr = PetscViewerFlush(viewer);CHKERRQ(ierr);
      ierr = PetscViewerASCIIPopSynchronized(viewer);CHKERRQ(ierr);
    }
  }
  PetscFunctionReturn(0);
}

PetscErrorCode DMDestroy_AKKT(DM dm) {
  DM_AKKT *kkt = (DM_AKKT*)(dm->data);
  PetscInt i;
  PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = MatDestroy(&(kkt->Aff));       CHKERRQ(ierr);
  ierr = DMDestroy(&(kkt->dm));         CHKERRQ(ierr);
  for(i = 0; i < 2; ++i) {
    ierr = DMDestroy(&(kkt->dmf[i]));   CHKERRQ(ierr);
    ierr = ISDestroy(&(kkt->isf[i]));   CHKERRQ(ierr);
    ierr = PetscFree(kkt->names[i]);    CHKERRQ(ierr);
  }
  ierr = DMDestroy(&(kkt->cdm));        CHKERRQ(ierr);
  ierr = MatDestroy(&(kkt->Pfc));       CHKERRQ(ierr);
  ierr = PetscFree(kkt); CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

static PetscErrorCode PCDestroy_Redistribute(PC pc)
{
  PC_Redistribute *red = (PC_Redistribute*)pc->data;
  PetscErrorCode  ierr;

  PetscFunctionBegin;
  ierr = VecScatterDestroy(&red->scatter);CHKERRQ(ierr);
  ierr = ISDestroy(&red->is);CHKERRQ(ierr);
  ierr = VecDestroy(&red->b);CHKERRQ(ierr);
  ierr = VecDestroy(&red->x);CHKERRQ(ierr);
  ierr = KSPDestroy(&red->ksp);CHKERRQ(ierr);
  ierr = VecDestroy(&red->work);CHKERRQ(ierr);
  ierr = PetscFree(red->drows);CHKERRQ(ierr);
  ierr = PetscFree(red->diag);CHKERRQ(ierr);
  ierr = PetscFree(pc->data);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

PetscErrorCode DMDestroy_SNESVI(DM_SNESVI *dmsnesvi)
{
  PetscErrorCode ierr;

  PetscFunctionBegin;
  /* reset the base methods in the DM object that were changed when the DM_SNESVI was reset */
  dmsnesvi->dm->ops->createinterpolation   = dmsnesvi->createinterpolation;
  dmsnesvi->dm->ops->coarsen            = dmsnesvi->coarsen;
  dmsnesvi->dm->ops->createglobalvector = dmsnesvi->createglobalvector;
  /* need to clear out this vectors because some of them may not have a reference to the DM
    but they are counted as having references to the DM in DMDestroy() */
  ierr = DMClearGlobalVectors(dmsnesvi->dm);CHKERRQ(ierr);

  ierr = ISDestroy(&dmsnesvi->inactive);CHKERRQ(ierr);
  ierr = PetscFree(dmsnesvi);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

/*
  There is a problem here with uninterpolated meshes. The index in numDof[] is not dimension in this case,
  but sieve depth.
*/
PetscErrorCode SetupSection(DM dm, AppCtx *user)
{
  PetscSection    section;
  const PetscInt  numFields           = NUM_FIELDS;
  PetscInt        dim                 = user->dim;
  PetscInt        numBC               = 0;
  PetscInt        bcFields[1]         = {0};
  IS              bcPoints[1]         = {NULL};
  PetscInt        numComp[NUM_FIELDS];
  const PetscInt *numFieldDof[NUM_FIELDS];
  PetscInt       *numDof;
  PetscInt        f, d;
  PetscErrorCode  ierr;

  PetscFunctionBeginUser;
  ierr = PetscFEGetNumComponents(user->fe[0], &numComp[0]);CHKERRQ(ierr);
  ierr = PetscFEGetNumComponents(user->fe[1], &numComp[1]);CHKERRQ(ierr);
  ierr = PetscFEGetNumDof(user->fe[0], &numFieldDof[0]);CHKERRQ(ierr);
  ierr = PetscFEGetNumDof(user->fe[1], &numFieldDof[1]);CHKERRQ(ierr);
  ierr = PetscMalloc(NUM_FIELDS*(dim+1) * sizeof(PetscInt), &numDof);CHKERRQ(ierr);
  for (f = 0; f < NUM_FIELDS; ++f) {
    for (d = 0; d <= dim; ++d) {
      numDof[f*(dim+1)+d] = numFieldDof[f][d];
    }
  }
  for (f = 0; f < numFields; ++f) {
    for (d = 1; d < dim; ++d) {
      if ((numDof[f*(dim+1)+d] > 0) && !user->interpolate) SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONG, "Mesh must be interpolated when unknowns are specified on edges or faces.");
    }
  }
  if (user->bcType == DIRICHLET) {
    numBC = 1;
    ierr  = DMPlexGetStratumIS(dm, "marker", 1, &bcPoints[0]);CHKERRQ(ierr);
  }
  ierr = DMPlexCreateSection(dm, dim, numFields, numComp, numDof, numBC, bcFields, bcPoints, &section);CHKERRQ(ierr);
  ierr = PetscSectionSetFieldName(section, 0, "velocity");CHKERRQ(ierr);
  ierr = PetscSectionSetFieldName(section, 1, "pressure");CHKERRQ(ierr);
  ierr = DMSetDefaultSection(dm, section);CHKERRQ(ierr);
  ierr = PetscSectionDestroy(&section);CHKERRQ(ierr);
  if (user->bcType == DIRICHLET) {
    ierr = ISDestroy(&bcPoints[0]);CHKERRQ(ierr);
  }
  ierr = PetscFree(numDof);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

/*
  DMLabelMakeInvalid_Private - Transfer stratum data from the sorted list format to the hash format

  Input parameter:
+ label - The DMLabel
- v - The stratum value

  Output parameter:
. label - The DMLabel with stratum in hash format

  Level: developer

.seealso: DMLabelCreate()
*/
static PetscErrorCode DMLabelMakeInvalid_Private(DMLabel label, PetscInt v)
{
  PETSC_UNUSED PetscHashIIter ret, iter;
  PetscInt                    p;
  const PetscInt              *points;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  if (!label->validIS[v]) PetscFunctionReturn(0);
  if (label->points[v]) {
    ierr = ISGetIndices(label->points[v],&points);CHKERRQ(ierr);
    for (p = 0; p < label->stratumSizes[v]; ++p) PetscHashIPut(label->ht[v], points[p], ret, iter);
    ierr = ISRestoreIndices(label->points[v],&points);CHKERRQ(ierr);
    ierr = ISDestroy(&(label->points[v]));CHKERRQ(ierr);
  }
  label->validIS[v] = PETSC_FALSE;
  PetscFunctionReturn(0);
}

static PetscErrorCode test5(DM dm, AppCtx *options)
{
  IS             cells;
  Vec            locX, locX_t, locA;
  PetscScalar    *u, *u_t, *a;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  locX_t = NULL;
  locA = NULL;
  ierr = ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &cells);CHKERRQ(ierr);
  ierr = DMGetLocalVector(dm, &locX);CHKERRQ(ierr);
  ierr = DMPlexGetCellFields(    dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr);
  ierr = DMPlexRestoreCellFields(dm, cells, locX, locX_t, locA, &u, &u_t, &a);CHKERRQ(ierr);
  ierr = DMRestoreLocalVector(dm, &locX);CHKERRQ(ierr);
  ierr = ISDestroy(&cells);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

int main(int argc,char **argv)
{
  PetscErrorCode ierr;
  PetscInt       step = 2;
  IS             is;

  ierr = PetscInitialize(&argc,&argv,(char*)0,help);CHKERRQ(ierr);

  ierr = PetscOptionsGetInt(NULL,"-step",&step,NULL);CHKERRQ(ierr);
  ierr = ISCreateStride(PETSC_COMM_SELF,10,0,step,&is);CHKERRQ(ierr);

  ierr = ISToGeneral(is);CHKERRQ(ierr);

  ierr = ISDestroy(&is);CHKERRQ(ierr);

  ierr = PetscFinalize();
  return 0;
}