C++ ProblemDomain::contains方法代码示例

void
SlabService::fillGraph(BaseFab<int>&      a_regIrregCovered,
                       Vector<IrregNode>& a_nodes,
                       const Box&         a_validRegion,
                       const Box&         a_ghostRegion,
                       const ProblemDomain&         a_domain,
                       const RealVect&    a_origin,
                       const Real&        a_dx) const
{
  Box grownCovBox = grow(m_coveredRegion, 1);
  grownCovBox &= a_ghostRegion;
  IntVectSet ivsIrreg(grownCovBox);
  ivsIrreg -= m_coveredRegion;
  ivsIrreg &= a_domain;

  //set regirregcoverred flags
  CH_assert(a_regIrregCovered.box().contains(a_ghostRegion));

  //set every cell to regular
  a_regIrregCovered.setVal(1);
  for (BoxIterator bit(a_ghostRegion); bit.ok(); ++bit)
    {
      if (m_coveredRegion.contains(bit()))
        {
          //set covered cells to -1
          a_regIrregCovered(bit(), 0) = -1;
        }
      else if (ivsIrreg.contains(bit()))
        {
          //set irreg cells to 0
          a_regIrregCovered(bit(), 0) = 0;
        }
    }

  //now loop through irreg cells and make Nodes for them
  a_nodes.resize(0);
  for (IVSIterator ivsit(ivsIrreg); ivsit.ok(); ++ivsit)
    {
      const IntVect& iv = ivsit();
      if (a_validRegion.contains(iv))
        {
          IrregNode node;
          //first the obvious
          node.m_cell = iv;
          node.m_volFrac = 1.0;
          node.m_cellIndex = 0;
          node.m_volCentroid    = RealVect::Zero;
          //any time the next cell over is in the covered
          //region, there is no face.  If there is a cell
          //but it is outside the domain, the arc=-1 to signify
          //a boundary face.  Otherwise, the arc is -2 if it
          //is to a regular cell and 0 if it is to an irregular cell
          for (int idir = 0; idir < SpaceDim; idir++)
            {
              for (SideIterator sit; sit.ok(); ++sit)
                {
                  int arcIndex = node.index(idir, sit());
                  Vector<int>&      arcs     = node.m_arc[arcIndex];
                  Vector<Real>&     areaFracs= node.m_areaFrac[arcIndex];
                  Vector<RealVect>& faceCents= node.m_faceCentroid[arcIndex];
                  IntVect otherIV = iv + sign(sit())*BASISV(idir);
                  if (m_coveredRegion.contains(otherIV))
                    {
                      //do nothing, covered face. leave the vector empty
                    }
                  else
                    {
                      int otherCellIndex;
                      if (ivsIrreg.contains(otherIV))
                        {
                          //arc irregular cell inside the domain
                          otherCellIndex = 0;
                        }
                      else if (!a_domain.contains(otherIV))
                        {
                          //boundary face
                          otherCellIndex = -1;
                        }
                      else
                        {
                          //arc to regular cell
                          otherCellIndex = -2;
                        }
                      arcs.push_back(otherCellIndex);
                      Real     areaFrac = 1.0;
                      RealVect faceCent = RealVect::Zero;

                      areaFracs.push_back(areaFrac);
                      faceCents.push_back(faceCent);
                    } //end otherIV not covered
                }
            }
          //the boundary centroid and normal
          //depend on which side is covered

          for (int idir = 0; idir < SpaceDim; idir++)
            {
              for (SideIterator sit; sit.ok(); ++sit)
                {
                  int arcIndex = node.index(idir, sit());
//.........这里部分代码省略.........

void
EBMGInterp::fillGhostCellsPWC(LevelData<EBCellFAB>& a_data,
                              const EBISLayout&     a_ebisl,
                              const ProblemDomain&  a_dom)
{
  for (int idir = 0; idir < SpaceDim; idir++)
    {
      if (m_ghost[idir] < 1)
        {
          MayDay::Error("EBMGInterp:I need a ghost cell for linear interpolation");
        }
    }
  //extrapolate to every ghost cell
  const DisjointBoxLayout& dbl = a_data.disjointBoxLayout();
  for (DataIterator dit = dbl.dataIterator(); dit.ok(); ++dit)
    {
      const Box& grid = dbl.get(dit());
      const EBGraph& ebgraph = a_ebisl[dit()].getEBGraph();
      for (int idir = 0; idir < SpaceDim; idir++)
        {
          for (SideIterator sit; sit.ok(); ++sit)
            {
              Side::LoHiSide flipSide = flip(sit());
              Box ghostBox = adjCellBox(grid, idir, sit(), 1);

              for (BoxIterator boxit(ghostBox); boxit.ok(); ++boxit)
                {
                  if (a_dom.contains(boxit()))
                    {
                      Vector<VolIndex> vofs = ebgraph.getVoFs(boxit());
                      for (int ivof = 0; ivof < vofs.size(); ivof++)
                        {
                          Real extrapVal = 0;
                          Vector<FaceIndex> faces = ebgraph.getFaces(vofs[ivof], idir, flipSide);
                          for (int ivar = 0; ivar < a_data.nComp(); ivar++)
                            {
                              for (int iface = 0; iface < faces.size(); iface++)
                                {
                                  const VolIndex& flipVoF = faces[iface].getVoF(flipSide);
                                  extrapVal += a_data[dit()](flipVoF, ivar);
                                }
                              if (faces.size() > 1) extrapVal /= faces.size();
                              a_data[dit()](vofs[ivof], ivar) = extrapVal;
                            }
                        }

                    }
                  else
                    {
                      //just put in the single valued part of the data holder something sensible
                      //if not inside domain
                      int isign = sign(flipSide);
                      BaseFab<Real>& bfdata = a_data[dit()].getSingleValuedFAB();
                      IntVect otherIV = boxit() + isign*BASISV(idir);
                      for (int ivar = 0; ivar < a_data.nComp(); ivar++)
                        {
                          bfdata(boxit(), ivar) = bfdata(otherIV, ivar);
                        }
                    }
                }
            }
        }
    }
  //do exchange to overwrite ghost cells where there exists real data
  a_data.exchange();
}

void
EBCompositeMACProjector::
correctTangentialVelocity(EBFaceFAB&        a_velocity,
                          const EBFaceFAB&  a_gradient,
                          const Box&        a_grid,
                          const EBISBox&    a_ebisBox,
                          const IntVectSet& a_cfivs)
{
  CH_TIME("EBCompositeMACProjector::correctTangentialVelocity");
  int velDir = a_velocity.direction();
  int gradDir = a_gradient.direction();
  //veldir is the face on which the velocity lives
  //graddir is the direction of the component
  //and the face on which the mac gradient lives
  CH_assert(velDir != gradDir);
  CH_assert(a_velocity.nComp() == 1);
  CH_assert(a_gradient.nComp() == 1);

  //updating in place in fortran so we have to save a acopy
  EBFaceFAB velSave(a_ebisBox, a_velocity.getCellRegion(),  velDir, 1);
  velSave.copy(a_velocity);

  //interior box is the box where all of the stencil can be reached
  //without going out of the domain
  ProblemDomain domainBox = a_ebisBox.getDomain();
  Box interiorBox = a_grid;
  interiorBox.grow(1);
  interiorBox &= domainBox;
  interiorBox.grow(-1);
  Box interiorFaceBox = surroundingNodes(interiorBox, velDir);

  if (!interiorBox.isEmpty())
    {
      BaseFab<Real>&       regVel  = a_velocity.getSingleValuedFAB();
      const BaseFab<Real>& regGrad = a_gradient.getSingleValuedFAB();
      FORT_REGCORRECTTANVEL(CHF_FRA1(regVel,0),
                            CHF_CONST_FRA1(regGrad,0),
                            CHF_BOX(interiorFaceBox),
                            CHF_INT(velDir),
                            CHF_INT(gradDir));
    }

  //do only irregular and boundary cells pointwise
  IntVectSet ivsIrreg = a_ebisBox.getIrregIVS(a_grid);
  IntVectSet ivsGrid(a_grid);
  ivsGrid -= interiorBox;
  ivsGrid |= ivsIrreg;
  FaceIterator faceit(ivsGrid, a_ebisBox.getEBGraph(), velDir, FaceStop::SurroundingWithBoundary);
  for (faceit.reset(); faceit.ok(); ++faceit)
    {
      //average neighboring grads in grad dir direction
      int numGrads = 0;
      Real gradAve = 0.0;
      const FaceIndex& velFace = faceit();
      for (SideIterator sitVel; sitVel.ok(); ++sitVel)
        {
          const VolIndex& vofSide = velFace.getVoF(sitVel());
          const IntVect&   ivSide = vofSide.gridIndex();
          //do not include stuff over coarse-fine interface and inside the domain
          //cfivs includes cells just outside domain
          if (!a_cfivs.contains(ivSide) && domainBox.contains(ivSide))
            {
              for (SideIterator sitGrad; sitGrad.ok(); ++sitGrad)
                {
                  Vector<FaceIndex> gradFaces = a_ebisBox.getFaces(vofSide, gradDir, sitGrad());
                  for (int iface = 0; iface < gradFaces.size(); iface++)
                    {
                      if (!gradFaces[iface].isBoundary())
                        {
                          numGrads++;
                          gradAve += a_gradient(gradFaces[iface], 0);
                        }
                    }
                }//end loop over gradient sides
            }//end cfivs check
          else
            {
              // inside coarse/fine interface or at domain boundary. extrapolate from neighboring vofs to get the gradient
              const Side::LoHiSide inSide    = flip(sitVel());
              const VolIndex&      inSideVof = velFace.getVoF(inSide);
              const IntVect&       inSideIV  = inSideVof.gridIndex();

              IntVect inSideFarIV = inSideIV;
              inSideFarIV[velDir] += sign(inSide);
              if (domainBox.contains(inSideIV) && domainBox.contains(inSideFarIV))
                {
                  Vector<VolIndex> farVofs = a_ebisBox.getVoFs(inSideFarIV);
                  if (farVofs.size()  == 1)
                    {
                      const VolIndex& inSideFarVof = farVofs[0];
                      for (SideIterator sitGrad; sitGrad.ok(); ++sitGrad)
                        {
                          //get the grad for the face adjoining inSideVof on the sitGrad side, in the gradDir direction
                          Vector<FaceIndex> gradFaces    = a_ebisBox.getFaces(inSideVof   , gradDir, sitGrad());
                          Vector<FaceIndex> gradFarFaces = a_ebisBox.getFaces(inSideFarVof, gradDir, sitGrad());
                          if ( (gradFaces.size() == 1) && (gradFarFaces.size() == 1) )
                            {
                              // if ( (!gradFaces[0].isBoundary()) && (!gradFarFaces[0].isBoundary()) )
                              //   {
                                  const Real& inSideGrad    = a_gradient(gradFaces[0], 0);
//.........这里部分代码省略.........

void oldLoHiCenter(Box&                 a_loBox,
                   int&                 a_hasLo,
                   Box&                 a_hiBox,
                   int&                 a_hasHi,
                   Box&                 a_centerBox,
                   Box&                 a_entireBox,
                   const Box&           a_inBox,
                   const ProblemDomain& a_domain,
                   const int&           a_dir)
{
    // Make a copy of the input box which can be modified
    Box inBox = a_inBox;
    inBox &= a_domain;

    // The centered difference box is always one smaller in a_dir
    a_centerBox = inBox;
    a_centerBox.grow(a_dir,-1);

    // The union of all the output boxes start off equal to the center
    // difference portion on the input box (intersected with the domain)
    a_entireBox = a_centerBox;

    // See if this chops off the high side of the input box
    Box tmp = a_inBox;
    tmp.shift(a_dir,-1);
    tmp &= a_domain;
    tmp.shift(a_dir,1);

    // If so, set up the high, one-sided difference box, a_hiBox, and expand
    // the entire box to include it
    if (!a_domain.contains(tmp))
    {
        a_hasHi = 1;
        tmp.shift(a_dir,-2);
        a_hiBox = adjCellHi(tmp,a_dir);
        a_entireBox.growHi(a_dir,1);
    }
    else
    {
        a_hasHi = 0;
        a_hiBox = Box();
    }

    // See if this chops off the low side of the input box
    tmp = a_inBox;
    tmp.shift(a_dir,1);
    tmp &= a_domain;
    tmp.shift(a_dir,-1);

    // If so, set up the low, one-sided difference box, a_loBox, and expand
    // the entire box to include it
    if (!a_domain.contains(tmp))
    {
        a_hasLo = 1;
        tmp.shift(a_dir,2);
        a_loBox = adjCellLo(tmp,a_dir);
        a_entireBox.growLo(a_dir,1);
    }
    else
    {
        a_hasLo = 0;
        a_loBox = Box();
    }

    // Make some simple sanity checks
    CH_assert(a_entireBox.contains(a_centerBox));

    if (a_hasLo == 1)
    {
        CH_assert(a_entireBox.contains(a_loBox));
    }

    if (a_hasHi == 1)
    {
        CH_assert(a_entireBox.contains(a_hiBox));
    }
}