* \date June 2019
* \brief Stationary heat diffusion equation solved with either finite elements or finite volume method.
* -\lambda\Delta T=\Phi + \lambda_{sf} (T_{fluid}-T)
+ * Dirichlet (imposed temperature) or Neumann (imposed normal flux) boundary conditions.
* */
//============================================================================
* \details see \ref StationaryDiffusionEqPage for more details
* -\lambda\Delta T=\Phi(T) + \lambda_{sf} (T_{fluid}-T)
*/
+
#ifndef StationaryDiffusionEquation_HXX_
#define StationaryDiffusionEquation_HXX_
#include "ProblemCoreFlows.hxx"
-#include "Node.hxx"
+
+/* For the laplacian spectrum */
+#include <slepceps.h>
+#include <slepcsvd.h>
using namespace std;
+/*! Boundary condition type */
+enum BoundaryTypeStationaryDiffusion { NeumannStationaryDiffusion, DirichletStationaryDiffusion, NoneBCStationaryDiffusion};
+
+/** \struct LimitField
+ * \brief value of some fields on the boundary */
+struct LimitFieldStationaryDiffusion{
+ LimitFieldStationaryDiffusion(){bcType=NoneBCStationaryDiffusion; T=0; normalFlux=0;}
+ LimitFieldStationaryDiffusion(BoundaryTypeStationaryDiffusion _bcType, double _T, double _normalFlux){
+ bcType=_bcType; T=_T; normalFlux=_normalFlux;
+ }
+
+ BoundaryTypeStationaryDiffusion bcType;
+ double T; //for Dirichlet
+ double normalFlux; //for Neumann
+};
+
class StationaryDiffusionEquation
{
/** \fn StationaryDiffusionEquation
* \brief Constructor for the temperature diffusion in a solid
* \param [in] int : space dimension
- * \param [in] double : solid density
- * \param [in] double : solid specific heat at constant pressure
+ * \param [in] bool : numerical method
* \param [in] double : solid conductivity
* */
StationaryDiffusionEquation( int dim,bool FECalculation=true,double lambda=1);
void setMesh(const Mesh &M);
- void setLinearSolver(linearSolver kspType, preconditioner pcType);
void setFileName(string fileName){
_fileName = fileName;
}
bool solveStationaryProblem();
Field getOutputTemperatureField();
+ //Linear system and spectrum
+ void setLinearSolver(linearSolver kspType, preconditioner pcType);
+ double getConditionNumber(bool isSingular=false, double tol=1e-6) const;
+ std::vector< double > getEigenvalues (int nev, EPSWhich which=EPS_SMALLEST_MAGNITUDE, double tol=1e-6) const;
+ std::vector< Vector > getEigenvectors(int nev, EPSWhich which=EPS_SMALLEST_MAGNITUDE, double tol=1e-6) const;
+ Field getEigenvectorsField(int nev, EPSWhich which=EPS_SMALLEST_MAGNITUDE, double tol=1e-6) const;
+
//Gestion du calcul
void initialize();
void terminate();//vide la mémoire et enregistre le résultat final
void save();
/* Boundary conditions */
- void setBoundaryFields(map<string, LimitField> boundaryFields){
+ void setBoundaryFields(map<string, LimitFieldStationaryDiffusion> boundaryFields){
_limitField = boundaryFields;
};
/** \fn setDirichletBoundaryCondition
* \param [out] void
* */
void setDirichletBoundaryCondition(string groupName,double Temperature){
- _limitField[groupName]=LimitField(Dirichlet,-1, vector<double>(_Ndim,0),vector<double>(_Ndim,0),
- vector<double>(_Ndim,0),Temperature,-1,-1,-1);
+ _limitField[groupName]=LimitFieldStationaryDiffusion(DirichletStationaryDiffusion,Temperature,-1);
};
/** \fn setNeumannBoundaryCondition
* \brief adds a new boundary condition of type Neumann
* \details
* \param [in] string : the name of the boundary
+ * \param [in] double : outward normal flux
* \param [out] void
* */
- void setNeumannBoundaryCondition(string groupName){
- _limitField[groupName]=LimitField(Neumann,-1, vector<double>(0),vector<double>(0),
- vector<double>(0),-1,-1,-1,-1);
+ void setNeumannBoundaryCondition(string groupName, double normalFlux=0){
+ _limitField[groupName]=LimitFieldStationaryDiffusion(NeumannStationaryDiffusion,-1, normalFlux);
};
void setDirichletValues(map< int, double> dirichletBoundaryValues);
+ void setNeumannValues (map< int, double> neumannBoundaryValues);
void setConductivity(double conductivite){
_conductivity=conductivite;
double _MaxIterLinearSolver;//nombre maximum d'iteration gmres obtenu au cours par les resolution de systemes lineaires au cours d'un pas de tmeps
bool _conditionNumber;//computes an estimate of the condition number
- map<string, LimitField> _limitField;
+ map<string, LimitFieldStationaryDiffusion> _limitField;
bool _onlyNeumannBC;//if true then the linear system is singular and should be solved up to a constant vector
bool _diffusionMatrixSet;
Vector gradientNodal(Matrix M, vector< double > v);//gradient of nodal shape functions
double computeDiffusionMatrixFE(bool & stop);
int fact(int n);
- int unknownNodeIndex(int globalIndex, std::vector< int > dirichletNodes);
- int globalNodeIndex(int unknownIndex, std::vector< int > dirichletNodes);
+ int unknownNodeIndex(int globalIndex, std::vector< int > dirichletNodes) const;
+ int globalNodeIndex(int unknownIndex, std::vector< int > dirichletNodes) const;
- /********* Possibility to set a boundary field as Dirichlet boundary condition *********/
+ /********* Possibility to set a boundary field as DirichletNeumann boundary condition *********/
bool _dirichletValuesSet;
+ bool _neumannValuesSet;
std::map< int, double> _dirichletBoundaryValues;
+ std::map< int, double> _neumannBoundaryValues;
};
#endif /* StationaryDiffusionEquation_HXX_ */