flashhenrylaw.cpp

Go to the documentation of this file.

#include "flashhenrylaw.h"
#include "hdf5orthowriter.h"

FlashHenryLaw::FlashHenryLaw(OrthoMesh &mesh,double Kh,double vn,double R,double T,double mc,double mw,double vc,double vw) 
   :FlashCompositional(2,2),
    m_invKh(1.0/Kh),
    m_c(-vn/(R*T)),
    m_vc(vc),
    m_vw(vw),
    m_co2MV(mc),
    m_H2OMV(mw),
    m_mesh(mesh)
{
  m_components_molar_mass.reinit(2);
  m_components_molar_mass(WATER)=18.01528e-3;
  m_components_molar_mass(CO2)=44.01e-3;

  //set the molar volumes
  //m_co2MV=6.1957868649318467E-5; //m3/moles
  //m_H2OMV=1.7157409523809526E-5;
}



 void FlashHenryLaw::getPhasesViscosities(double P,const FlashData &data,VecDouble &visc) 
{
  assert(visc.size() == 2);
  visc(0)=m_vw;
  visc(1)=m_vc;
}

 const VecDouble& FlashHenryLaw::getComponentsMolarMass() const 
{
  return m_components_molar_mass;
}


 void FlashHenryLaw::getPhasesVolume(double P, const FlashData &data, VecDouble &phasesVol) 
{
  assert(data.getNPhases() == 2);
  assert(data.getNComponents() ==2);
  assert(phasesVol.size() ==2);

  phasesVol(1)=(data.getMoles(AQUEOUS,CO2) + data.getMoles(CO2_RICH,CO2))*m_co2MV;
  phasesVol(0)=data.getMoles(AQUEOUS,WATER)*m_H2OMV;

  
  
}

 void FlashHenryLaw::flash(double P, const VecDouble &compTotalMoles, FlashData &data) 
{
  //Calculate the maximum allowed mole fraction of dissolved co2
  double XD = P*m_invKh*exp(m_c*P);

  //Calculate the maximum co2 moles that can dissolve.
  double mdCO2=XD*compTotalMoles(WATER)/(1-XD);

  //Get the total moles of co2;
  double totalCO2=compTotalMoles(CO2);

  if (totalCO2>mdCO2)
  {
    //water always in the aqueous phase
    data.setMoles(AQUEOUS,WATER,compTotalMoles(WATER));
    data.setMoles(CO2_RICH,WATER,0);

    //Maximum dissolved co2 goes to water
    data.setMoles(AQUEOUS,CO2,mdCO2);
    //the remaining co2 stays in co2 rich phase
    data.setMoles(CO2_RICH,CO2,totalCO2-mdCO2);
    
  }
  else
  {
    //water always in the aqueous phase
    data.setMoles(AQUEOUS,WATER,compTotalMoles(WATER));
    data.setMoles(CO2_RICH,WATER,0);

    //All the CO2 goes into water
    data.setMoles(AQUEOUS,CO2,totalCO2);
    //nothing to be in the co2  rich phase
    data.setMoles(CO2_RICH,CO2,0);
  }
}


 FlashHenryLaw::~FlashHenryLaw() 
{

}


 void FlashHenryLaw::printOutput() 
{
  HDF5OrthoWriter &hdf5 = HDF5OrthoWriter::getHDF5OrthoWriter(); 
  unsigned nCells = m_mesh.numCells();
  VecDouble vc(nCells);
  VecDouble vcSatC(nCells);
  FlashData data(numPhases(),numComponents(),NULL);
  VecDouble vv(m_mesh.numVertices());
  VecDouble phasesVol(2);
  for (unsigned i=0;i<nCells;i++)
  {
    FlashCompositional::flash(i,data);
    getPhasesVolume(getDynamicModule().getPressureAtCells()(i),data,phasesVol);
    vcSatC(i)=phasesVol(CO2_RICH)/(phasesVol(AQUEOUS)+phasesVol(CO2_RICH));
  }
  m_mesh.projectCentralValuesAtVertices(vcSatC,vv);
  hdf5.writeScalarField(vv,"SatC");

}

---