1 //============================================================================
3 * \file StationaryDiffusionEquation.hxx
4 * \author Michael NDJINGA
7 * \brief Stationary heat diffusion equation solved with either finite elements or finite volume method.
8 * -\lambda\Delta T=\Phi + \lambda_{sf} (T_{fluid}-T)
10 //============================================================================
12 /*! \class StationaryDiffusionEquation StationaryDiffusionEquation.hxx "StationaryDiffusionEquation.hxx"
13 * \brief Scalar stationary heat equation solved with either finite elements or finite volume method.
14 * \details see \ref StationaryDiffusionEqPage for more details
15 * -\lambda\Delta T=\Phi(T) + \lambda_{sf} (T_{fluid}-T)
17 #ifndef StationaryDiffusionEquation_HXX_
18 #define StationaryDiffusionEquation_HXX_
20 #include "ProblemCoreFlows.hxx"
25 class StationaryDiffusionEquation
29 /** \fn StationaryDiffusionEquation
30 * \brief Constructor for the temperature diffusion in a solid
31 * \param [in] int : space dimension
32 * \param [in] double : solid density
33 * \param [in] double : solid specific heat at constant pressure
34 * \param [in] double : solid conductivity
37 StationaryDiffusionEquation( int dim,bool FECalculation=true,double lambda=1);
39 void setMesh(const Mesh &M);
40 void setLinearSolver(linearSolver kspType, preconditioner pcType);
41 void setFileName(string fileName){
44 bool solveStationaryProblem();
45 Field getOutputTemperatureField();
49 void terminate();//vide la mémoire et enregistre le résultat final
50 double computeDiffusionMatrix(bool & stop);
51 double computeTimeStep(bool & stop);//For coupling calculations
52 bool iterateNewtonStep(bool &ok);
55 /* Boundary conditions */
56 void setBoundaryFields(map<string, LimitField> boundaryFields){
57 _limitField = boundaryFields;
59 /** \fn setDirichletBoundaryCondition
60 * \brief adds a new boundary condition of type Dirichlet
62 * \param [in] string : the name of the boundary
63 * \param [in] double : the value of the temperature at the boundary
66 void setDirichletBoundaryCondition(string groupName,double Temperature){
67 _limitField[groupName]=LimitField(Dirichlet,-1, vector<double>(_Ndim,0),vector<double>(_Ndim,0),
68 vector<double>(_Ndim,0),Temperature,-1,-1,-1);
71 /** \fn setNeumannBoundaryCondition
72 * \brief adds a new boundary condition of type Neumann
74 * \param [in] string : the name of the boundary
77 void setNeumannBoundaryCondition(string groupName){
78 _limitField[groupName]=LimitField(Neumann,-1, vector<double>(0),vector<double>(0),
79 vector<double>(0),-1,-1,-1,-1);
82 void setDirichletValues(map< int, double> dirichletBoundaryValues);
84 void setConductivity(double conductivite){
85 _conductivity=conductivite;
87 void setFluidTemperatureField(Field coupledTemperatureField){
88 _fluidTemperatureField=coupledTemperatureField;
89 _fluidTemperatureFieldSet=true;
91 void setFluidTemperature(double fluidTemperature){
92 _fluidTemperature=fluidTemperature;
94 Field& getRodTemperatureField(){
97 Field& getFluidTemperatureField(){
98 return _fluidTemperatureField;
100 void setDiffusiontensor(Matrix DiffusionTensor){
101 _DiffusionTensor=DiffusionTensor;
104 /** \fn setHeatPowerField
105 * \brief set the heat power field (variable in space)
110 void setHeatPowerField(Field heatPower){
111 _heatPowerField=heatPower;
112 _heatPowerFieldSet=true;
115 /** \fn setHeatPowerField
116 * \brief set the heat power field (variable in space)
118 * \param [in] string fileName (including file path)
119 * \param [in] string fieldName
122 void setHeatPowerField(string fileName, string fieldName){
123 _heatPowerField=Field(fileName, CELLS,fieldName);
124 _heatPowerFieldSet=true;
131 int _Ndim;//space dimension
132 int _nVar;//Number of equations to solve=1
137 bool _initializedMemory;
138 int _Nmailles;//number of cells for FV calculation
139 int _neibMaxNbCells;//maximum number of cells around a cell
142 double _precision_Newton;
143 double _erreur_rel;//norme(Uk+1-Uk)
144 bool _computationCompletedSuccessfully;
146 //Linear solver and petsc
152 int _maxPetscIts;//nombre maximum d'iteration gmres autorisé au cours d'une resolution de système lineaire
153 int _maxNewtonIts;//nombre maximum d'iteration de Newton autorise au cours de la resolution d'un pas de temps
154 int _PetscIts;//the number of iterations of the linear solver
156 Mat _A;//Linear system matrix
157 Vec _b;//Linear system right hand side
158 double _MaxIterLinearSolver;//nombre maximum d'iteration gmres obtenu au cours par les resolution de systemes lineaires au cours d'un pas de tmeps
159 bool _conditionNumber;//computes an estimate of the condition number
161 map<string, LimitField> _limitField;
162 bool _onlyNeumannBC;//if true then the linear system is singular and should be solved up to a constant vector
164 bool _diffusionMatrixSet;
166 Matrix _DiffusionTensor;
167 Vec _deltaT, _Tk, _Tkm1, _b0;
170 //Heat transfert variables
171 double _conductivity, _fluidTemperature;
172 Field _heatPowerField, _fluidTemperatureField;
173 bool _heatPowerFieldSet, _fluidTemperatureFieldSet;
174 double _heatTransfertCoeff, _heatSource;
177 bool _verbose, _system;
178 ofstream * _runLogFile;//for creation of a log file to save the history of the simulation
180 string _fileName;//name of the calculation
181 string _path;//path to execution directory used for saving results
182 saveFormat _saveFormat;//file saving format : MED, VTK or CSV
184 double computeRHS(bool & stop);
185 double computeDiffusionMatrixFV(bool & stop);
187 /************ Data for FE calculation *************/
189 int _Nnodes;/* number of nodes for FE calculation */
190 int _neibMaxNbNodes;/* maximum number of nodes around a node */
191 int _NunknownNodes;/* number of unknown nodes for FE calculation */
192 int _NboundaryNodes;/* total number of boundary nodes */
193 int _NdirichletNodes;/* number of boundary nodes with Dirichlet BC for FE calculation */
194 std::vector< int > _boundaryNodeIds;/* List of boundary nodes */
195 std::vector< int > _dirichletNodeIds;/* List of boundary nodes with Dirichlet BC */
197 /*********** Functions for finite element method ***********/
198 Vector gradientNodal(Matrix M, vector< double > v);//gradient of nodal shape functions
199 double computeDiffusionMatrixFE(bool & stop);
201 int unknownNodeIndex(int globalIndex, std::vector< int > dirichletNodes);
202 int globalNodeIndex(int unknownIndex, std::vector< int > dirichletNodes);
204 /********* Possibility to set a boundary field as Dirichlet boundary condition *********/
205 bool _dirichletValuesSet;
206 std::map< int, double> _dirichletBoundaryValues;
209 #endif /* StationaryDiffusionEquation_HXX_ */