LaxFriedrichsForSystem Class Reference

#include <laxfriedrichsforsystem.h>

Inheritance diagram for LaxFriedrichsForSystem:
Inheritance graph
[legend]
Collaboration diagram for LaxFriedrichsForSystem:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 LaxFriedrichsForSystem (OrthoMesh &mesh, Function3D &fInitU, const VecDouble &cPor, Function3D &fPrescribedU, FaceFluxFunction &flux, FixedValueCondition &fixedC, double CFL)
virtual ArrayOfVecDoublegetSolutionAtCells ()
virtual void updateVelocities (DynamicBase &dynMod)
const VecDoublegetPorosity ()
virtual void iterateN (unsigned nSteps, double dt)
virtual double getDt (double t, double tEnd)
virtual ~LaxFriedrichsForSystem ()
virtual void printOutput ()
virtual void updateFixedPressureBC (Function3D *fDP, FlashCompositional &flash)
void getAdjCellOfBCFaces (OrthoMesh &mesh, VecIncBC &vbc, VecIndex &vCells)
void updateBCForCompressibility (VecIncBC &vbc, DynamicBase &dyn, FlashCompositional &flash)

Protected Attributes

VecIncBC m_vbc
VecIndex m_cellsOfBCFaces
OrthoMeshm_mesh
const VecDoublem_cPor
FaceFluxFunctionm_flux
const unsigned n_components
ArrayOfVecDouble m_sol
ArrayOfVecDouble m_solAnt
FixedValueCondition m_fixedC
double m_CFL

Detailed Description

Definition at line 14 of file laxfriedrichsforsystem.h.


Constructor & Destructor Documentation

LaxFriedrichsForSystem::LaxFriedrichsForSystem ( OrthoMesh mesh,
Function3D fInitU,
const VecDouble cPor,
Function3D fPrescribedU,
FaceFluxFunction flux,
FixedValueCondition fixedC,
double  CFL 
)

Definition at line 8 of file laxfriedrichsforsystem.cpp.

00009   :m_mesh(mesh),m_cPor(cPor),m_flux(flux),n_components(flux.numComponents()),m_fixedC(fixedC)
00010     
00011 
00012  {
00013    assert(fInitU.n_components() == n_components);
00014    assert(fPrescribedU.n_components() == n_components);
00015   assert(cPor.size() == mesh.numCells());
00016   m_CFL = CFL;
00017   
00018 
00019   //Alloc Vectors
00020   m_sol.reinit(mesh.numCells(),n_components);
00021   m_solAnt.reinit(mesh.numCells(),n_components);
00022 
00023   
00024   
00025   mesh.setCentralValuesFromFunction(fInitU,m_sol);
00026     
00027 
00028   //Initialize the boundary condtions
00029   m_vbc.set_bc_from_function(mesh,fPrescribedU);
00030 
00031   //  std::cout << m_vbc << std::endl;
00032   //exit(0);
00033   //Initial conditions
00034   //m_vSol.reinit(mesh.numCells(),n_components);
00035   //m_vSolAnt.reinit(mesh.numCells(),n_components);
00036   //mesh.setCentralValuesFromFunction(fInitU,m_vSol);
00037   
00038  }

virtual LaxFriedrichsForSystem::~LaxFriedrichsForSystem (  )  [inline, virtual]

Definition at line 42 of file laxfriedrichsforsystem.h.

00042 {}


Member Function Documentation

void LaxFriedrichsForSystem::getAdjCellOfBCFaces ( OrthoMesh mesh,
VecIncBC vbc,
VecIndex vCells 
)
double LaxFriedrichsForSystem::getDt ( double  t,
double  tEnd 
) [virtual]

Implements TransportBase.

Reimplemented in LaxFriedrichsSystemMPI.

Definition at line 65 of file laxfriedrichsforsystem.cpp.

00066 {
00067   assert(tEnd >= t);
00068   double timeInterval = tEnd-t;
00069   double dt = this->calculateTimeStep(m_mesh,timeInterval,m_CFL,getPorosity(),m_flux);
00070   return dt;
00071 }

const VecDouble& LaxFriedrichsForSystem::getPorosity (  )  [inline]

Definition at line 39 of file laxfriedrichsforsystem.h.

00039 {return m_cPor;} 

ArrayOfVecDouble & LaxFriedrichsForSystem::getSolutionAtCells (  )  [virtual]

Implements TransportBase.

Definition at line 42 of file laxfriedrichsforsystem.cpp.

00043 {
00044   return m_sol;
00045 }

void LaxFriedrichsForSystem::iterateN ( unsigned  nSteps,
double  dt 
) [virtual]

