C++ Interval::begin方法代碼示例

 virtual void linearIn(EBCellFAB& arg,  void* buf, const Box& R,
                       const Interval& comps) const
 {
   EBCellFAB tmp;
   tmp.clone(arg);
   tmp.linearIn(buf, R, comps);
   arg.plus(tmp, R, comps.begin(), comps.begin(), comps.size());
 }

void
EBMGAverage::averageFAB(EBCellFAB&       a_coar,
                        const Box&       a_boxCoar,
                        const EBCellFAB& a_refCoar,
                        const DataIndex& a_datInd,
                        const Interval&  a_variables) const
{
  CH_TIMERS("EBMGAverage::average");
  CH_TIMER("regular_average", t1);
  CH_TIMER("irregular_average", t2);
  CH_assert(isDefined());

  const Box& coarBox = a_boxCoar;

  //do all cells as if they were regular
  Box refBox(IntVect::Zero, IntVect::Zero);
  refBox.refine(m_refRat);
  int numFinePerCoar = refBox.numPts();

  BaseFab<Real>& coarRegFAB =             a_coar.getSingleValuedFAB();
  const BaseFab<Real>& refCoarRegFAB = a_refCoar.getSingleValuedFAB();

  //set to zero because the fortran is a bit simpleminded
  //and does stuff additively
  a_coar.setVal(0.);
  CH_START(t1);
  for (int comp = a_variables.begin();  comp <= a_variables.end(); comp++)
    {
      FORT_REGAVERAGE(CHF_FRA1(coarRegFAB,comp),
                      CHF_CONST_FRA1(refCoarRegFAB,comp),
                      CHF_BOX(coarBox),
                      CHF_BOX(refBox),
                      CHF_CONST_INT(numFinePerCoar),
                      CHF_CONST_INT(m_refRat));
    }
  CH_STOP(t1);

  //this is really volume-weighted averaging even though it does
  //not look that way.

  //so (in the traditional sense) we want to preserve
  //rhoc * volc = sum(rhof * volf)
  //this translates to
  //volfrac_C * rhoC = (1/numFinePerCoar)(sum(volFrac_F * rhoF))
  //but the data input to this routine is all kappa weigthed so
  //the volumefractions have already been multiplied
  //which means
  // rhoC = (1/numFinePerCoar)(sum(rhoF))
  //which is what this does

  CH_START(t2);
  for (int comp = a_variables.begin();  comp <= a_variables.end(); comp++)
    {
      m_averageEBStencil[a_datInd]->apply(a_coar, a_refCoar, false, comp);
    }
  CH_STOP(t2);

}

 void op(EBCellFAB& dest,
         const Box& RegionFrom,
         const Interval& Cdest,
         const Box& RegionTo,
         const EBCellFAB& src,
         const Interval& Csrc) const
 {
   dest.plus(src, RegionFrom, Csrc.begin(), Cdest.begin(), Cdest.size());
 }

