CoreFlows/examples/DiffusionEquation_1DHeatedRod.cxx

   1 #include "DiffusionEquation.hxx"\r
   2 \r
   3 using namespace std;\r
   4 \r
   5 #define PI 3.14159265\r
   6 \r
   7 void power_field_diffusionTest(Field & Phi){\r
   8         double L=4.2;\r
   9         double lambda=0.2;\r
  10         double phi=1e5;\r
  11         double x;\r
  12         Mesh M = Phi.getMesh();\r
  13         int nbCells = M.getNumberOfCells();\r
  14         for (int j = 0; j < nbCells; j++) {\r
  15                 x=M.getCell(j).x();\r
  16                 Phi(j) = phi*cos(PI*(x-L/2)/(L+lambda));\r
  17         }\r
  18 }\r
  19 \r
  20 int main(int argc, char** argv)\r
  21 {\r
  22         //Preprocessing: mesh and group creation\r
  23         double xinf=0.0;\r
  24         double xsup=4.2;\r
  25         int nx=10;\r
  26         cout << "Building of a 1D mesh with "<<nx<<" cells" << endl;\r
  27         Mesh M(xinf,xsup,nx);\r
  28         double eps=1.E-6;\r
  29         M.setGroupAtPlan(xsup,0,eps,"Neumann");\r
  30         M.setGroupAtPlan(xinf,0,eps,"Neumann");\r
  31         int spaceDim = M.getSpaceDimension();\r
  32 \r
  33 \r
  34         //Solid parameters\r
  35         double cp_ur=300;//Uranium specific heat\r
  36         double rho_ur=10000;//Uranium density\r
  37         double lambda_ur=5;\r
  38  \r
  39     bool FEcalculation=false;\r
  40         DiffusionEquation  myProblem(spaceDim,FEcalculation,rho_ur,cp_ur,lambda_ur);\r
  41         Field VV("Solid temperature", CELLS, M, 1);\r
  42 \r
  43         //Set fluid temperature (temperature du fluide)\r
  44         double fluidTemp=573;//fluid mean temperature\r
  45         double heatTransfertCoeff=1000;//fluid/solid exchange coefficient\r
  46         myProblem.setFluidTemperature(fluidTemp);\r
  47         myProblem.setHeatTransfertCoeff(heatTransfertCoeff);\r
  48         //Set heat source\r
  49         Field Phi("Heat power field", CELLS, M, 1);\r
  50         power_field_diffusionTest(Phi);\r
  51         myProblem.setHeatPowerField(Phi);\r
  52         Phi.writeVTK("1DheatPowerField");\r
  53 \r
  54         //Initial field creation\r
  55         Vector VV_Constant(1);\r
  56         VV_Constant(0) = 623;//Rod clad temperature\r
  57 \r
  58         cout << "Building initial data" << endl;\r
  59         myProblem.setInitialFieldConstant(M,VV_Constant);\r
  60 \r
  61         //set the boundary conditions\r
  62         myProblem.setNeumannBoundaryCondition("Neumann");\r
  63 \r
  64         // set the numerical method\r
  65         myProblem.setNumericalScheme(upwind, Explicit);\r
  66 \r
  67         // name result file\r
  68         string fileName = "1DRodTemperature_FV";\r
  69 \r
  70         // parameters calculation\r
  71         unsigned MaxNbOfTimeStep =3;\r
  72         int freqSave = 1;\r
  73         double cfl = 0.5;\r
  74         double maxTime = 1000000;\r
  75         double precision = 1e-6;\r
  76 \r
  77         myProblem.setCFL(cfl);\r
  78         myProblem.setPrecision(precision);\r
  79         myProblem.setMaxNbOfTimeStep(MaxNbOfTimeStep);\r
  80         myProblem.setTimeMax(maxTime);\r
  81         myProblem.setFreqSave(freqSave);\r
  82         myProblem.setFileName(fileName);\r
  83 \r
  84         // set display option to monitor the calculation\r
  85         myProblem.setVerbose( true);\r
  86         //set file saving format\r
  87         myProblem.setSaveFileFormat(CSV);\r
  88 \r
  89         // evolution\r
  90         myProblem.initialize();\r
  91         bool ok = myProblem.run();\r
  92         if (ok)\r
  93                 cout << "Simulation "<<fileName<<" is successful !" << endl;\r
  94         else\r
  95                 cout << "Simulation "<<fileName<<"  failed ! " << endl;\r
  96 \r
  97         cout << "------------ End of calculation -----------" << endl;\r
  98         myProblem.terminate();\r
  99 \r
 100         return EXIT_SUCCESS;\r
 101 }\r