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  * dT/dt - \lambda\Delta T = \Phi + \lambda_{sf} (T_{fluid}-T)
   9  * Dirichlet (imposed temperature) or Neumann (imposed normal flux) boundary conditions.
  10  * */
  11 //============================================================================
  12
  13 /*! \class DiffusionEquation DiffusionEquation.hxx "DiffusionEquation.hxx"
  14  *  \brief Scalar heat equation for the Uranium rods temperature
  15  *  \details see \ref DiffusionEqPage for more details
  16  * dT/dt - \lambda\Delta T = \Phi + \lambda_{sf} (T_{fluid}-T)
  17  */
  18 #ifndef DiffusionEquation_HXX_
  19 #define DiffusionEquation_HXX_
  20
  21 #include "ProblemCoreFlows.hxx"
  22 #include "Node.hxx"
  23
  24 using namespace std;
  25
  26 //! enumeration BoundaryType
  27 /*! Boundary condition type  */
  28 enum BoundaryTypeDiffusion      { NeumannDiffusion, DirichletDiffusion, NoneBCDiffusion};
  29
  30 /** \struct LimitField
  31  * \brief value of some fields on the boundary  */
  32 struct LimitFieldDiffusion{
  33         LimitFieldDiffusion(){bcType=NoneBCDiffusion; T=0; normalFlux=0;}
  34         LimitFieldDiffusion(BoundaryTypeDiffusion _bcType, double _T,   double _normalFlux){
  35                 bcType=_bcType; T=_T; normalFlux=_normalFlux;
  36         }
  37
  38         BoundaryTypeDiffusion bcType;
  39         double T; //for Dirichlet
  40         double normalFlux; //for Neumann
  41 };
  42
  43 class DiffusionEquation: public ProblemCoreFlows
  44 {
  45
  46 public :
  47         /** \fn DiffusionEquation
  48                          * \brief Constructor for the temperature diffusion in a solid
  49                          * \param [in] int : space dimension
  50                          * \param [in] double : solid density
  51                          * \param [in] double : solid specific heat at constant pressure
  52                          * \param [in] double : solid conductivity
  53                          *  */
  54
  55         DiffusionEquation( int dim,bool FECalculation=true,double rho=10000,double cp=300,double lambda=5);
  56
  57         //Gestion du calcul
  58         void initialize();
  59         void terminate();//vide la mémoire et enregistre le résultat final
  60         bool initTimeStep(double dt);
  61         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
  62         void abortTimeStep();//efface les inconnues calculées par solveTimeStep() et reinitialise dt à 0
  63         bool iterateTimeStep(bool &ok);
  64         void save();
  65         void validateTimeStep();
  66
  67     /* Boundary conditions */
  68         void setBoundaryFields(map<string, LimitFieldDiffusion> boundaryFields){
  69                 _limitField = boundaryFields;
  70     };
  71         /** \fn setDirichletBoundaryCondition
  72                          * \brief adds a new boundary condition of type Dirichlet
  73                          * \details
  74                          * \param [in] string : the name of the boundary
  75                          * \param [in] double : the value of the temperature at the boundary
  76                          * \param [out] void
  77                          *  */
  78         void setDirichletBoundaryCondition(string groupName,double Temperature=0){
  79                 _limitField[groupName]=LimitFieldDiffusion(DirichletDiffusion,Temperature,-1);
  80         };
  81         /** \fn setNeumannBoundaryCondition
  82                          * \brief adds a new boundary condition of type Neumann
  83                          * \details
  84                          * \param [in] string : the name of the boundary
  85                          * \param [out] void
  86                          *  */
  87         void setNeumannBoundaryCondition(string groupName, double normalFlux=0){
  88                 _limitField[groupName]=LimitFieldDiffusion(NeumannDiffusion,-1, normalFlux);
  89         };
  90
  91         void setDirichletValues(map< int, double> dirichletBoundaryValues);
  92         void setNeumannValues(map< int, double> neumannBoundaryValues);
  93
  94         void setRodDensity(double rho){
  95                 _rho=rho;
  96         };
  97         void setConductivity(double conductivite){
  98                 _conductivity=conductivite;
  99         };
 100         void setFluidTemperatureField(Field coupledTemperatureField){
 101                 _fluidTemperatureField=coupledTemperatureField;
 102                 _fluidTemperatureFieldSet=true;
 103         };
 104
 105         void setDiffusiontensor(Matrix DiffusionTensor){
 106                 _DiffusionTensor=DiffusionTensor;
 107         };
 108
 109         void setFluidTemperature(double fluidTemperature){
 110         _fluidTemperature=fluidTemperature;
 111         }
 112
 113         //get output fields for postprocessing or coupling
 114         vector<string> getOutputFieldsNames() ;//liste tous les champs que peut fournir le code pour le postraitement
 115         Field&         getOutputField(const string& nameField );//Renvoie un champs pour le postraitement
 116
 117         Field& getRodTemperatureField(){
 118                 return _VV;
 119         }
 120         Field& getFluidTemperatureField(){
 121                 return _fluidTemperatureField;
 122         }
 123
 124 protected :
 125         double computeDiffusionMatrix(bool & stop);
 126         double computeDiffusionMatrixFV(bool & stop);
 127         double computeRHS(bool & stop);
 128
 129         Field _fluidTemperatureField;
 130         bool _fluidTemperatureFieldSet, _diffusionMatrixSet;
 131         double _conductivity,_diffusivity, _fluidTemperature;
 132         double _rho;
 133         double _cp;
 134         Vector _normale;
 135         Matrix _DiffusionTensor;
 136         Vec _Tn, _deltaT, _Tk, _Tkm1, _b0;
 137         double _dt_diffusion, _dt_src;
 138
 139     /************ Data for FE calculation *************/
 140     bool _FECalculation;
 141         int _neibMaxNbNodes;/* maximum number of nodes around a node */
 142         int _NunknownNodes;/* number of unknown nodes for FE calculation */
 143         int _NboundaryNodes;/* total number of boundary nodes */
 144         int _NdirichletNodes;/* number of boundary nodes with Dirichlet BC for FE calculation */
 145     std::vector< int > _boundaryNodeIds;/* List of boundary nodes */
 146     std::vector< int > _dirichletNodeIds;/* List of boundary nodes with Dirichlet BC */
 147
 148     /*********** Functions for finite element method ***********/
 149     static Vector gradientNodal(Matrix M, vector< double > v);//gradient of nodal shape functions
 150         double computeDiffusionMatrixFE(bool & stop);
 151     static int fact(int n);
 152     static int unknownNodeIndex(int globalIndex, std::vector< int > dirichletNodes);
 153     static int globalNodeIndex(int unknownIndex, std::vector< int > dirichletNodes);
 154
 155         TimeScheme _timeScheme;
 156         map<string, LimitFieldDiffusion> _limitField;
 157
 158     /********* Possibility to set a boundary field as DirichletNeumann boundary condition *********/
 159     bool _dirichletValuesSet;
 160     bool _neumannValuesSet;
 161     std::map< int, double> _dirichletBoundaryValues;
 162     std::map< int, double> _neumannBoundaryValues;
 163     };
 164
 165 #endif /* DiffusionEquation_HXX_ */