void
MappedLevelFluxRegister::incrementCoarse(const FArrayBox& a_coarseFlux,
                                   Real a_scale,
                                   const DataIndex& a_coarseDataIndex,
                                   const Interval& a_srcInterval,
                                   const Interval& a_dstInterval,
                                   int a_dir,
                                   Side::LoHiSide a_sd)
{
    CH_assert(isDefined());
    if (!(m_isDefined & FluxRegCoarseDefined)) return;
    CH_TIME("MappedLevelFluxRegister::incrementCoarse");

    const Vector<Box>& intersect =
        m_coarseLocations[a_dir + a_sd * CH_SPACEDIM][a_coarseDataIndex];

    FArrayBox& coarse = m_coarFlux[a_coarseDataIndex];

    // We cast away the constness in a_coarseFlux for the scope of this function. This
    // should be acceptable, since at the end of the day there is no change to it. -JNJ
    FArrayBox& coarseFlux = const_cast<FArrayBox&>(a_coarseFlux); // Muhahaha.
    coarseFlux.shiftHalf(a_dir, sign(a_sd));
    Real scale = -sign(a_sd) * a_scale;
    int s = a_srcInterval.begin();
    int d = a_dstInterval.begin();
    int size = a_srcInterval.size();

    for (int b = 0; b < intersect.size(); ++b) {
        const Box& box = intersect[b];
        Vector<Real> regbefore(coarse.nComp());
        Vector<Real> regafter(coarse.nComp());
        if (s_verbose && (a_dir == debugdir) && box.contains(ivdebnoeb)) {
            for (int ivar = 0; ivar < coarse.nComp(); ivar++) {
                regbefore[ivar] = coarse(ivdebnoeb, ivar);
            }
        }

        coarse.plus(coarseFlux, box, box, scale, s, d, size);

        if (s_verbose && (a_dir == debugdir) && box.contains(ivdebnoeb)) {
            for (int ivar = 0; ivar < coarse.nComp(); ivar++) {
                regafter[ivar] = coarse(ivdebnoeb, ivar);
            }

            pout() << "levelfluxreg::incrementCoar: scale = " << scale << ", ";
            for (int ivar = 0; ivar < coarse.nComp(); ivar++) {
                pout() << " input flux = " << coarseFlux(ivdebnoeb, ivar) << ", ";
                pout() << " reg before = " <<             regbefore[ivar] << ", ";
                pout() << " reg after  = " <<                regafter[ivar] << ", ";
            }
            pout() << endl;

        }
    }

    coarseFlux.shiftHalf(a_dir, - sign(a_sd));
}

void
EBCoarsen::coarsenIrreg(EBCellFAB&       a_coar,
                        const EBCellFAB& a_fine,
                        const DataIndex& a_dit,
                        const Interval&  a_variables)
{
  const BaseIVFAB<VoFStencil>& stenBaseIVFAB = m_coarsenStencil[a_dit];
  VoFIterator& vofit = m_vofIt[a_dit];
  for (vofit.reset(); vofit.ok(); ++vofit)
    {
      const VolIndex&    vofCoar   = vofit();
      const VoFStencil&  stencil   = stenBaseIVFAB(vofCoar,0);

      for (int icomp=a_variables.begin();icomp ==a_variables.end();icomp++)
        {
          //coarsen irreg fine vofs to coarse
          //  compute the coarsening by using fine data
          Real phi = 0.0;
          for (int i = 0; i < stencil.size(); ++i )
            {
              const Real&      weight = stencil.weight(i);
              const VolIndex& vofFine = stencil.vof(i);
              const Real&        phiF = a_fine(vofFine,icomp);
              phi += weight * phiF;
            }

          //set the coarse value
          a_coar(vofCoar,icomp) = phi;
        }
    }
}

// ---------------------------------------------------------
// 28 March 2003:
// This is called by other functions, and should not be called directly.
Real
integral(const BoxLayoutData<NodeFArrayBox>& a_layout,
         const Real a_dx,
         const Interval& a_interval,
         bool a_verbose)
{
  Real integralTotal = 0.;

  for (DataIterator it = a_layout.dataIterator(); it.ok(); ++it)
    {
      const NodeFArrayBox& thisNfab = a_layout[it()];
      const Box& thisBox(a_layout.box(it())); // CELL-centered
      Real thisNfabIntegral =
        integral(thisNfab, a_dx, thisBox,
                 a_interval.begin(), a_interval.size());
      integralTotal += thisNfabIntegral;
    }
# ifdef CH_MPI
  Real recv;
  // add up
  int result = MPI_Allreduce(&integralTotal, &recv, 1, MPI_CH_REAL,
                             MPI_SUM, Chombo_MPI::comm);
  if (result != MPI_SUCCESS)
    { //bark!!!
      MayDay::Error("sorry, but I had a communication error on integral");
    }
  integralTotal = recv;
# endif
  return integralTotal;
}

