CO2INJECTION: facefluxsimpleblackoilmass.cpp Source File

00001 #include "facefluxsimpleblackoilmass.h"
00002 #include "numericmethods.h"
00003 
00004 
00005  FaceFluxSimpleBlackOilMass::~FaceFluxSimpleBlackOilMass() 
00006 {
00007 
00008 }
00009 
00010  FaceFluxSimpleBlackOilMass::FaceFluxSimpleBlackOilMass(OrthoMesh &mesh,FlashCompositional &flash) 
00011    :FaceFluxSimpleBlackOil(mesh,flash),m_vNormalMassFlux(mesh.numFaces())
00012 {
00013 }
00014 
00015  void FaceFluxSimpleBlackOilMass::updateDynamicData(DynamicBase &dynMod) 
00016 {
00017   FaceFluxSimpleBlackOil::updateDynamicData(dynMod);
00018   dynMod.getNormalMassFluxAtFaces(m_vNormalMassFlux);
00019   //  m_vNormalMassFlux.print(std::cout,8);
00020 }
00021 
00022 using NumericMethods::_div;
00023 
00024  void FaceFluxSimpleBlackOilMass::fluxAtFace(VecDouble &vFlux, OrthoMesh::Face_It &faceIt,int face,int cell1,int cell2,const VecDouble &Q1,const VecDouble &Q2) 
00025 {
00026   assert(vFlux.size()==2);
00027   assert((unsigned) face < m_vNormalVel.size());
00028   static FlashData data1(2,2),data2(2,2);
00029   static VecDouble phasesVolume1(2);
00030   static VecDouble phasesVolume2(2);
00031   static VecDouble vFlux1(2);
00032   static VecDouble vFlux2(2);
00033   const VecDouble &P = getPressure();
00034   m_bUnsaturated=false;
00035 
00036   m_flash.flash(P(cell1),Q1,data1);
00037   m_flash.flash(P(cell2),Q2,data2);
00038 
00039   m_flash.getPhasesViscosities(P(cell1),data1,m_vViscosities1);
00040   m_flash.getPhasesViscosities(P(cell2),data2,m_vViscosities2);
00041   
00042   m_flash.getPhasesVolume(P(cell1),data1,phasesVolume1);
00043   m_Sw1=phasesVolume1(0);
00044 
00045   
00046   m_mob.setParameters(m_vViscosities1);
00047   double mobO=m_mob(m_Sw1);
00048 
00049   vFlux(0)=data1.getMoles(0,0)*mobO/(data1.getMoles(0,0)*mobO + data1.getMoles(0,1)*mobO + _div(data1.getMoles(1,1),phasesVolume1(1))*(1-mobO)*phasesVolume1(0));
00050   vFlux(1)=1-vFlux(0);
00051   
00052 
00053 
00054   m_flash.getPhasesVolume(P(cell2),data2,phasesVolume2);
00055   m_Sw2=phasesVolume2(0);
00056   mobO=m_mob(m_Sw2);
00057   m_mob.setParameters(m_vViscosities2);
00058 
00059 
00060   vFlux2(0)=data2.getMoles(0,0)*mobO/(data2.getMoles(0,0)*mobO + data2.getMoles(0,1)*mobO + _div(data2.getMoles(1,1),phasesVolume2(1))*(1-mobO)*phasesVolume2(0));
00061   vFlux2(1)=1-vFlux2(0);
00062 
00063   
00064   vFlux+=vFlux2;
00065   vFlux*= m_vNormalMassFlux(face)/2.0;
00066 
00067   if (data1.getMoles(1,1) == 0.0 || data2.getMoles(1,1)==0.0)
00068   {
00069     //printf("Unsaturated at %g\n",faceIt->barycenter()[0]);
00070     m_bUnsaturated=true;
00071   }
00072 
00073 
00074 //Store the cells just to check at getLocalCharVelocity
00075   m_cell1=cell1;
00076   m_cell2=cell2;
00077 
00078   /*
00079   VecDouble vF(2);
00080   FaceFluxSimpleBlackOil::fluxAtFace(vF, faceIt, face, cell1, cell2,Q1,Q2);
00081   printf("vFlux(0)=%g %g\n",vFlux(0),vF(0));
00082   ********/
00083 }
00084 
00085 
00086  void FaceFluxSimpleBlackOilMass::maxGlobalCharVelocity(double vel[3]) 
00087 {
00088   double maxVel=NumericMethods::vectorMaxAbsValue(m_vNormalMassFlux);
00089   double maxVel2=NumericMethods::vectorMaxAbsValue(m_vNormalVel);
00090   printf("Max Vel: %g %g",maxVel,maxVel2);
00091   double P=NumericMethods::vectorMinAbsValue(getPressure());
00092   static FlashData data(2,2);
00093   data.allocateOwnMemory();
00094   double fact=pow((P/10000 + 1.0)/(0.05*P+ 1),2.0)*(0.06*P + 1);
00095 
00096   m_flash.getPhasesViscosities(P,data,m_vViscosities1);
00097   m_Dmob.setParameters(m_vViscosities1);
00098   m_Dmob.updateInflectionPoint();
00099   vel[X]= fact*maxVel*m_Dmob.getMaxNorm(0,1);
00100   vel[Y]=0;
00101   vel[Z]=0;
00102   return;
00103 
00104 }
00105