Implements TransportBase.

Reimplemented in GodunovMethod, LaxFriedrichsSystemMPI, RussanovForSystem, RussanovSystemMPI, and UpwindForSystem.

Definition at line 74 of file laxfriedrichsforsystem.cpp.

00075 {
00076   unsigned posCell,negCell;
00077 
00078   //We iterate the method along faces, so 
00079   //alloc the vectors to contain the solutions in neg
00080   //and pos cells of the faces.  Remeber that the solution
00081   //in a cell is vector with a value for each component
00082   static VecDouble  vBC(n_components),vFlux(n_components);
00083   static VecDoubleRef vQPos(n_components),vQNeg(n_components);
00084   //For each time step
00085   for (unsigned count = 0; count < nSteps; count++)
00086   {
00087     //Get the face iterator and the number of faces.
00088     OrthoMesh::Cash_Face_It face = m_mesh.begin_cash_face();
00089     OrthoMesh::Cash_Face_It endf = m_mesh.end_cash_face();
00090     m_vbc.rewind();
00091     bool hasbc;
00092 
00093     m_solAnt=m_sol;
00094     
00095     //For each face
00096     for (;face!=endf;face++)
00097     {
00098       face->getAdjCellIndices(negCell,posCell);
00099       //Get the values of the previous solution in the right and left of the face
00100       hasbc=this->getValuesOfTheFaceCells(m_mesh,m_vbc,m_solAnt,face->index(),negCell,posCell,vQNeg,vQPos);
00101 
00102       //printf("face %d : Values %g %g %g | %g %g %g\n",face->index(),vQNeg(0),vQNeg(1),vQNeg(2),vQPos(0),vQPos(1),vQPos(2));
00103       
00104       //Get the indices of the cells that contain the face.
00105       //If the face is at boundary, it has just one cell.
00106       //In this case the method getAdjCells() return posCell == negCell
00107       if (posCell == OrthoMesh::INVALID_INDEX)
00108         posCell=negCell;
00109       else if (negCell == OrthoMesh::INVALID_INDEX)
00110         negCell=posCell;
00111         
00112 
00113 
00114       double posPor = getPorosity()(posCell);
00115       double negPor = getPorosity()(negCell);
00116       double mPor = (posPor + negPor)/2.0;
00117 
00118 
00119       //get Flux
00120       m_flux.fluxAtFace(vFlux,face,vQNeg,vQPos);
00121       
00122       //Now update the cell values component by component
00123       for (unsigned cmp=0;cmp<n_components;cmp++)
00124       {
00125         
00126         //See if the face has prescribed boundary condition BC      
00127         double flux;
00128         if (!hasbc)    
00129         {
00130           /*
00131             The face doesnt have boundary condition BC
00132             proceed to the usual computation of the flux
00133             to calculate the amount of the mass passed through
00134             the face and divide it by the cell volume
00135           */
00136            
00137           flux = dt*face->areaPerCellVol() * vFlux(cmp) - mPor*(vQPos(cmp)-vQNeg(cmp))/2.0;
00138         }
00139         else //the face has BC, so dont add the stability term
00140         {
00141           flux = dt*face->areaPerCellVol() * vFlux(cmp);
00142         }
00143 
00144       
00145       
00146         //printf("face %d = %g,Pos: %d,Neg %d, SwPos = %g,SwNeg = %g \n",faceIndex,flux,face->getPosCellIndex(),face->getNegCellIndex(),SwPos,SwNeg);
00147         if (face->hasPosCell())
00148           m_sol(posCell,cmp) +=  + flux/posPor;
00149 
00150         if (face->hasNegCell())
00151           m_sol(negCell,cmp) +=  - flux/negPor;
00152         
00153       } //end for each component
00154     } //end for each face
00155   } //end for each time step
00156 }

void LaxFriedrichsForSystem::printOutput (  )  [virtual]

Implements TransportBase.

Definition at line 163 of file laxfriedrichsforsystem.cpp.

00164 {
00165   VecDouble vValues(m_mesh.numVertices());
00166   VecDouble cValues(m_mesh.numCells());
00167 
00168   double acum=0.0;
00169   for (unsigned cmp=0;cmp<n_components;cmp++)
00170   {
00171     char str[100];
00172     m_sol.getComponent(cValues,cmp);
00173     this->m_mesh.projectCentralValuesAtVertices(cValues,vValues);
00174     //printf("Mass S%d: %g\n",cmp+1,m_mesh.getIntegralAtCells(cValues);
00175     acum+=m_mesh.getIntegralAtCells(cValues);
00176     HDF5OrthoWriter &hdf5 = HDF5OrthoWriter::getHDF5OrthoWriter(); 
00177     sprintf(str,"S%d",cmp+1);
00178     hdf5.writeScalarField(vValues,str);
00179 
00180   }
00181 }