// ---------------------------------------------------------
// 27 March 2003:
// This is called by other functions, and should not be called directly.
Real
maxnorm(const BoxLayoutData<NodeFArrayBox>& a_layout,
        const Interval& a_interval,
        bool a_verbose)
{
  Real normTotal = 0.;
  // a_p == 0:  max norm
  for (DataIterator it = a_layout.dataIterator(); it.ok(); ++it)
    {
      const Box& thisBox(a_layout.box(it())); // CELL-centered
      const NodeFArrayBox& thisNfab = a_layout[it()];
      Real thisNfabNorm =
        maxnorm(thisNfab, thisBox, a_interval.begin(), a_interval.size());
      if (a_verbose)
        cout << "maxnorm(" << thisBox << ") = " << thisNfabNorm << endl;
      normTotal = Max(normTotal, thisNfabNorm);
    }
# ifdef CH_MPI
  Real recv;
  int result = MPI_Allreduce(&normTotal, &recv, 1, MPI_CH_REAL,
                             MPI_MAX, Chombo_MPI::comm);
  if (result != MPI_SUCCESS)
    { //bark!!!
      MayDay::Error("sorry, but I had a communication error on maxnorm");
    }
  normTotal = recv;
# endif
  return normTotal;
}

void
LevelFluxRegisterEdge::incrementFine(
                                     FArrayBox& a_fineFlux,
                                     Real a_scale,
                                     const DataIndex& a_fineDataIndex,
                                     const Interval& a_srcInterval,
                                     const Interval& a_dstInterval)
{
  CH_assert(isDefined());
  CH_assert(!a_fineFlux.box().isEmpty());
  CH_assert(a_srcInterval.size() == a_dstInterval.size());
  CH_assert(a_srcInterval.begin() >= 0);
  CH_assert(a_srcInterval.end() < a_fineFlux.nComp());
  CH_assert(a_dstInterval.begin() >= 0);
  CH_assert(a_dstInterval.end() < m_nComp);

  int edgeDir = -1;
  for (int sideDir = 0; sideDir<SpaceDim; sideDir++)
    {
      if (a_fineFlux.box().type(sideDir) == IndexType::CELL)
        {
          edgeDir = sideDir;
        }
    }
  CH_assert(edgeDir >= 0);
  CH_assert(edgeDir < SpaceDim);

  for (int faceDir=0; faceDir<SpaceDim; faceDir++)
    {
      if (faceDir != edgeDir)
        {

          SideIterator sit;
          for (sit.begin(); sit.ok(); ++sit)
            {
              incrementFine(a_fineFlux,
                            a_scale,
                            a_fineDataIndex,
                            a_srcInterval,
                            a_dstInterval,
                            faceDir,
                            sit());
            }
        }
    }
}

// ---------------------------------------------------------
// 7 Dec 2005
Real
maxnorm(const BoxLayoutData<NodeFArrayBox>& a_layout,
        const LevelData<NodeFArrayBox>& a_mask,
        const ProblemDomain& a_domain,
        const Interval& a_interval,
        bool a_verbose)
{
  Real normTotal = 0.;
  // a_p == 0:  max norm
  int ncomp = a_interval.size();
  for (DataIterator it = a_layout.dataIterator(); it.ok(); ++it)
    {
      const NodeFArrayBox& thisNfab = a_layout[it()];
      const FArrayBox& dataFab = thisNfab.getFab();
      const FArrayBox& maskFab = a_mask[it()].getFab();
      const Box& thisBox(a_layout.box(it())); // CELL-centered
      NodeFArrayBox dataMasked(thisBox, ncomp);
      FArrayBox& dataMaskedFab = dataMasked.getFab();
      dataMaskedFab.copy(dataFab);
      // dataMaskedFab *= maskFab;
      for (int comp = a_interval.begin(); comp <= a_interval.end(); comp++)
        {
          // Set dataMaskedFab[comp] *= maskFab[0].
          dataMaskedFab.mult(maskFab, 0, comp);
        }
      Real thisNfabNorm =
        maxnorm(dataMasked, thisBox, a_interval.begin(), a_interval.size());
      if (a_verbose)
        cout << "maxnorm(" << thisBox << ") = " << thisNfabNorm << endl;
       normTotal = Max(normTotal, thisNfabNorm);
    }
# ifdef CH_MPI
  Real recv;
  // add up (a_p is not 0)
  int result = MPI_Allreduce(&normTotal, &recv, 1, MPI_CH_REAL,
                             MPI_MAX, Chombo_MPI::comm);
  if (result != MPI_SUCCESS)
    { //bark!!!
      MayDay::Error("sorry, but I had a communication error on norm");
    }
  normTotal = recv;
# endif
  return normTotal;
}

