CoreFlows/examples/C/FiveEqsTwoFluid_1DDepressurisation.cxx

   1 #include "FiveEqsTwoFluid.hxx"
   2
   3 using namespace std;
   4
   5 int main(int argc, char** argv)
   6 {
   7         //Preprocessing: mesh and group creation
   8         cout << "Building cartesian mesh" << endl;
   9         double xinf=0.0;
  10         double xsup=4.2;
  11         int nx=50;
  12         Mesh M(xinf,xsup,nx);
  13         double eps=1.E-8;
  14         M.setGroupAtPlan(xsup,0,eps,"Outlet");//Neumann
  15         M.setGroupAtPlan(xinf,0,eps,"Wall");//
  16         int spaceDim = M.getSpaceDimension();
  17
  18         // set the limit field for each boundary
  19         LimitField limitOutlet, limitWall;
  20         map<string, LimitField> boundaryFields;
  21         limitOutlet.bcType=Outlet;
  22         limitOutlet.p = 100e5;
  23         boundaryFields["Outlet"] = limitOutlet;
  24
  25         limitWall.bcType=Wall;
  26         limitWall.T = 600;
  27         limitWall.v_x = vector<double>(2,0);
  28         boundaryFields["Wall"]= limitWall;
  29
  30         // physical constants
  31         double latentHeat=1e6;
  32         double satTemp=618;
  33         double dHsatl_over_dp=0.05;
  34         double Psat=85e5;
  35
  36         FiveEqsTwoFluid  myProblem(around155bars600K,spaceDim);
  37         int nbPhase = myProblem.getNumberOfPhases();
  38         int nVar = myProblem.getNumberOfVariables();
  39
  40         //Initial field creation
  41         Vector VV_Constant(nVar);
  42         VV_Constant(0) = 0.;
  43         VV_Constant(1) = 155e5;
  44         for (int idim=0; idim<spaceDim;idim++){
  45                 VV_Constant(2+idim) = 0;
  46                 VV_Constant(2+idim +spaceDim) =0;
  47         }
  48         VV_Constant(2+spaceDim*nbPhase) = 600;
  49
  50         cout << "Number of Phases = " << nbPhase << endl;
  51         cout << "Building initial data " << endl;
  52
  53         // generate initial condition
  54         myProblem.setInitialFieldConstant(M,VV_Constant);
  55
  56         //set the boundary conditions
  57         myProblem.setBoundaryFields(boundaryFields);
  58         /* set physical parameters*/
  59 //      myProblem.setLatentHeat(latentHeat);
  60 //      myProblem.setSatPressure( Psat, dHsatl_over_dp);
  61
  62         // set the numerical method
  63         myProblem.setNumericalScheme(upwind, Explicit);
  64         myProblem.setEntropicCorrection(true);
  65
  66         // name file save
  67         string fileName = "1DDepressurisation";
  68
  69         // set numerical parameters
  70         unsigned MaxNbOfTimeStep =3;
  71         int freqSave = 1;
  72         double cfl = 0.5;
  73         double maxTime = 5;
  74         double precision = 1e-6;
  75
  76         myProblem.setCFL(cfl);
  77         myProblem.setPrecision(precision);
  78         myProblem.setMaxNbOfTimeStep(MaxNbOfTimeStep);
  79         myProblem.setTimeMax(maxTime);
  80         myProblem.setFreqSave(freqSave);
  81         myProblem.setFileName(fileName);
  82         bool ok;
  83
  84         // evolution
  85         myProblem.initialize();
  86
  87         ok = myProblem.run();
  88         if (ok)
  89                 cout << "Simulation "<<fileName<<" is successful !" << endl;
  90         else
  91                 cout << "Simulation "<<fileName<<"  failed ! " << endl;
  92
  93         cout << "------------ End of calculation -----------" << endl;
  94         myProblem.terminate();
  95
  96         return EXIT_SUCCESS;
  97 }