---

BaseIFFABI.H

Go to the documentation of this file.

/* _______              __
  / ___/ /  ___  __ _  / /  ___
 / /__/ _ \/ _ \/  ' \/ _ \/ _ \
 \___/_//_/\___/_/_/_/_.__/\___/ 
*/
//
// This software is copyright (C) by the Lawrence Berkeley
// National Laboratory.  Permission is granted to reproduce
// this software for non-commercial purposes provided that
// this notice is left intact.
// 
// It is acknowledged that the U.S. Government has rights to
// this software under Contract DE-AC03-765F00098 between
// the U.S.  Department of Energy and the University of
// California.
//
// This software is provided as a professional and academic
// contribution for joint exchange. Thus it is experimental,
// is provided ``as is'', with no warranties of any kind
// whatsoever, no support, no promise of updates, or printed
// documentation. By using this software, you acknowledge
// that the Lawrence Berkeley National Laboratory and
// Regents of the University of California shall have no
// liability with respect to the infringement of other
// copyrights by any part of this software.
//

#ifndef _BaseIFFABI_H_
#define _BaseIFFABI_H_
#include "MayDay.H"
#include "FaceIterator.H"
/******************/
template <class T> inline 
const EBGraph& 
BaseIFFAB<T>::getEBGraph() const
{
  return m_ebgraph;
}
/******************/
template <class T> inline 
BaseIFFAB<T>::BaseIFFAB()
{
  setDefaultValues();
}
/******************/
template <class T> inline 
BaseIFFAB<T>::~BaseIFFAB()
{
  clear();
}
/******************/
template <class T> inline 
BaseIFFAB<T>::BaseIFFAB(const IntVectSet& a_ivsin, 
                        const EBGraph& a_ebgraph,
                        const int& a_direction,
                        const int& a_nvarin)
{
  setDefaultValues();
  define(a_ivsin, a_ebgraph, a_direction, a_nvarin);
}
/******************/
template <class T> inline 
void
BaseIFFAB<T>::define(const IntVectSet& a_ivsin, 
                     const EBGraph& a_ebgraph,
                     const int& a_direction,
                     const int& a_nvarin)
{
  clear();
  m_isDefined = true;
  assert(a_nvarin > 0);
  assert((a_direction >= 0) && (a_direction < SpaceDim));
  const Box& domain = a_ebgraph.getDomain();
  m_direction = a_direction;
  m_ivs = a_ivsin;
  m_ebgraph = a_ebgraph;
  m_nComp = a_nvarin;
  if(!a_ivsin.isEmpty())
    {
      Box minbox = a_ivsin.minBox();
      m_ivmap.resize(surroundingNodes(minbox,a_direction), 1);
    }
  m_nFaces = 0;
  for(IVSIterator ivsit(a_ivsin); ivsit.ok(); ++ivsit)
    {
      //iterate through sides 
      //so we get the surrounding nodes
      //of the irregular set.
      const IntVect thisIV = ivsit();
      Vector<VolIndex> theseVoFs = m_ebgraph.getVoFs(thisIV);
      for(SideIterator sit; sit.ok(); ++sit)
        {
          int iside = sign(sit());
          IntVect otherIV = thisIV + iside*BASISV(m_direction);
          //need combinatoric number of faces because we are using
          //lexographic ordering to index into the vector.
          //on boundary faces this means that only one vof has a number
          //of vofs.
          int numfaces = 0;
          if(domain.contains(otherIV))
            {
              numfaces = (m_ebgraph.numVoFs(thisIV))*(m_ebgraph.numVoFs(otherIV));
#ifndef NDEBUG
              Vector<FaceIndex> facesHi = m_ebgraph.getAllFaces(thisIV, m_direction, sit());
              Vector<FaceIndex> facesLo = m_ebgraph.getAllFaces(otherIV, m_direction, flip(sit()));
              if(facesHi.size() != facesLo.size())
                {
                  MayDay::Error("internally inconsistent geometry");
                }
#endif
            }
          else
            {
              numfaces = m_ebgraph.numVoFs(thisIV);
            }
          //by convention, the high vof of the face is where
          // we store the indexing information 
          IntVect ivface;
          if(sit() == Side::Lo)
            ivface = thisIV;
          else
            ivface = otherIV;
          
          Vector<int>& vIndex = m_ivmap(ivface, 0);
          vIndex.resize(numfaces);
          //now loop and find the actual faces
          //and label them.
          for(int ivof = 0; ivof < theseVoFs.size(); ivof++)
            {
              Vector<FaceIndex> sideFaces = 
                m_ebgraph.getFaces(theseVoFs[ivof],m_direction, sit());
              for(int iface = 0; iface < sideFaces.size(); iface++)
                {
                  vIndex[iface] = m_nFaces;
                  m_nFaces++;
                }
            } //end loop over vofs in the cell
        } //end loop over sides
    } //end loop over cells in the irregular set

  if(m_nFaces > 0)
    {
      //set up face low and high vectors
      for (int dir = 0; dir < SpaceDim; ++dir)
        {
          m_loVect[dir] = 1;
          m_hiVect[dir] = 1;
        }
      m_hiVect[0] = m_nFaces;
      m_dataPtr = new T[m_nComp*m_nFaces];
    }
  else
    {
      m_loVect = IntVect::Unit;
      m_hiVect = IntVect::Zero;
      m_dataPtr = NULL;
    }
}
/******************/
template <class T> inline 
const IntVectSet& 
BaseIFFAB<T>::getIVS() const
{
  return m_ivs;
}
/******************/
template <class T> inline 
int 
BaseIFFAB<T>::direction() const
{
  return m_direction;
}
/******************/
template <class T> inline 
void
BaseIFFAB<T>::setVal(const T& a_value)
{
  assert(isDefined());
  for(int ivec = 0; ivec < m_nFaces*m_nComp; ivec++)
    m_dataPtr[ivec] = a_value;
}
/******************/
template <class T> inline 
void
BaseIFFAB<T>::copy(const Box& a_fromBox,       
                   const Interval& a_dstInterval,
                   const Box& a_toBox,       
                   const BaseIFFAB<T>& a_src, 
                   const Interval& a_srcInterval)
{
  assert(isDefined());
  assert(a_src.isDefined());
  assert(a_srcInterval.size() == a_dstInterval.size());
  assert(a_dstInterval.begin() >= 0);
  assert(a_srcInterval.begin() >= 0);
  assert(a_dstInterval.end()   < m_nComp);
  assert(a_srcInterval.end()   < a_src.m_nComp);
  assert(a_fromBox == a_toBox);

  Box intBox = a_toBox;
  IntVectSet ivsIntersect = m_ivs;
  ivsIntersect &=  a_src.m_ivs;
  ivsIntersect &=  intBox;
  int compSize = a_srcInterval.size();
  FaceStop::WhichFaces stopCrit = FaceStop::SurroundingWithBoundary;
  for(FaceIterator faceit(ivsIntersect, m_ebgraph, m_direction, stopCrit); 
      faceit.ok(); ++faceit)
    {
      const FaceIndex& face = faceit();
      for(int icomp = 0; icomp < compSize; icomp++)
        {
          int isrccomp = a_srcInterval.begin() + icomp;
          int idstcomp = a_dstInterval.begin() + icomp;
          int isrc = a_src.getIndex(face, isrccomp);
          int idst =       getIndex(face, idstcomp);
          m_dataPtr[idst] = a_src.m_dataPtr[isrc];
        } //end loop over components
    }//end loop over faces
} 
/*********/
template <class T> inline 
int
BaseIFFAB<T>::getLocalVecIndex(const FaceIndex& a_face) const
{
  int retval;
  if(!a_face.isBoundary())
    {
      //this checks that both vofs are defined
      assert(a_face.cellIndex(Side::Lo) >= 0);
      assert(a_face.cellIndex(Side::Hi) >= 0);
      int xlen = m_ebgraph.numVoFs(a_face.gridIndex(Side::Lo));
      int loCell = a_face.cellIndex(Side::Lo);
      int hiCell = a_face.cellIndex(Side::Hi);
      retval =  loCell + xlen*hiCell;
    }
  else
    {
      int loCell = a_face.cellIndex(Side::Lo);
      int hiCell = a_face.cellIndex(Side::Hi);
      //one should be -1, the other should be larger
      assert(((loCell == -1)&&(hiCell > -1))||
             ((hiCell == -1)&&(loCell > -1)));
      //return the one that is not -1
      retval = Max(loCell, hiCell);
    }
  return retval;
}
/*********/
template <class T> inline 
int
BaseIFFAB<T>::getIndex(const FaceIndex& a_face, const int& a_comp) const
{
  assert((a_comp >= 0) && (a_comp < m_nComp));

  //the index is stored in the high vof gridindex by convention
  //(surroundingnodes(minbox) kind of forces it too)
  const Vector<int>& vecOffset = m_ivmap(a_face.gridIndex(Side::Hi), 0);
  int ivec   = getLocalVecIndex(a_face);
  //if this fails on a size() thing, then the
  //basefab<vector<int> > was constructed incorrectly
  //(the assert should catch whether the vof is outside the orig ivs)
  int ioffset = vecOffset[ivec];
  assert(ioffset >= 0);
  assert(ioffset < m_nFaces);
  //now add offset from componentnitude
  ioffset += m_nFaces*a_comp;
  return ioffset;
}
/********************/
template <class T> inline 
void
BaseIFFAB<T>::clear()
{
  m_nComp = 0;
  m_nFaces = 0;
  m_direction = -1;
  m_ivs.makeEmpty();
  m_ivmap.clear();
  if(m_dataPtr != NULL)
    {
      delete[] m_dataPtr;
      m_dataPtr = NULL;
    }
  m_isDefined = false;
}
/*************************/
template <class T> inline 
bool 
BaseIFFAB<T>::isDefined() const
{
  return (m_isDefined);
}
/*************************/
template <class T> inline 
int
BaseIFFAB<T>::numFaces() const 
{
  return m_nFaces;
}
/*************************/
template <class T> inline 
int
BaseIFFAB<T>::nComp() const 
{
  return m_nComp;
}
/*************************/
template <class T> inline 
T&
BaseIFFAB<T>::operator() (const FaceIndex& a_ndin,
                          const int& a_comp)
{ 
  assert(isDefined());
  assert(a_comp >= 0);
  assert(a_comp < m_nComp);
  assert((m_ivs.contains(a_ndin.gridIndex(Side::Lo)) || 
          m_ivs.contains(a_ndin.gridIndex(Side::Hi))));
  assert(a_ndin.direction() == m_direction);
  int iloc = getIndex(a_ndin, a_comp);
  return(m_dataPtr[iloc]);
}
/**************************/
template <class T> inline 
const T&
BaseIFFAB<T>::operator() (const FaceIndex& a_ndin,
                          const int& a_comp) const
{ 
  assert(isDefined());
  assert(a_comp >= 0);
  assert(a_comp < m_nComp);
  assert((m_ivs.contains(a_ndin.gridIndex(Side::Lo)) || 
          m_ivs.contains(a_ndin.gridIndex(Side::Hi))));
  assert(a_ndin.direction() == m_direction);
  int iloc = getIndex(a_ndin, a_comp);
  return(m_dataPtr[iloc]);
}
/******************/
template <class T> inline 
const T*
BaseIFFAB<T>::dataPtr(const int& a_comp) const
{
  assert(isDefined());
  assert(a_comp >= 0);
  assert(a_comp < m_nComp);
  return m_dataPtr + a_comp*m_nFaces;
}
/******************/
template <class T> inline 
T* 
BaseIFFAB<T>::dataPtr(const int& a_comp) 
{
  assert(isDefined());
  return m_dataPtr + a_comp*m_nFaces;
}
/******************/
template <class T> inline 
const int*
BaseIFFAB<T>::loVect() const 
{
  return m_loVect.getVect();
}
/******************/
template <class T> inline 
const int*
BaseIFFAB<T>::hiVect() const 
{
  return m_hiVect.getVect();
}
/******************/
template <class T> inline 
void
BaseIFFAB<T>::setDefaultValues()
{
  m_isDefined = false;
  m_dataPtr = NULL;
  m_nFaces = 0;
  m_nComp = 0;
  m_direction = -1;
  m_loVect = IntVect::Unit;
  m_hiVect = IntVect::Zero;
}
/******************/
template <class T> inline 
BaseIFFAB<T>& 
BaseIFFAB<T>::operator= (const BaseIFFAB<T>& a_input)
{
  MayDay::Error("BaseIFFAB operator = not defined");
  return *this;
}
/******************/
template <class T> inline 
BaseIFFAB<T>::BaseIFFAB (const BaseIFFAB<T>& a_input)
{
  MayDay::Error("BaseIFFAB copy constructor not defined");
}
/******************/
template <class T> inline 
int BaseIFFAB<T>::size(const Box& a_region, 
                       const Interval& a_comps) const
{
  assert(isDefined());
  //create set of cells in fab that are also in the input region
  IntVectSet subIVS = m_ivs;
  subIVS &= a_region;

  //account for number of points in the intersection
  int retval = 0;
  FaceStop::WhichFaces stopCrit = FaceStop::SurroundingWithBoundary;
  for(FaceIterator faceit(subIVS, m_ebgraph, m_direction, stopCrit); 
      faceit.ok(); ++faceit)
    {
      retval++;
    }
  //account for number of components
  retval *= a_comps.size();
  //account for the size of the template object
  retval *= sizeof(T);
  return retval;
}
/********************/
template <class T> inline 
void BaseIFFAB<T>::linearOut(void* a_buf, 
                             const Box& a_region, 
                             const Interval& a_comps) const
{
  assert(isDefined());
  //cast buffer into a pointer of the correct type
  //this way incrementing the pointer goes by the correct amount
  T* buffer = (T*)a_buf;
  //create set of cells in fab that are also in the input region
  IntVectSet subIVS = m_ivs;
  subIVS &= a_region;
  const BaseIFFAB<T>& thisFAB = *this;

  FaceStop::WhichFaces stopCrit = FaceStop::SurroundingWithBoundary;
  for(FaceIterator faceit(subIVS, m_ebgraph, m_direction, stopCrit); 
      faceit.ok(); ++faceit)
    {
      const FaceIndex& face = faceit();
      for(int icomp = a_comps.begin(); icomp <= a_comps.end(); icomp++)
        {
          const T& thisVal = thisFAB(face, icomp);
          //put the value in the buffer
          *buffer = thisVal;
          //increment the buffer offset
          ++buffer;
        }
    }
}
/********************/
template <class T> inline 
void BaseIFFAB<T>::linearIn(void* a_buf, const Box& a_region, const Interval& a_comps)
{
  assert(isDefined());
  //cast buffer into a pointer of the correct type
  //this way incrementing the pointer goes by the correct amount
  T* buffer = (T*)a_buf;
  //create set of cells in fab that are also in the input region
  IntVectSet subIVS = m_ivs;
  subIVS &= a_region;
  BaseIFFAB<T>& thisFAB = *this;

  FaceStop::WhichFaces stopCrit = FaceStop::SurroundingWithBoundary;
  for(FaceIterator faceit(subIVS, m_ebgraph, m_direction, stopCrit); 
      faceit.ok(); ++faceit)
    {
      const FaceIndex& face = faceit();
      for(int icomp = a_comps.begin(); icomp <= a_comps.end(); icomp++)
        {
          T& thisVal = thisFAB(face, icomp);
          //put the value in the buffer
          thisVal = *buffer;
          //increment the buffer offset
          ++buffer;
        }
    }
}
/********************/
/******************/

#endif

---