void
LevelFluxRegisterEdge::incrementCoarse(FArrayBox& a_coarseFlux,
                                       Real a_scale,
                                       const DataIndex& a_coarseDataIndex,
                                       const Interval& a_srcInterval,
                                       const Interval& a_dstInterval)
{
  CH_assert(isDefined());
  CH_assert(!a_coarseFlux.box().isEmpty());

  CH_assert(a_srcInterval.size() == a_dstInterval.size());
  CH_assert(a_srcInterval.begin() >= 0);
  CH_assert(a_srcInterval.end() < a_coarseFlux.nComp());
  CH_assert(a_dstInterval.begin() >= 0);
  CH_assert(a_dstInterval.end() < m_nComp);

  // get edge-centering of coarseFlux
  const Box& edgeBox = a_coarseFlux.box();
  int edgeDir = -1;
  for (int dir=0; dir<SpaceDim; dir++)
    {
      if (edgeBox.type(dir) == IndexType::CELL)
        {
          if (edgeDir == -1)
            {
              edgeDir = dir;
            }
          else
            {
              // already found a cell-centered direction (should only be
              // one for edge-centering)
              MayDay::Error("LevelFluxRegisterEdge::incrementCoarse -- e-field not edge-centered");
            }
        }
    } // end loop over directions
  CH_assert(edgeDir != -1);

  FArrayBox& thisCrseReg = m_regCoarse[a_coarseDataIndex][edgeDir];

  thisCrseReg.plus(a_coarseFlux, -a_scale, a_srcInterval.begin(),
                   a_dstInterval.begin(), a_srcInterval.size());

}

// ---------------------------------------------------------
// 27 March 2003:
// This is called by other functions, and should not be called directly.
Real
norm(const BoxLayoutData<NodeFArrayBox>& a_layout,
     const Real a_dx,
     const int a_p,
     const Interval& a_interval,
     bool a_verbose)
{
  if (a_p == 0)
    return maxnorm(a_layout, a_interval, a_verbose);

  Real normTotal = 0.;

  for (DataIterator it = a_layout.dataIterator(); it.ok(); ++it)
    {
      const NodeFArrayBox& thisNfab = a_layout[it()];
      const Box& thisBox(a_layout.box(it())); // CELL-centered
      Real thisNfabNorm = norm(thisNfab, a_dx, thisBox, a_p,
                               a_interval.begin(), a_interval.size());
      if (a_verbose)
        cout << a_p << "norm(" << thisBox << ") = " << thisNfabNorm << endl;

      if (a_p == 1)
        {
          normTotal += thisNfabNorm;
        }
      else if (a_p == 2)
        {
          normTotal += thisNfabNorm * thisNfabNorm;
        }
      else
        {
          normTotal += pow(thisNfabNorm, Real(a_p));
        }
    }
# ifdef CH_MPI
  Real recv;
  // add up (a_p is not 0)
  int result = MPI_Allreduce(&normTotal, &recv, 1, MPI_CH_REAL,
                             MPI_SUM, Chombo_MPI::comm);
  if (result != MPI_SUCCESS)
    { //bark!!!
      MayDay::Error("sorry, but I had a communication error on norm");
    }
  normTotal = recv;
# endif
  // now do sqrt, etc
  if (a_p == 2)
    normTotal = sqrt(normTotal);
  else
    if ((a_p != 0) && (a_p != 1))
      normTotal = pow(normTotal, (Real)1.0/Real(a_p));

  return normTotal;
}

