1 //============================================================================
3 * \file TransportEquation.hxx
4 * \author Michael NDJINGA
7 * \brief Fluid enthalpy transport equation
9 //============================================================================
11 /*! \class TransportEquation TransportEquation.hxx "TransportEquation.hxx"
12 * \brief Scalar advection equation for a fluid enthalpy
13 * \details see \ref TransportEqPage for more details
15 #ifndef TransportEquation_HXX_
16 #define TransportEquation_HXX_
18 #include "ProblemCoreFlows.hxx"
23 //! enumeration BoundaryType
24 /*! Boundary condition type */
25 enum BoundaryTypeTransport {InletTransport, OutletTransport, NeumannTransport, DirichletTransport, NoneBCTransport};//Actually Inlet=Dirichlet and Outlet=Neumann
27 //! enumeration phaseType
28 /*! The material phase can be Solid, Gas or liquid */
31 Air,/**< Material considered is air */
32 Water/**< Material considered is water */
35 /** \struct LimitField
36 * \brief value of some fields on the boundary */
37 struct LimitFieldTransport{
38 LimitFieldTransport(){bcType=NoneBCTransport; T=0; h=0; flux=0; }
39 LimitFieldTransport(BoundaryTypeTransport _bcType, double _T, double _h,double _flux ){
40 bcType=_bcType; T=_T; h=_h; flux=_flux;
43 BoundaryTypeTransport bcType;
44 double T; //for inlet or Dirichlet
45 double h; //for inlet or Dirichlet
46 double flux; //for Neumann or outlet
49 class TransportEquation: public ProblemCoreFlows
53 /** \fn TransportEquation
54 * \brief Constructor for the enthalpy transport in a fluid
55 * \param [in] phase : \ref Liquid or \ref Gas
56 * \param [in] pressureEstimate : \ref around1bar or \ref around155bars
57 * \param [in] vector<double> : fluid velocity (assumed constant)
59 TransportEquation(FluidMaterial fluid, pressureEstimate pEstimate,vector<double> vitesseTransport, MPI_Comm comm = MPI_COMM_WORLD);
62 virtual void initialize();
63 virtual void terminate();//vide la mémoire et enregistre le résultat final
64 bool initTimeStep(double dt);
65 double computeTimeStep(bool & stop);//propose un pas de temps pour le calcul. Celà nécessite de discrétiser les opérateur (convection, diffusion, sources) et pour chacun d'employer la condition cfl. En cas de problème durant ce calcul (exemple t=tmax), renvoie stop=true
66 void abortTimeStep();//efface les inconnues calculées par solveTimeStep() et reinitialise dt à 0
67 bool iterateTimeStep(bool &ok);
69 virtual void validateTimeStep();
71 /* Boundary conditions */
72 /** \fn setIntletBoundaryCondition
73 * \brief adds a new boundary condition of type Inlet
74 * \details same as setDirichletBoundaryCondition
75 * \param [in] string : the name of the boundary
76 * \param [in] double : the value of the enthalpy at the boundary
79 void setInletBoundaryCondition(string groupName,double enthalpy){
80 _limitField[groupName]=LimitFieldTransport(InletTransport,-1,enthalpy,-1);
82 /** \fn setDirichletBoundaryCondition
83 * \brief adds a new boundary condition of type Dirichlet
84 * \details same as setInletBoundaryCondition
85 * \param [in] string : the name of the boundary
86 * \param [in] double : the value of the enthalpy at the boundary
89 void setDirichletBoundaryCondition(string groupName,double enthalpy){
90 _limitField[groupName]=LimitFieldTransport(DirichletTransport,-1,enthalpy,-1);
92 /** \fn setDirichletBoundaryCondition
93 * \brief adds a new boundary condition of type Inlet
94 * \details Reads the boundary field in a med file
95 * \param [in] string : the name of the boundary
96 * \param [in] string : the file name
97 * \param [in] string : the field name
98 * \param [in] int : the time step number
99 * \param [in] int : int corresponding to the enum CELLS or NODES
102 void setDirichletBoundaryCondition(string groupName, string fileName, string fieldName, int timeStepNumber, int order, int meshLevel, int field_support_type);
103 void setDirichletBoundaryCondition(string groupName, Field bc_field){
104 _limitField[groupName]=LimitFieldTransport(DirichletTransport, -1, 0, -1);
107 /** \fn setNeumannBoundaryCondition
108 * \brief adds a new boundary condition of type Neumann
109 * \details same as setOutletBoundaryCondition
110 * \param [in] string the name of the boundary
111 * \param [in] double : the value of the enthalpy flux at the boundary
114 void setNeumannBoundaryCondition(string groupName, double flux=0){
115 _limitField[groupName]=LimitFieldTransport(NeumannTransport,-1,-1,flux);
117 /** \fn setOutletBoundaryCondition
118 * \brief adds a new boundary condition of type Outlet
119 * \details same as setNeumannBoundaryCondition
120 * \param [in] string the name of the boundary
121 * \param [in] double : the value of the enthalpy flux at the boundary
124 void setOutletBoundaryCondition(string groupName, double flux=0){
125 _limitField[groupName]=LimitFieldTransport(OutletTransport,-1,-1,flux);
127 /** \fn setNeumannBoundaryCondition
128 * \brief adds a new boundary condition of type Neumann
129 * \details Reads the boundary field in a med file
130 * \param [in] string : the name of the boundary
131 * \param [in] string : the file name
132 * \param [in] string : the field name
133 * \param [in] int : the time step number
134 * \param [in] int : int corresponding to the enum CELLS or NODES
137 void setNeumannBoundaryCondition(string groupName, string fileName, string fieldName, int timeStepNumber, int order, int meshLevel, int field_support_type);
138 void setNeumannBoundaryCondition(string groupName, Field bc_field){
139 _limitField[groupName]=LimitFieldTransport(NeumannTransport,-1,-1, 0);
142 /** \fn setBoundaryFields
143 * \brief met à jour _limitField ( le type de condition limite )
148 void setBoundaryFields(map<string, LimitFieldTransport> boundaryFields){
149 _limitField = boundaryFields;
153 /*Physical parameters*/
154 void setLiqSatEnthalpy(double hsatl){
157 void setVapSatEnthalpy(double hsatv){
160 void setLiqSatDensity(double rhosatl){
163 void setVapSatDensity(double rhosatv){
166 void setTransportVelocity(Vector v){
170 /* set input fields to prepare the simulation */
171 vector<string> getInputFieldsNames();
172 void setInputField(const string& nameField, Field& inputField );//supply of a required input field
174 /** \fn setRodTemperatureField
175 * \brief Set the rod temperature field
180 void setRodTemperatureField(Field rodTemperature);
182 /** \fn setRodTemperature
183 * \brief Set a constant rod temperature field
188 void setRodTemperature(double rodTemp){
189 _rodTemperature=rodTemp;
190 _isStationary=false;//Source term may be changed after previously reaching a stationary state
193 /** \fn getRodTemperatureField
199 Field& getRodTemperatureField(){ // ?? je ne retrouve pas cet attribut dans le file.cxx
200 return _rodTemperatureField;
203 /* get output fields for postprocessing or coupling */
204 vector<string> getOutputFieldsNames() ;//liste tous les champs que peut fournir le code pour le postraitement
205 Field& getOutputField(const string& nameField );//Renvoie un champs pour le postraitement
207 Field& getFluidTemperatureField(){
211 Field& getEnthalpyField(){
215 Field& getVoidFractionField(){
219 Field& getDensityField(){
224 double computeTransportMatrix();
226 void updatePrimitives();
228 /* Postprocessing fields */
229 Field _TT, _Alpha, _Rho;//Fields of temperature, void fraction, density. Unknown field is enthalpy (_VV)
230 double _rhosatv, _rhosatl;
231 double _Tref, _href, _cpref;
233 double temperature(double h){
234 return _Tref+(h-_href)/_cpref;
236 double voidFraction(double h){
237 double titre=(h-_href)/(_hsatv-_hsatl);
242 else return titre*_rhosatl/(titre*_rhosatl+(1-titre)*_rhosatv);
244 double density(double alpha){
245 return alpha*_rhosatv+(1-alpha)*_rhosatl;
248 Vector _vitesseTransport, _normale;
249 bool _transportMatrixSet;
250 Vec _Hn, _deltaH, _Hk, _Hkm1, _b0;
251 Vec _Hn_seq; // Local sequential copy of the parallel vector _Hn, used for saving result files
252 double _dt_transport, _dt_src;
254 map<string, LimitFieldTransport> _limitField;
257 bool _rodTemperatureFieldSet;
258 Field _rodTemperatureField;
259 double _rodTemperature;
262 #endif /* TransportEquation_HXX_ */