void LaxFriedrichsForSystem::updateBCForCompressibility ( VecIncBC vbc,
DynamicBase dyn,
FlashCompositional flash 
)

Definition at line 234 of file laxfriedrichsforsystem.cpp.

00235 {
00236   const VecDouble &P = dyn.getPressureAtCells();
00237   VecIncBC::data_iterator itBC = vbc.begin_data();
00238   VecDouble phasesVol(flash.numPhases());
00239   FlashData flashData(flash.numPhases(),flash.numComponents());
00240   VecDoubleRef m(vbc.bc_dim());
00241   for (;itBC!= vbc.end_data();itBC++)
00242   {
00243     //get the bc
00244     m.setData(itBC->data.p);
00245     //appy the flash
00246     flash.flash(P(itBC->data.nearCell),m,flashData);
00247     //calculate the volumes
00248     flash.getPhasesVolume(P(itBC->data.nearCell),flashData,phasesVol);
00249     double factor=1.0/NumericMethods::vectorSum(phasesVol);
00250     //Multiply the components by such factor
00251     m*=factor;
00252   }
00253 }

void LaxFriedrichsForSystem::updateFixedPressureBC ( Function3D fDP,
FlashCompositional flash 
) [virtual]

Implements ConservativeMethodForSystem.

Definition at line 186 of file laxfriedrichsforsystem.cpp.

00187 {
00188   FlashData flashData(flash.numPhases(),flash.numComponents());
00189   VecDouble phasesVol(flash.numPhases());
00190   assert(n_components == flash.numComponents());
00191   OrthoMesh::Face_It fend = m_mesh.end_face();
00192   flashData.allocateOwnMemory();
00193   Point3D pt;
00194   VecDoubleRef bc(n_components);
00195   m_vbc.rewind();
00196   //for all faces
00197   //m_fBC[0].print(std::cout);
00198   for (OrthoMesh::Face_It face=m_mesh.begin_face(); face!=fend; face++)
00199     {
00200       //If face is at boundary and there is boundary conditions for the transport
00201       unsigned i = face->index();
00202       if (m_vbc.get_ref_bc(i,bc))
00203       {
00204         face->barycenter(pt);
00205         if (fDP->isInDomain(pt))
00206         {
00207           //Get components for face i
00208 
00209           double P = (*fDP)(pt);
00210 
00211           //Make the flash and get the phase volumes
00212           flash.flash(P,bc,flashData);
00213           flash.getPhasesVolume(P,flashData,phasesVol);
00214 
00215           //calculate factor
00216           double factor=1.0/NumericMethods::vectorSum(phasesVol);
00217           //printf("%g %g %g\n",factor,phasesVol(0),phasesVol(1));
00218           //flashData.print();
00219           bc*=factor;
00220         }
00221       }
00222     }
00223   std::cout << m_vbc;
00224 }

void LaxFriedrichsForSystem::updateVelocities ( DynamicBase dynMod  )  [virtual]

Update the velocities. The method pass the dynamic module to the FaceFluxFunction class so that the flux class can retrieve the necessary information about the velocity fields.

Parameters:
dynMod 

Implements TransportBase.

Definition at line 53 of file laxfriedrichsforsystem.cpp.

00054 {
00055   m_flux.updateDynamicData(dynMod);
00056   FlashCompositional *flash = getFlash();
00057 
00058   if (flash)
00059     this->updateBCForCompressibility(m_vbc,dynMod,*flash);
00060   
00061 }


Member Data Documentation

Definition at line 18 of file laxfriedrichsforsystem.h.

double LaxFriedrichsForSystem::m_CFL [protected]

Definition at line 29 of file laxfriedrichsforsystem.h.

Vector containg the porosities in each cell

Definition at line 20 of file laxfriedrichsforsystem.h.

Definition at line 24 of file laxfriedrichsforsystem.h.

Definition at line 21 of file laxfriedrichsforsystem.h.

The mesh

Definition at line 19 of file laxfriedrichsforsystem.h.

Definition at line 23 of file laxfriedrichsforsystem.h.

Definition at line 23 of file laxfriedrichsforsystem.h.

Definition at line 17 of file laxfriedrichsforsystem.h.

const unsigned LaxFriedrichsForSystem::n_components [protected]

Definition at line 22 of file laxfriedrichsforsystem.h.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Sun Apr 8 23:13:14 2012 for CO2INJECTION by  doxygen 1.6.3