void
EBFineToCoarRedist::
redistribute(LevelData<EBCellFAB>& a_coarSolution,
             const Interval& a_variables)
{
  CH_TIME("EBFineToCoarRedist::redistribute");
  Real nrefD = 1.0;
  for (int idir = 0; idir < SpaceDim; idir++)
    nrefD *= m_refRat;
  //copy the buffer to the coarse layout
  m_regsFine.copyTo(a_variables, m_regsRefCoar, a_variables);
  //redistribute the refined coarse registers to the coarse solution
  int ibox = 0;
  for (DataIterator dit = m_gridsCoar.dataIterator(); dit.ok(); ++dit)
    {
      const BaseIVFAB<Real>& regRefCoar = m_regsRefCoar[dit()];
      const IntVectSet& ivsRefCoar = m_setsRefCoar[dit()];
      const EBISBox& ebisBoxRefCoar = m_ebislRefCoar[dit()];
      const EBISBox& ebisBoxCoar = m_ebislCoar[dit()];
      const BaseIVFAB<VoFStencil>& stenFAB = m_stenRefCoar[dit()];

      EBCellFAB& solFAB = a_coarSolution[dit()];

      for (VoFIterator vofit(ivsRefCoar, ebisBoxRefCoar.getEBGraph());
          vofit.ok(); ++vofit)
        {
          const VolIndex& srcVoFFine = vofit();
          const VoFStencil& vofsten = stenFAB(srcVoFFine, 0);
          for (int isten = 0; isten < vofsten.size(); isten++)
            {
              const Real& weight = vofsten.weight(isten);
              const VolIndex& dstVoFFine = vofsten.vof(isten);
              VolIndex dstVoFCoar =
                m_ebislRefCoar.coarsen(dstVoFFine,m_refRat, dit());
              //by construction...
              CH_assert(m_gridsCoar.get(dit()).contains(dstVoFCoar.gridIndex()));
              Real dstVolFracFine = ebisBoxRefCoar.volFrac(dstVoFFine);
              Real dstVolFracCoar = ebisBoxCoar.volFrac(dstVoFCoar);
              Real denom = dstVolFracCoar*nrefD;
              for (int ivar = a_variables.begin();
                  ivar <= a_variables.end();  ivar++)
                {
                  Real dmFine = regRefCoar(srcVoFFine, ivar);
                  //ucoar+= massfine/volcoar, ie.
                  //ucoar+= (wcoar*dmCoar*volFracfine/volfraccoar)=massfine/volcoar
                  Real dUCoar = dmFine*weight*dstVolFracFine/denom;
                  solFAB(dstVoFCoar, ivar) += dUCoar;
                }
            }
        }
      ibox++;
    }
}

void
EBCoarToFineRedist::
redistribute(LevelData<EBCellFAB>& a_fineSolution,
             const Interval& a_variables)
{
  CH_TIME("EBCoarToFineRedist::redistribute");
  //copy the buffer to the fine layout
  m_regsCoar.copyTo(a_variables, m_regsCedFine, a_variables);
  //redistribute the coarsened fine registers to the fine solution
  for (DataIterator dit = m_gridsFine.dataIterator(); dit.ok(); ++dit)
    {
      const BaseIVFAB<Real>& regCoar = m_regsCedFine[dit()];
      const IntVectSet& ivsCoar = m_setsCedFine[dit()];
      const EBISBox& ebisBoxCoar = m_ebislCedFine[dit()];
      const BaseIVFAB<VoFStencil>& stenFAB = m_stenCedFine[dit()];

      EBCellFAB& solFAB = a_fineSolution[dit()];

      for (VoFIterator vofitCoar(ivsCoar, ebisBoxCoar.getEBGraph());
          vofitCoar.ok(); ++vofitCoar)
        {
          const VolIndex& srcVoFCoar = vofitCoar();
          const VoFStencil& vofsten = stenFAB(srcVoFCoar, 0);
          for (int isten = 0; isten < vofsten.size(); isten++)
            {
              const Real& weight = vofsten.weight(isten);
              const VolIndex& dstVoFCoar = vofsten.vof(isten);
              Vector<VolIndex> vofsFine =
                m_ebislCedFine.refine(dstVoFCoar,m_refRat, dit());

              for (int ivar = a_variables.begin();
                  ivar <= a_variables.end();  ivar++)
                {
                  Real dmCoar = regCoar(srcVoFCoar, ivar);
                  for (int ifine = 0; ifine < vofsFine.size(); ifine++)
                    {

                      const VolIndex& dstVoFFine = vofsFine[ifine];
                      //ufine += (wcoar*dmCoar) (piecewise constant density diff)
                      Real dUFine = dmCoar*weight;
                      solFAB(dstVoFFine, ivar) += dUFine;
                    }
                }
            }
        }
    }
}

