CoreFlows/Models/inc/DiffusionEquation.hxx

   1 //============================================================================
   2 /**
   3  * \file DiffusionEquation.hxx
   4  * \author Michael NDJINGA
   5  * \version 1.0
   6  * \date 24 March 2015
   7  * \brief Heat diffusion equation
   8  * */
   9 //============================================================================
  10
  11 /*! \class DiffusionEquation DiffusionEquation.hxx "DiffusionEquation.hxx"
  12  *  \brief Scalar heat equation for the Uranium rods temperature
  13  *  \details see \ref DiffusionEqPage for more details
  14  */
  15 #ifndef DiffusionEquation_HXX_
  16 #define DiffusionEquation_HXX_
  17
  18 #include "ProblemCoreFlows.hxx"
  19 #include "Node.hxx"
  20
  21 using namespace std;
  22
  23 class DiffusionEquation: public ProblemCoreFlows
  24 {
  25
  26 public :
  27         /** \fn DiffusionEquation
  28                          * \brief Constructor for the temperature diffusion in a solid
  29                          * \param [in] int : space dimension
  30                          * \param [in] double : solid density
  31                          * \param [in] double : solid specific heat at constant pressure
  32                          * \param [in] double : solid conductivity
  33                          *  */
  34
  35         DiffusionEquation( int dim,bool FECalculation=true,double rho=10000,double cp=300,double lambda=5);
  36
  37         //Gestion du calcul
  38         void initialize();
  39         void terminate();//vide la mémoire et enregistre le résultat final
  40         bool initTimeStep(double dt);
  41         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
  42         void abortTimeStep();//efface les inconnues calculées par solveTimeStep() et reinitialise dt à 0
  43         bool iterateTimeStep(bool &ok);
  44         void save();
  45         void validateTimeStep();
  46
  47     /* Boundary conditions */
  48         void setBoundaryFields(map<string, LimitField> boundaryFields){
  49                 _limitField = boundaryFields;
  50     };
  51         /** \fn setDirichletBoundaryCondition
  52                          * \brief adds a new boundary condition of type Dirichlet
  53                          * \details
  54                          * \param [in] string : the name of the boundary
  55                          * \param [in] double : the value of the temperature at the boundary
  56                          * \param [out] void
  57                          *  */
  58         void setDirichletBoundaryCondition(string groupName,double Temperature){
  59                 _limitField[groupName]=LimitField(Dirichlet,-1,vector<double>(_Ndim,0),vector<double>(_Ndim,0),vector<double>(_Ndim,0),Temperature,-1,-1,-1);
  60         };
  61         /** \fn setNeumannBoundaryCondition
  62                          * \brief adds a new boundary condition of type Neumann
  63                          * \details
  64                          * \param [in] string : the name of the boundary
  65                          * \param [out] void
  66                          *  */
  67         void setNeumannBoundaryCondition(string groupName){
  68                 _limitField[groupName]=LimitField(Neumann,-1, vector<double>(0),vector<double>(0),
  69                                                       vector<double>(0),-1,-1,-1,-1);
  70         };
  71
  72         void setRodDensity(double rho){
  73                 _rho=rho;
  74         };
  75         void setConductivity(double conductivite){
  76                 _conductivity=conductivite;
  77         };
  78         void setFluidTemperatureField(Field coupledTemperatureField){
  79                 _fluidTemperatureField=coupledTemperatureField;
  80                 _fluidTemperatureFieldSet=true;
  81         };
  82         void setFluidTemperature(double fluidTemperature){
  83         _fluidTemperature=fluidTemperature;
  84         }
  85         Field& getRodTemperatureField(){
  86                 return _VV;
  87         }
  88         Field& getFluidTemperatureField(){
  89                 return _fluidTemperatureField;
  90         }
  91         void setDiffusiontensor(Matrix DiffusionTensor){
  92                 _DiffusionTensor=DiffusionTensor;
  93         };
  94 protected :
  95         double computeDiffusionMatrix(bool & stop);
  96         double computeDiffusionMatrixFV(bool & stop);
  97         double computeRHS(bool & stop);
  98
  99         Field _fluidTemperatureField;
 100         bool _fluidTemperatureFieldSet, _diffusionMatrixSet;
 101         double _conductivity,_diffusivity, _fluidTemperature;
 102         double _rho;
 103         double _cp;
 104         Vector _normale;
 105         Matrix _DiffusionTensor;
 106         Vec _Tn, _deltaT, _Tk, _Tkm1, _b0;
 107         double _dt_diffusion, _dt_src;
 108
 109     /************ Data for FE calculation *************/
 110     bool _FECalculation;
 111         int _NunknownNodes;/* number of unknown nodes for FE calculation */
 112         int _NboundaryNodes;/* total number of boundary nodes */
 113         int _NdirichletNodes;/* number of boundary nodes with Dirichlet BC for FE calculation */
 114     std::vector< int > _boundaryNodeIds;/* List of boundary nodes */
 115     std::vector< int > _dirichletNodeIds;/* List of boundary nodes with Dirichlet BC */
 116
 117     /*********** Functions for finite element method ***********/
 118     Vector gradientNodal(Matrix M, vector< double > v);//gradient of nodal shape functions
 119         double computeDiffusionMatrixFE(bool & stop);
 120     int fact(int n);
 121     int unknownNodeIndex(int globalIndex, std::vector< int > dirichletNodes);
 122     int globalNodeIndex(int unknownIndex, std::vector< int > dirichletNodes);
 123
 124 };
 125
 126 #endif /* DiffusionEquation_HXX_ */