russanovmpioverlap.cpp

Go to the documentation of this file.

#include "russanovmpioverlap.h"
 RussanovMPIOverlap::RussanovMPIOverlap(OrthoMesh &mesh,Function3D &fInitU,const VecDouble &cPor,Function3D &fPrescribedU,FaceFluxFunction &flux, double CFL,int mesh_overlap) 
   :LaxFriedrichsMPIOverlap(mesh,fInitU,cPor,fPrescribedU,flux, CFL, mesh_overlap)
{

}




void RussanovMPIOverlap::iterate(double t, double tEnd) 
{
  
  double timeInterval = tEnd-t;
  double localDt = this->calculateTimeStep(m_mesh,timeInterval,m_CFL,getPorosity(),m_flux);
  double dt=negotiateTimeStep(localDt);
  printf("Russ Negotiating dt: Sent %g, got %g\n",localDt,dt); 
  

  int nIteraction = (int) round(timeInterval/dt); //number of iteractions


  //For each time step
  for (int count = 0; count < nIteraction; count++)
  {
    //Get the face iterator and the number of faces.
    OrthoMesh::Face_It face = m_mesh.begin_face();
    unsigned nFaces = m_mesh.numFaces();
    
    m_cValuesPrev=m_cValues;

    //For each face
    VecTag::iterator itDZ = fDeadZone.begin();
    for (unsigned faceIndex=0;faceIndex < nFaces; faceIndex++,face++,itDZ++)
    {
      if (*itDZ != ACTIVE_ZONE)
        continue;
      //Get the values of the previous solution in the right and left of the face
      double SwPos,SwNeg;
      unsigned cell1,cell2;
      this->getValuesOfTheFaceCells(m_mesh,m_fValues,m_cValuesPrev,face,SwNeg,SwPos);

      //Get the values of the cells that contain the face.
      //If the face is at boundary, it has just one cell.
      //In this case the method getAdjCells() return cell1 == cell2
      face->getValidCellIndices(cell1,cell2);


      
      double posPor = getPorosity()(cell1);
      double negPor = getPorosity()(cell2);
      
      //calculate the right coefficient of diffusion 
      double a;
      if ( SwPos >= SwNeg)
        a = fabs(m_flux.maxLocalCharVelocity(face,faceIndex,cell1,cell2,SwNeg,SwPos));
      else
        a = fabs(m_flux.maxLocalCharVelocity(face,faceIndex,cell1,cell2,SwPos,SwNeg));

      //printf("%g ",a);
       //Calculate the amount of the mass passed through the face and divide it by the cell volume
      double flux = dt*face->areaPerCellVol()* ( m_flux.fluxAtFace(face,faceIndex,cell1,cell2,SwPos,SwNeg) - a*(SwNeg-SwPos)/2.0);
      //printf("face %d = %g,Pos: %d,Neg %d, SwPos = %g,SwNeg = %g \n",faceIndex,flux,face->getPosCellIndex(),face->getNegCellIndex(),SwPos,SwNeg);
      if (face->hasPosCell())
        m_cValues(face->getPosCell()) +=  - flux/posPor;

      if (face->hasNegCell())
        m_cValues(face->getNegCell()) +=  + flux/negPor;
    }
    updateGhostCells();
  }
  updateDeadZone();
}



RussanovMPIOverlap::~RussanovMPIOverlap() 
{

}

---