void
EBCoarseAverage::averageFAB(BaseIVFAB<Real>&       a_coar,
                            const BaseIVFAB<Real>& a_fine,
                            const DataIndex&       a_datInd,
                            const Interval&        a_variables) const
{
  CH_assert(isDefined());
  //recall that datInd is from the fine layout.
  const EBISBox& ebisBoxCoar = m_eblgCoFi.getEBISL()[a_datInd];
  const EBISBox& ebisBoxFine = m_eblgFine.getEBISL()[a_datInd];
  const IntVectSet& coarIrregIVS = a_coar.getIVS();
  const IntVectSet& fineIrregIVS = a_fine.getIVS();

  for (VoFIterator vofitCoar(coarIrregIVS, ebisBoxCoar.getEBGraph());
      vofitCoar.ok(); ++vofitCoar)
    {
      const VolIndex& coarVoF = vofitCoar();
      Vector<VolIndex> fineVoFs =
        m_eblgCoFi.getEBISL().refine(coarVoF, m_refRat, a_datInd);

      for (int ivar = a_variables.begin(); ivar <= a_variables.end(); ivar++)
        {
          int  numVoFs = 0;
          Real areaTot = 0;
          Real dataVal = 0;
          for (int ifine = 0; ifine < fineVoFs.size(); ifine++)
            {
              const VolIndex& fineVoF = fineVoFs[ifine];
              if (fineIrregIVS.contains(fineVoF.gridIndex()))
                {
                  Real bndryArea = ebisBoxFine.bndryArea(fineVoF);
                  if (bndryArea > 0)
                    {
                      areaTot += bndryArea;
                      numVoFs++;
                      dataVal += a_fine(fineVoF, ivar);
                    }
                }
            }
          if (numVoFs > 1)
            {
              dataVal /= Real(numVoFs);
            }
          a_coar(coarVoF, ivar) = dataVal;
        }
    }
}

void BlockTab::findIntervals( void )
{
	vector< bool > visited( curBlockNum + 1, false );
	list< BasicBlock * > left;
	left.push_back( blockList[ 1 ] );        //entry node
	visited[ blockList[ 1 ] -> no ] = true;
	
	while ( !left.empty( ) ) { 
		Interval current;	
		current.Head( *left.begin( ) );
		bool added = true;
		while ( added == true ) {
		added = false;
		for ( list< BasicBlock * >::iterator i = left.begin( );
							i != left.end( ); ) {
//includes:check if first set includes the second set 

			if ( current.Head( ) == *i || 
				includes( current.begin( ), current.end( ), 
				( *i ) -> predecessors.begin( ), 
						(*i) -> predecessors.end())) {

				current.insert(*i);
				(*i) -> head = current.Head();
               			if ((*i) -> getTakenPtr() && 
					!visited[(*i) -> getTakenPtr() -> no]) {
					
				    visited[(*i) -> getTakenPtr() -> no] = true;
					left.push_back((*i) -> getTakenPtr());
					added = true;
				}
	       			if ((*i) -> getNTakenPtr() && 
					!visited[(*i)-> getNTakenPtr() -> no]) {
				   visited[(*i) -> getNTakenPtr() -> no] = true;
					left.push_back((*i) -> getNTakenPtr());
					added = true;
				}
				left.erase(i++); 
	       		}
			else
				i++;
		}
		}
		all.push_back(current);
	}
}