CHOLESKY/**< preconditioner is actually a direct solver for symmetric matrices (CHOLESKY factorisation)*/
};
+//! enumeration nonLinearSolver
+/*! the nonLinearSolver can be Newton_SOLVERLAB or using PETSc, Newton_PETSC_LINESEARCH, Newton_PETSC_TRUSTREGION (see Petsc documentation) */
+enum nonLinearSolver
+{
+ Newton_SOLVERLAB,/**< nonLinearSolver is Newton_SOLVERLAB */
+ Newton_PETSC_LINESEARCH,/**< nonLinearSolver is Newton_PETSC_LINESEARCH */
+ Newton_PETSC_LINESEARCH_BASIC,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_BASIC */
+ Newton_PETSC_LINESEARCH_BT,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_BT */
+ Newton_PETSC_LINESEARCH_SECANT,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_SECANT */
+ Newton_PETSC_LINESEARCH_NLEQERR,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_LEQERR */
+ Newton_PETSC_TRUSTREGION,/**< nonLinearSolver is Newton_PETSC_TRUSTREGION */
+ Newton_PETSC_NGMRES,/**< nonLinearSolver is Newton_PETSC_NGMRES */
+ Newton_PETSC_ASPIN/**< nonLinearSolver is Newton_PETSC_ASPIN */
+};
+
//! enumeration saveFormat
/*! the numerical results are saved using MED, VTK or CSV format */
enum saveFormat
Implicit/**< Implicit numerical scheme */
};
+//! enumeration pressureEstimate
+/*! the pressure estimate needed to fit physical parameters */
+enum pressureEstimate
+{
+ around1bar300K,/**< pressure is around 1 bar and temperature around 300K (for TransportEquation, SinglePhase and IsothermalTwoFluid) or 373 K (saturation for DriftModel and FiveEqsTwoFluid) */
+ around155bars600K/**< pressure is around 155 bars and temperature around 618 K (saturation) */
+};
+
class ProblemCoreFlows
{
public :
//! Constructeur par défaut
- ProblemCoreFlows();
+ ProblemCoreFlows(MPI_Comm comm = MPI_COMM_WORLD);
virtual ~ProblemCoreFlows();
- // -*-*-*- Gestion du calcul (interface ICoCo -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
+ // -*-*-*- Gestion du calcul (interface ICoCo) -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
/** \fn initialize
* \brief Alloue la mémoire et vérifie que le maillage et les conditions limites/initiales sont bien définis
* \param [in] bool, passage par reférence.
* \param [out] bool
* */
- virtual bool iterateTimeStep(bool &converged) = 0; //??
+ virtual bool iterateTimeStep(bool &converged) = 0;
/** \fn isStationary
* \brief vérifie la stationnairité du problème .
* */
virtual double presentTime() const;
+ /** \fn setStationaryMode
+ * \brief Perform the search of a stationary regime
+ * \details
+ * \param [in] bool
+ * \param [out]
+ * */
+ virtual void setStationaryMode(bool stationaryMode){ _stationaryMode=stationaryMode;};
+
+ /** \fn getStationaryMode
+ * \brief Tells if we are seeking a stationary regime
+ * \details
+ * \param [in]
+ * \param [out] bool
+ * */
+ virtual bool getStationaryMode(){return _stationaryMode;};
+
+ /** \fn resetTime
+ * \brief sets the current time (typically to start a new calculation)
+ * \details
+ * \param [in] double
+ * \param [out] void
+ * */
+ void resetTime (double time);
+
/*
//Coupling interface
- virtual vector<string> getInputFieldsNames()=0 ;//Renvoie les noms des champs dont le problème a besoin (données initiales)
- virtual Field& getInputFieldTemplate(const string& name)=0;//Renvoie le format de champs attendu (maillage, composantes etc)
- virtual void setInputField(const string& name, const Field& afield)=0;//enregistre les valeurs d'une donnée initiale
+ virtual void getInputMEDDoubleFieldTemplate(const std::string& name, MEDDoubleField& afield) const;//Renvoie le format de champs attendu (maillage, composantes etc)
virtual vector<string> getOutputFieldsNames()=0 ;//liste tous les champs que peut fournir le code pour le postraitement
- virtual Field& getOutputField(const string& nameField )=0;//Renvoie un champs pour le postraitement
+ virtual void getOutputMEDDoubleField(const std::string& name, MEDDoubleField& afield)//Renvoie un champs pour le postraitement ou le couplage
*/
+ /** Set input fields to prepare the simulation or coupling **/
+ virtual vector<string> getInputFieldsNames()=0;
+ virtual void setInputField(const string& nameField, Field& inputField )=0;//supply of a required input field
+
+ /*! @brief (Optional) Provide the code with a scalar double data.
+ *
+ * See Problem documentation for more details on the time semantic of a scalar value.
+ *
+ * @param[in] name name of the scalar value that is given to the code.
+ * @param[in] val value passed to the code.
+ * @throws ICoCo::WrongArgument exception if the scalar name ('name' parameter) is invalid.
+ */
+ /* virtual void setInputDoubleValue(const std::string& name, const double& val); */
+
+ /*! @brief (Optional) Retrieve a scalar double value from the code.
+ *
+ * See Problem documentation for more details on the time semantic of a scalar value.
+ *
+ * @param[in] name name of the scalar value to be read from the code.
+ * @return the double value read from the code.
+ * @throws ICoCo::WrongArgument exception if the scalar name ('name' parameter) is invalid.
+ */
+ /* virtual double getOutputDoubleValue(const std::string& name) const; */
+
+ /*! @brief (Optional) Similar to setInputDoubleValue() but for an int value.
+ * @sa setInputDoubleValue()
+ */
+ /* virtual void setInputIntValue(const std::string& name, const int& val); */
+
+ /*! @brief (Optional) Similar to getOutputDoubleValue() but for an int value.
+ * @sa getOutputDoubleValue()
+ */
+ /* virtual int getOutputIntValue(const std::string& name) const; */
+
+ /*! @brief (Optional) Similar to setInputDoubleValue() but for an string value.
+ * @sa setInputDoubleValue()
+ */
+ /* virtual void setInputStringValue(const std::string& name, const std::string& val); */
+
+ /*! @brief (Optional) Similar to getOutputDoubleValue() but for an string value.
+ * @sa getOutputDoubleValue()
+ */
+ /* virtual std::string getOutputStringValue(const std::string& name) const; */
+
+ /*! @brief Return ICoCo interface major version number.
+ * @return ICoCo interface major version number (2 at present)
+ */
+ static int GetICoCoMajorVersion() { return 2; }
+
+ /*! @brief (Optional) Get MEDCoupling major version, if the code was built with MEDCoupling support.
+ *
+ * This can be used to assess compatibility between codes when coupling them.
+ *
+ * @return the MEDCoupling major version number (typically 7, 8, 9, ...)
+ */
+ virtual int getMEDCouplingMajorVersion() const{ return MEDCOUPLING_VERSION_MAJOR; };
+
+ /*! @brief (Optional) Indicate whether the code was built with a 64-bits version of MEDCoupling.
+ *
+ * Implemented if the code was built with MEDCoupling support.
+ * This can be used to assess compatibility between codes when coupling them.
+ *
+ * @return the MEDCoupling major version number
+ */
+ virtual bool isMEDCoupling64Bits() const;
+
+ /*! @brief (Optional) Get the list of input scalars accepted by the code.
+ *
+ * @return the list of scalar names that represent inputs of the code
+ * @throws ICoCo::WrongContext exception if called before initialize() or after terminate().
+ */
+ /* virtual std::vector<std::string> getInputValuesNames() const; */
+
+ /*! @brief (Optional) Get the list of output scalars that can be provided by the code.
+ *
+ * @return the list of scalar names that can be returned by the code
+ * @throws ICoCo::WrongContext exception if called before initialize() or after terminate().
+ */
+ /* virtual std::vector<std::string> getOutputValuesNames() const; */
+
Field getUnknownField() const;
//paramètres du calcul -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
- /** \fn setPresentTime
- * \brief met à jour _time (le temps courant du calcul)
- * \details
- * \param [in] double
- * \param [out] void
- * */
- void setPresentTime (double time);
-
/** \fn setMaxNbOfTimeStep
* \brief met à jour _maxNbOfTimeStep ( le nombre maximum d'itération du calcul )
* \details
* */
void setInitialField(const Field &VV);
+ /** \fn setInitialField
+ * \brief sets the initial field from a field in a med file.
+ * \details This function is added because we have not been able yet to swig properly the enum EntityType. It is replaced by an integer.
+ * \param [in] string : the file name
+ * \param [in] string : the field name
+ * \param [in] int : the time step number
+ * \param [in] int : int corresponding to the enum CELLS, NODES or FACES
+ * \param [out] void
+ * */
+ void setInitialField(string fileName, string fieldName, int timeStepNumber, int order, int meshLevel, int field_support_type);
+
/** \fn setInitialField
* \brief sets the initial field from a field in a med file
* \details
* \param [in] string : the file name
* \param [in] string : the field name
* \param [in] int : the time step number
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
- void setInitialField(string fileName, string fieldName, int timeStepNumber);
+ void setInitialField(string fileName, string fieldName, int timeStepNumber, int order = 0, int meshLevel=0, EntityType typeField = CELLS);
/** \fn setInitialFieldConstant
* \brief sets a constant initial field on a mesh stored in a med file
* \details
* \param [in] string : the file name
* \param [in] vector<double> : the value in each cell
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
- void setInitialFieldConstant(string fileName, const vector<double> Vconstant);
+ void setInitialFieldConstant(string fileName, const vector<double> Vconstant, EntityType typeField = CELLS);
/** \fn setInitialFieldConstant
* \brief sets a constant initial field
* \details
* \param [in] Mesh
* \param [in] Vector
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
- void setInitialFieldConstant(const Mesh& M, const Vector Vconstant);
+ void setInitialFieldConstant(const Mesh& M, const Vector Vconstant, EntityType typeField = CELLS);
/** \fn setInitialFieldConstant
* \brief sets a constant initial field
* \details
* \param [in] Mesh
* \param [in] vector<double>
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
- void setInitialFieldConstant(const Mesh& M, const vector<double> Vconstant);
+ void setInitialFieldConstant(const Mesh& M, const vector<double> Vconstant, EntityType typeField = CELLS);
/** \fn setInitialFieldConstant
* \brief sets a constant initial field
* \param [in] double the highest value in the z direction
* \param [in] string name of the bottom boundary
* \param [in] string name of the top boundary
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
void setInitialFieldConstant( int nDim, const vector<double> Vconstant, double xmin, double xmax,int nx, string leftSide, string rightSide,
double ymin=0, double ymax=0, int ny=0, string backSide="", string frontSide="",
- double zmin=0, double zmax=0, int nz=0, string bottomSide="", string topSide="");
+ double zmin=0, double zmax=0, int nz=0, string bottomSide="", string topSide="", EntityType typeField = CELLS);
+
+ /** \fn setInitialFieldConstant
+ * \brief sets a constant initial field
+ * \details This function is added because we have not been able yet to swig roperly the enum EntityType. It is replaced by an integer.
+ * \param [in] int the space dimension
+ * \param [in] vector<double> the value in each cell
+ * \param [in] double the lowest value in the x direction
+ * \param [in] double the highest value in the x direction
+ * \param [in] string name of the left boundary
+ * \param [in] string name of the right boundary
+ * \param [in] double the lowest value in the y direction
+ * \param [in] double the highest value in the y direction
+ * \param [in] string name of the back boundary
+ * \param [in] string name of the front boundary
+ * \param [in] double the lowest value in the z direction
+ * \param [in] double the highest value in the z direction
+ * \param [in] string name of the bottom boundary
+ * \param [in] string name of the top boundary
+ * \param [in] integer corresponding to the field support enum : CELLS, NODES or FACES
+ * \param [out] void
+ * */
+ void setInitialFieldConstant( int nDim, const vector<double> Vconstant, double xmin, double xmax,int nx, string leftSide, string rightSide,
+ double ymin, double ymax, int ny, string backSide, string frontSide,
+ double zmin, double zmax, int nz, string bottomSide, string topSide, int type_of_field );
/** \fn setInitialFieldStepFunction
* \brief sets a step function initial field (Riemann problem)
* \param [in] Vector
* \param [in] double position of the discontinuity on one of the three axis
* \param [in] int direction (axis carrying the discontinuity) : 0 for x, 1 for y, 2 for z
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
- void setInitialFieldStepFunction(const Mesh M, const Vector Vleft, const Vector Vright, double disc_pos, int direction=0);
+ void setInitialFieldStepFunction(const Mesh M, const Vector Vleft, const Vector Vright, double disc_pos, int direction=0, EntityType typeField = CELLS);
/** \fn setInitialFieldStepFunction
* \brief sets a constant initial field
* \param [in] double the highest value in the z direction
* \param [in] string name of the bottom boundary
* \param [in] string name of the top boundary
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
void setInitialFieldStepFunction( int nDim, const vector<double> VV_Left, vector<double> VV_Right, double xstep,
double xmin, double xmax,int nx, string leftSide, string rightSide,
double ymin=0, double ymax=0, int ny=0, string backSide="", string frontSide="",
- double zmin=0, double zmax=0, int nz=0, string bottomSide="", string topSide="");
+ double zmin=0, double zmax=0, int nz=0, string bottomSide="", string topSide="", EntityType typeField = CELLS);
/** \fn setInitialFieldSphericalStepFunction
* \brief sets a step function initial field with value Vin inside the ball with radius Radius and Vout outside
* \param [in] Vector Vout, value outside the ball
* \param [in] double radius of the ball
* \param [in] Vector Center, coordinates of the ball center
+ * \param [in] EntityType : CELLS, NODES or FACES
* \param [out] void
* */
- void setInitialFieldSphericalStepFunction(const Mesh M, const Vector Vin, const Vector Vout, double Radius, Vector Center);
+ void setInitialFieldSphericalStepFunction(const Mesh M, const Vector Vin, const Vector Vout, double Radius, Vector Center, EntityType typeField = CELLS);
/** \fn getTime
* \brief renvoie _time (le temps courant du calcul)
return _ksptype;
};
+ /** \fn getNonLinearSolver
+ * \brief renvoie _nonLinearSolver (le type du solveur de Newton utilisé)
+ * \details
+ * \param [in] void
+ * \param [out] string
+ * */
+ nonLinearSolver getNonLinearSolver() {
+ return _nonLinearSolver;
+ };
+
/** \fn getNumberOfVariables
* \brief le nombre d'inconnues du problème
* \details
* @param kspType linear solver type (GMRES or BICGSTAB)
* @param pcType preconditioner (ILU,LU or NONE)
*/
- void setLinearSolver(linearSolver solverName, preconditioner pcType);
+ void setLinearSolver(linearSolver solverName, preconditioner pcType, double maxIts=50);
/** \fn setNewtonSolver
- * \brief set the Newton algorithm parameters
- * \param [in] int maximum number of newton iterations
- * \param [in] double precision required for the convergence of the newton scheme
+ * \brief sets the Newton type algorithm for solving the nonlinear algebraic system arising from the discretisation of the PDE
+ * \param [in] double : precision required for the convergence of the newton scheme
+ * \param [in] int : maximum number of newton iterations
+ * \param [in] nonLinearSolver : the algorithm to be used to solve the nonlinear system
* \param [out] void
* */
- void setNewtonSolver(double precision,int iterations=20)
- {
- _maxNewtonIts=iterations;
- _precision_Newton=precision;
- };
+ void setNewtonSolver(double precision, int iterations=20, nonLinearSolver solverName=Newton_SOLVERLAB);
/** \fn displayConditionNumber
* \brief display the condition number of the preconditioned linear systems
* \param [in] Field
* \param [out] void
* */
- void setHeatPowerField(Field heatPower){
- _heatPowerField=heatPower;
- _heatPowerFieldSet=true;
- }
+ void setHeatPowerField(Field heatPower);
/** \fn setHeatPowerField
* \brief set the heat power field (variable in space)
* \param [in] string fieldName
* \param [out] void
* */
- void setHeatPowerField(string fileName, string fieldName){
- _heatPowerField=Field(fileName, CELLS,fieldName);
- _heatPowerFieldSet=true;
- }
+ void setHeatPowerField(string fileName, string fieldName, int iteration = 0, int order = 0, int meshLevel=0);
/** \fn setHeatSource
- * \brief met à jour la puissance thermique ( _phi )
+ * \brief sets a constant heat power field
* \details
* \param [in] double
* \param [out] void
* */
void setHeatSource(double phi){
_heatSource=phi;
+ _isStationary=false;//Source term may be changed after previously reaching a stationary state
}
/** \fn getHeatPowerField
- * \brief renvoie le champs ?? ( _heatPowerField )
+ * \brief returns the heat power field
* \details
* \param [in] void
* \param [out] Field
return _heatPowerField;
}
- /** \fn setRodTemperatureField ??
- * \brief
- * \details
- * \param [in] Field
- * \param [out] void
- * */
- void setRodTemperatureField(Field rodTemperature){
- _rodTemperatureField=rodTemperature;
- _rodTemperatureFieldSet=true;
- }
-
- /** \fn setRodTemperature ??
- * \brief
- * \details
- * \param [in] double
- * \param [out] void
- * */
- void setRodTemperature(double rodTemp){
- _rodTemperature=rodTemp;
- }
-
- /** \fn getRodTemperatureField
- * \brief
- * \details
- * \param [in] void
- * \param [out] Field
- * */
- virtual Field& getRodTemperatureField(){ // ?? je ne retrouve pas cet attribut dans le file.cxx
- return _rodTemperatureField;
- }
-
/** \fn setHeatTransfertCoeff
* \brief set the heat transfert coefficient for heat exchange between fluid and solid
* \details
* */
void setHeatTransfertCoeff(double heatTransfertCoeff){
_heatTransfertCoeff=heatTransfertCoeff;
+ _isStationary=false;//Source term may be changed after previously reaching a stationary state
}
- /** \fn setDISPLAY
- * \brief met à jour les paramètres de l'affichage
+ /** \fn setVerbose
+ * \brief Updates display options
* \details
* \param [in] bool
* \param [in] bool
int _Nmailles;//number of cells
int _Nnodes;//number of nodes
int _Nfaces;//number of faces
- int _neibMaxNb;//maximum number of neighbours around a cell
+ int _neibMaxNbCells;//maximum number of neighbours around a cell
int _neibMaxNbNodes;/* maximum number of nodes around a node */
Mesh _mesh;
Field _perimeters;
PC _pc;
PCType _pctype;
string _pc_hypre;
+ nonLinearSolver _nonLinearSolver;
int _maxPetscIts;//nombre maximum d'iteration gmres autorise au cours d'une resolution de systeme lineaire
int _PetscIts;//the number of iterations of the linear solver
int _maxNewtonIts;//nombre maximum d'iteration de Newton autorise au cours de la resolution d'un pas de temps
int _NEWTON_its;
- Mat _A;//Linear system matrix
+ Mat _A;//Linear system matrix
Vec _b;//Linear system right hand side
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
+ /** \fn createKSP
+ * \brief Create PETSc solver and preconditioner structures
+ * */
+ void createKSP();
//simulation monitoring variables
bool _isStationary;
bool _initialDataSet;
bool _initializedMemory;
+ bool _stationaryMode;//ICoCo V2
bool _restartWithNewTimeScheme;
bool _restartWithNewFileName;
double _timeMax,_time;
ofstream * _runLogFile;//for creation of a log file to save the history of the simulation
//Heat transfert variables
- Field _heatPowerField, _rodTemperatureField;
- bool _heatPowerFieldSet, _rodTemperatureFieldSet;
+ Field _heatPowerField;
+ bool _heatPowerFieldSet;
double _heatTransfertCoeff;
- double _heatSource, _rodTemperature;
+ double _heatSource;
double _hsatv, _hsatl;//all models appart from DiffusionEquation will need this
//Display variables
string _path;//path to execution directory used for saving results
saveFormat _saveFormat;//file saving format : MED, VTK or CSV
+
+ //MPI related variables
+ PetscMPIInt _mpi_size; /* size of communicator */
+ PetscMPIInt _mpi_rank; /* processor rank */
+ VecScatter _scat; /* For the distribution of a local vector */
+ int _globalNbUnknowns, _localNbUnknowns;
+ int _d_nnz, _o_nnz; /* local and "non local" numbers of non zeros corfficients */
};
#endif /* PROBLEMCOREFLOWS_HXX_ */
