1 //============================================================================
2 // Name : ProblemCoreFlows
5 // Copyright : CEA Saclay 2014
6 // Description : Generic class for PDEs problems
7 //============================================================================
8 /* A ProblemCoreFlows class */
10 /*! \class ProblemCoreFlows ProblemCoreFlows.hxx "ProblemCoreFlows.hxx"
11 * \brief Generic class for thermal hydraulics problems
12 * \details Common functions to CoreFlows models
15 #ifndef PROBLEMCOREFLOWS_HXX_
16 #define PROBLEMCOREFLOWS_HXX_
33 #include "CdmathException.hxx"
37 //! enumeration linearSolver
38 /*! the linearSolver can be GMRES or BiCGStab (see Petsc documentation) */
41 GMRES,/**< linearSolver is GMRES */
42 BCGS,/**< linearSolver is BiCGSstab */
43 CG/**< linearSolver is CG */
46 //! enumeration preconditioner
47 /*! the preconditioner can be ILU or LU (see Petsc documentation) */
50 ILU,/**< preconditioner is ILU(0) */
51 LU,/**< preconditioner is actually a direct solver (LU factorisation)*/
52 NONE,/**< no preconditioner used */
53 ICC,/**< preconditioner is ICC(0) */
54 CHOLESKY/**< preconditioner is actually a direct solver for symmetric matrices (CHOLESKY factorisation)*/
57 //! enumeration nonLinearSolver
58 /*! the nonLinearSolver can be Newton_SOLVERLAB or using PETSc, Newton_PETSC_LINESEARCH, Newton_PETSC_TRUSTREGION (see Petsc documentation) */
61 Newton_SOLVERLAB,/**< nonLinearSolver is Newton_SOLVERLAB */
62 Newton_PETSC_LINESEARCH,/**< nonLinearSolver is Newton_PETSC_LINESEARCH */
63 Newton_PETSC_LINESEARCH_BASIC,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_BASIC */
64 Newton_PETSC_LINESEARCH_BT,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_BT */
65 Newton_PETSC_LINESEARCH_SECANT,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_SECANT */
66 Newton_PETSC_LINESEARCH_NLEQERR,/**< nonLinearSolver is Newton_PETSC_LINESEARCH_LEQERR */
67 Newton_PETSC_TRUSTREGION,/**< nonLinearSolver is Newton_PETSC_TRUSTREGION */
68 Newton_PETSC_NGMRES,/**< nonLinearSolver is Newton_PETSC_NGMRES */
69 Newton_PETSC_ASPIN/**< nonLinearSolver is Newton_PETSC_ASPIN */
72 //! enumeration saveFormat
73 /*! the numerical results are saved using MED, VTK or CSV format */
76 MED,/**< MED format is used */
77 VTK,/**< VTK format is used */
78 CSV/**< CSV format is used */
81 //! enumeration TimeScheme
82 /*! The numerical method can be Explicit or Implicit */
85 Explicit,/**< Explicit numerical scheme */
86 Implicit/**< Implicit numerical scheme */
89 class ProblemCoreFlows
92 //! Constructeur par défaut
93 ProblemCoreFlows(MPI_Comm comm = MPI_COMM_WORLD);
94 virtual ~ProblemCoreFlows();
96 // -*-*-*- Gestion du calcul (interface ICoCo) -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
99 * \brief Alloue la mémoire et vérifie que le maillage et les conditions limites/initiales sont bien définis
100 * \Details c'est une fonction virtuelle pure, on la surcharge dans le problem Fluide .
103 virtual void initialize()=0;
106 * \brief vide la mémoire et enregistre le résultat final
107 * \Details on la surcharge dans le problem fluid
110 virtual void terminate()=0;
112 /** \fn computeTimeStep
113 * \brief Propose un pas de temps pour le calcul
114 * \details Pour proposer un pas de temps la fonction nécessite de discrétiser les opérateurs
115 * convection, diffusion et sources. Pour chacun on emploie la condition cfl.
116 * En cas de problème durant ce calcul (exemple t=tmax), la fonction renvoie stop=true sinon stop = false
117 * Cest une fonction virtuelle pure, on la surcharge dans Problème Fulide
118 * @param Stop booléen correspond a True si ya un problème dans le code (exemple t=tmax) False sinon
119 * \return dt le pas de temps
121 virtual double computeTimeStep(bool & stop)=0;
124 * \brief Enregistre le nouveau pas de temps dt et l'intègre aux opérateurs
125 * discrets (convection, diffusion, sources)
126 * \details c'est une fonction virtuelle pure, on la surcharge dans le problemfluid
127 * @param dt est le nouvel pas de temps (double)
128 * \return false si dt <0 et True sinon
130 virtual bool initTimeStep(double dt)=0;
132 /** \fn solveTimeStep
133 * \brief calcule les valeurs inconnues au pas de temps +1 .
134 * \details c'est une fonction virtuelle
136 * \return Renvoie false en cas de problème durant le calcul (valeurs non physiques..)
138 virtual bool solveTimeStep();//
140 /** \fn validateTimeStep
141 * \brief Valide le calcule au temps courant
142 * \details met à jour le temps présent t=t+dt, sauvegarde les champs inconnus
143 * et reinitialise dt à 0, teste la stationnarité .
144 * c'est une fonction virtuel , on la surchage dans chacun des modèles
146 * \return Renvoie false en cas de problème durant le calcul
148 virtual void validateTimeStep()=0;//
150 /** \fn abortTimeStep
151 * \brief efface les inconnues calculées par solveTimeStep() et reinitialise dt à 0
152 * \details c'est une fonction virtuelle pure, elle est surchargée dans ProblemFluid, TransportEquation et DiffusionEquation
154 virtual void abortTimeStep()=0;
157 * \brief Vérifie tous les paramètres et lance le code en supposant que la cfl ne changera pas
158 * \details Renvoie "false" en cas de problème durant le calcul
159 * C'est une fonction virtuelle pure, on la surcharge dans problème Fluide
161 * \return Renvoie "false" en cas de problème durant le calcul
163 virtual bool run();//
165 /** \fn iterateTimeStep
166 * \brief Calcul d'une sous-itération du pas de temps en cours, typiquement une itération de Newton pour un schéma implicite
167 * \details Deux paramètres booléen (converged et ok) sont retournés
168 * converged, Vaut true si on peut passer au pas de temps suivant (shéma explicite ou convergence du schéma de Newton dans le schéma implicite)
169 * ok vaut true si le calcul n'a pas rencontré d'erreur ou de problème particulier dans la mise à jour des variables physiques
170 * \param [in] bool, passage par reférence.
173 virtual bool iterateTimeStep(bool &converged) = 0;
176 * \brief vérifie la stationnairité du problème .
177 * \details Renvoie "true" si le problème atteint l'état stationnaire et "false" sinon
181 virtual bool isStationary() const;
184 * \brief Calcule la valeur du temps courant .
187 * \param [out] double
189 virtual double presentTime() const;
191 /** \fn setStationaryMode
192 * \brief Perform the search of a stationary regime
197 virtual void setStationaryMode(bool stationaryMode){ _stationaryMode=stationaryMode;};
199 /** \fn getStationaryMode
200 * \brief Tells if we are seeking a stationary regime
205 virtual bool getStationaryMode(){return _stationaryMode;};
208 * \brief sets the current time (typically to start a new calculation)
213 void resetTime (double time);
217 virtual void getInputMEDDoubleFieldTemplate(const std::string& name, MEDDoubleField& afield) const;//Renvoie le format de champs attendu (maillage, composantes etc)
218 virtual vector<string> getOutputFieldsNames()=0 ;//liste tous les champs que peut fournir le code pour le postraitement
219 virtual void getOutputMEDDoubleField(const std::string& name, MEDDoubleField& afield)//Renvoie un champs pour le postraitement ou le couplage
222 /** Set input fields to prepare the simulation or coupling **/
223 virtual vector<string> getInputFieldsNames()=0;
224 virtual void setInputField(const string& nameField, Field& inputField )=0;//supply of a required input field
226 /*! @brief (Optional) Provide the code with a scalar double data.
228 * See Problem documentation for more details on the time semantic of a scalar value.
230 * @param[in] name name of the scalar value that is given to the code.
231 * @param[in] val value passed to the code.
232 * @throws ICoCo::WrongArgument exception if the scalar name ('name' parameter) is invalid.
234 /* virtual void setInputDoubleValue(const std::string& name, const double& val); */
236 /*! @brief (Optional) Retrieve a scalar double value from the code.
238 * See Problem documentation for more details on the time semantic of a scalar value.
240 * @param[in] name name of the scalar value to be read from the code.
241 * @return the double value read from the code.
242 * @throws ICoCo::WrongArgument exception if the scalar name ('name' parameter) is invalid.
244 /* virtual double getOutputDoubleValue(const std::string& name) const; */
246 /*! @brief (Optional) Similar to setInputDoubleValue() but for an int value.
247 * @sa setInputDoubleValue()
249 /* virtual void setInputIntValue(const std::string& name, const int& val); */
251 /*! @brief (Optional) Similar to getOutputDoubleValue() but for an int value.
252 * @sa getOutputDoubleValue()
254 /* virtual int getOutputIntValue(const std::string& name) const; */
256 /*! @brief (Optional) Similar to setInputDoubleValue() but for an string value.
257 * @sa setInputDoubleValue()
259 /* virtual void setInputStringValue(const std::string& name, const std::string& val); */
261 /*! @brief (Optional) Similar to getOutputDoubleValue() but for an string value.
262 * @sa getOutputDoubleValue()
264 /* virtual std::string getOutputStringValue(const std::string& name) const; */
266 /*! @brief Return ICoCo interface major version number.
267 * @return ICoCo interface major version number (2 at present)
269 static int GetICoCoMajorVersion() { return 2; }
271 /*! @brief (Optional) Get MEDCoupling major version, if the code was built with MEDCoupling support.
273 * This can be used to assess compatibility between codes when coupling them.
275 * @return the MEDCoupling major version number (typically 7, 8, 9, ...)
277 virtual int getMEDCouplingMajorVersion() const{ return MEDCOUPLING_VERSION_MAJOR; };
279 /*! @brief (Optional) Indicate whether the code was built with a 64-bits version of MEDCoupling.
281 * Implemented if the code was built with MEDCoupling support.
282 * This can be used to assess compatibility between codes when coupling them.
284 * @return the MEDCoupling major version number
286 virtual bool isMEDCoupling64Bits() const;
288 /*! @brief (Optional) Get the list of input scalars accepted by the code.
290 * @return the list of scalar names that represent inputs of the code
291 * @throws ICoCo::WrongContext exception if called before initialize() or after terminate().
293 /* virtual std::vector<std::string> getInputValuesNames() const; */
295 /*! @brief (Optional) Get the list of output scalars that can be provided by the code.
297 * @return the list of scalar names that can be returned by the code
298 * @throws ICoCo::WrongContext exception if called before initialize() or after terminate().
300 /* virtual std::vector<std::string> getOutputValuesNames() const; */
302 Field getUnknownField() const;
304 //paramètres du calcul -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
306 /** \fn setMaxNbOfTimeStep
307 * \brief met à jour _maxNbOfTimeStep ( le nombre maximum d'itération du calcul )
312 void setMaxNbOfTimeStep(int maxNbOfTimeStep);
315 * \brief met à jour _timeMax (Le temps maximum du calcul)
320 void setTimeMax(double timeMax);
323 * \brief met à jour la _CFL
328 void setCFL(double cfl);
331 * \brief met à jour _precision (la précision du calcule)
336 void setPrecision(double precision);
338 /** \fn setInitialField
339 * \brief sets the initial field
344 void setInitialField(const Field &VV);
346 /** \fn setInitialField
347 * \brief sets the initial field from a field in a med file.
348 * \details This function is added because we have not been able yet to swig properly the enum EntityType. It is replaced by an integer.
349 * \param [in] string : the file name
350 * \param [in] string : the field name
351 * \param [in] int : the time step number
352 * \param [in] int : int corresponding to the enum CELLS, NODES or FACES
355 void setInitialField(string fileName, string fieldName, int timeStepNumber, int order, int meshLevel, int field_support_type);
357 /** \fn setInitialField
358 * \brief sets the initial field from a field in a med file
360 * \param [in] string : the file name
361 * \param [in] string : the field name
362 * \param [in] int : the time step number
363 * \param [in] EntityType : CELLS, NODES or FACES
366 void setInitialField(string fileName, string fieldName, int timeStepNumber, int order = 0, int meshLevel=0, EntityType typeField = CELLS);
368 /** \fn setInitialFieldConstant
369 * \brief sets a constant initial field on a mesh stored in a med file
371 * \param [in] string : the file name
372 * \param [in] vector<double> : the value in each cell
373 * \param [in] EntityType : CELLS, NODES or FACES
376 void setInitialFieldConstant(string fileName, const vector<double> Vconstant, EntityType typeField = CELLS);
378 /** \fn setInitialFieldConstant
379 * \brief sets a constant initial field
383 * \param [in] EntityType : CELLS, NODES or FACES
386 void setInitialFieldConstant(const Mesh& M, const Vector Vconstant, EntityType typeField = CELLS);
388 /** \fn setInitialFieldConstant
389 * \brief sets a constant initial field
392 * \param [in] vector<double>
393 * \param [in] EntityType : CELLS, NODES or FACES
396 void setInitialFieldConstant(const Mesh& M, const vector<double> Vconstant, EntityType typeField = CELLS);
398 /** \fn setInitialFieldConstant
399 * \brief sets a constant initial field
401 * \param [in] int the space dimension
402 * \param [in] vector<double> the value in each cell
403 * \param [in] double the lowest value in the x direction
404 * \param [in] double the highest value in the x direction
405 * \param [in] string name of the left boundary
406 * \param [in] string name of the right boundary
407 * \param [in] double the lowest value in the y direction
408 * \param [in] double the highest value in the y direction
409 * \param [in] string name of the back boundary
410 * \param [in] string name of the front boundary
411 * \param [in] double the lowest value in the z direction
412 * \param [in] double the highest value in the z direction
413 * \param [in] string name of the bottom boundary
414 * \param [in] string name of the top boundary
415 * \param [in] EntityType : CELLS, NODES or FACES
418 void setInitialFieldConstant( int nDim, const vector<double> Vconstant, double xmin, double xmax,int nx, string leftSide, string rightSide,
419 double ymin=0, double ymax=0, int ny=0, string backSide="", string frontSide="",
420 double zmin=0, double zmax=0, int nz=0, string bottomSide="", string topSide="", EntityType typeField = CELLS);
422 /** \fn setInitialFieldConstant
423 * \brief sets a constant initial field
424 * \details This function is added because we have not been able yet to swig roperly the enum EntityType. It is replaced by an integer.
425 * \param [in] int the space dimension
426 * \param [in] vector<double> the value in each cell
427 * \param [in] double the lowest value in the x direction
428 * \param [in] double the highest value in the x direction
429 * \param [in] string name of the left boundary
430 * \param [in] string name of the right boundary
431 * \param [in] double the lowest value in the y direction
432 * \param [in] double the highest value in the y direction
433 * \param [in] string name of the back boundary
434 * \param [in] string name of the front boundary
435 * \param [in] double the lowest value in the z direction
436 * \param [in] double the highest value in the z direction
437 * \param [in] string name of the bottom boundary
438 * \param [in] string name of the top boundary
439 * \param [in] integer corresponding to the field support enum : CELLS, NODES or FACES
442 void setInitialFieldConstant( int nDim, const vector<double> Vconstant, double xmin, double xmax,int nx, string leftSide, string rightSide,
443 double ymin, double ymax, int ny, string backSide, string frontSide,
444 double zmin, double zmax, int nz, string bottomSide, string topSide, int type_of_field );
446 /** \fn setInitialFieldStepFunction
447 * \brief sets a step function initial field (Riemann problem)
452 * \param [in] double position of the discontinuity on one of the three axis
453 * \param [in] int direction (axis carrying the discontinuity) : 0 for x, 1 for y, 2 for z
454 * \param [in] EntityType : CELLS, NODES or FACES
457 void setInitialFieldStepFunction(const Mesh M, const Vector Vleft, const Vector Vright, double disc_pos, int direction=0, EntityType typeField = CELLS);
459 /** \fn setInitialFieldStepFunction
460 * \brief sets a constant initial field
462 * \param [in] int the space dimension
463 * \param [in] vector<double> the value left of the discontinuity
464 * \param [in] vector<double> the value right of the discontinuity
465 * \param [in] double the position of the discontinuity in the x direction
466 * \param [in] double the lowest value in the x direction
467 * \param [in] double the highest value in the x direction
468 * \param [in] string name of the left boundary
469 * \param [in] string name of the right boundary
470 * \param [in] double the lowest value in the y direction
471 * \param [in] double the highest value in the y direction
472 * \param [in] string name of the back boundary
473 * \param [in] string name of the front boundary
474 * \param [in] double the lowest value in the z direction
475 * \param [in] double the highest value in the z direction
476 * \param [in] string name of the bottom boundary
477 * \param [in] string name of the top boundary
478 * \param [in] EntityType : CELLS, NODES or FACES
481 void setInitialFieldStepFunction( int nDim, const vector<double> VV_Left, vector<double> VV_Right, double xstep,
482 double xmin, double xmax,int nx, string leftSide, string rightSide,
483 double ymin=0, double ymax=0, int ny=0, string backSide="", string frontSide="",
484 double zmin=0, double zmax=0, int nz=0, string bottomSide="", string topSide="", EntityType typeField = CELLS);
486 /** \fn setInitialFieldSphericalStepFunction
487 * \brief sets a step function initial field with value Vin inside the ball with radius Radius and Vout outside
490 * \param [in] Vector Vin, value inside the ball
491 * \param [in] Vector Vout, value outside the ball
492 * \param [in] double radius of the ball
493 * \param [in] Vector Center, coordinates of the ball center
494 * \param [in] EntityType : CELLS, NODES or FACES
497 void setInitialFieldSphericalStepFunction(const Mesh M, const Vector Vin, const Vector Vout, double Radius, Vector Center, EntityType typeField = CELLS);
500 * \brief renvoie _time (le temps courant du calcul)
503 * \param [out] double
507 /** \fn getNbTimeStep
508 * \brief renvoie _nbTimeStep le Numéro d'itération courant
511 * \param [out] unsigned
513 unsigned getNbTimeStep();
516 * \brief renvoie la _CFL
519 * \param [out] double
524 * \brief renvoie _precision (la précision du calcul)
527 * \param [out] double
529 double getPrecision();
532 * \brief renvoie _Mesh (le maillage du problème)
540 * \brief met à jour _fileName le nom du fichier
545 void setFileName(string fileName);
548 * \brief met à jour _FreqSave (la fréquence du sauvgarde de la solution)
553 void setFreqSave(int freqSave);
556 * \brief sauvgarde les données dans des fichiers MED ou VTK
557 * \details c'est une fonction virtuelle pure , on la surcharge
558 * dans chacun des modèles
561 virtual void save() = 0;
563 /** \fn getLinearSolver
564 * \brief renvoie _ksptype (le type du solveur linéaire utilisé)
567 * \param [out] string
569 string getLinearSolver() {
573 /** \fn getNonLinearSolver
574 * \brief renvoie _nonLinearSolver (le type du solveur de Newton utilisé)
577 * \param [out] string
579 nonLinearSolver getNonLinearSolver() {
580 return _nonLinearSolver;
583 /** \fn getNumberOfVariables
584 * \brief le nombre d'inconnues du problème
587 * \return renvoie _nVar (le nombre d'inconnues du problème)
589 int getNumberOfVariables(){
593 /** \fn setWellBalancedCorrection
594 * \brief include a well balanced correction to treat stiff source terms
595 * @param boolean that is true if a well balanced correction should be applied
597 void setWellBalancedCorrection(bool wellBalancedCorr){
598 _wellBalancedCorrection=wellBalancedCorr;
601 /** \fn setLinearSolver
602 * \brief sets the linear solver and preconditioner
603 * \details virtual function overloaded by intanciable classes
604 * @param kspType linear solver type (GMRES or BICGSTAB)
605 * @param pcType preconditioner (ILU,LU or NONE)
607 void setLinearSolver(linearSolver solverName, preconditioner pcType, double maxIts=50);
609 /** \fn setNewtonSolver
610 * \brief sets the Newton type algorithm for solving the nonlinear algebraic system arising from the discretisation of the PDE
611 * \param [in] double : precision required for the convergence of the newton scheme
612 * \param [in] int : maximum number of newton iterations
613 * \param [in] nonLinearSolver : the algorithm to be used to solve the nonlinear system
616 void setNewtonSolver(double precision, int iterations=20, nonLinearSolver solverName=Newton_SOLVERLAB);
618 /** \fn displayConditionNumber
619 * \brief display the condition number of the preconditioned linear systems
621 void displayConditionNumber(bool display=true){
622 _conditionNumber=display;
625 /** \fn setSaveFileFormat
626 * \brief sets the numerical results file format (MED, VTK or CSV)
628 * \param [in] saveFormat
631 void setSaveFileFormat(saveFormat saveFileFormat){
632 _saveFormat=saveFileFormat;
635 /** \fn setResultDirectory
636 * \brief sets the directory where the results will be saved
638 * \param [in] resultsPath
641 void setResultDirectory(string resultsPath){
645 /** \fn getTimeScheme
646 * \brief returns the time scheme name
648 * \param [out] enum TimeScheme (explicit or implicit)
650 TimeScheme getTimeScheme();
652 /** \fn setNumericalScheme
653 * \brief sets the numerical method ( explicit vs implicit )
655 * \param [in] TimeScheme
658 void setTimeScheme( TimeScheme method);
660 //Couplages Thermohydraulique-thermique-neutronique *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
662 /** \fn setHeatPowerField
663 * \brief set the heat power field (variable in space)
668 void setHeatPowerField(Field heatPower);
670 /** \fn setHeatPowerField
671 * \brief set the heat power field (variable in space)
673 * \param [in] string fileName (including file path)
674 * \param [in] string fieldName
677 void setHeatPowerField(string fileName, string fieldName, int iteration = 0, int order = 0, int meshLevel=0);
679 /** \fn setHeatSource
680 * \brief sets a constant heat power field
685 void setHeatSource(double phi){
687 _isStationary=false;//Source term may be changed after previously reaching a stationary state
690 /** \fn getHeatPowerField
691 * \brief returns the heat power field
696 Field getHeatPowerField(){
697 return _heatPowerField;
700 /** \fn setHeatTransfertCoeff
701 * \brief set the heat transfert coefficient for heat exchange between fluid and solid
706 void setHeatTransfertCoeff(double heatTransfertCoeff){
707 _heatTransfertCoeff=heatTransfertCoeff;
708 _isStationary=false;//Source term may be changed after previously reaching a stationary state
712 * \brief Updates display options
718 void setVerbose(bool verbose, bool system=false)
725 double getConditionNumber(bool isSingular=false, double tol=1e-6) const;
726 std::vector< double > getEigenvalues (int nev, EPSWhich which=EPS_SMALLEST_MAGNITUDE, double tol=1e-6) const;
727 std::vector< Vector > getEigenvectors(int nev, EPSWhich which=EPS_SMALLEST_MAGNITUDE, double tol=1e-6) const;
728 Field getEigenvectorsField(int nev, EPSWhich which=EPS_SMALLEST_MAGNITUDE, double tol=1e-6) const;
729 std::vector< double > getSingularValues( int nsv, SVDWhich which=SVD_SMALLEST, double tol=1e-6) const;
730 std::vector< Vector > getSingularVectors(int nsv, SVDWhich which=SVD_SMALLEST, double tol=1e-6) const;
732 // some supplementary functions
734 /** \fn displayMatrix
735 * \brief displays a matrix of size "size x size" for profiling
736 * @param matrix is a pointer of size "size"
737 * @param size, size of the matrix
738 * @param name, string, name or description of the matrix
739 * @return displays the matrix on the terminal
741 static void displayMatrix(double *matrix, int size, string name="Vector coefficients :");
743 /** \fn displayMatrix
744 * \brief displays a vector of size "size" for profiling
745 * @param vector is a pointer of size "size"
746 * @param size, size of the vector
747 * @param name, string, name or description of the vector
748 * @return displays the vector on the terminal
750 static void displayVector(double *vector, int size, string name="Matrix coefficients :");
754 int _Ndim;//space dimension
755 int _nVar;//Number of equations to sole
756 int _Nmailles;//number of cells
757 int _Nnodes;//number of nodes
758 int _Nfaces;//number of faces
759 int _neibMaxNbCells;//maximum number of neighbours around a cell
760 int _neibMaxNbNodes;/* maximum number of nodes around a node */
770 double _maxvp;//valeur propre max pour calcul cfl
771 double _minl;//minimum cell diameter
773 /** boolean used to specify that a well balanced correction should be used */
774 bool _wellBalancedCorrection;
775 TimeScheme _timeScheme;
777 //Linear solver and petsc
783 nonLinearSolver _nonLinearSolver;
784 int _maxPetscIts;//nombre maximum d'iteration gmres autorise au cours d'une resolution de systeme lineaire
785 int _PetscIts;//the number of iterations of the linear solver
786 int _maxNewtonIts;//nombre maximum d'iteration de Newton autorise au cours de la resolution d'un pas de temps
788 Mat _A;//Linear system matrix
789 Vec _b;//Linear system right hand side
790 double _MaxIterLinearSolver;//nombre maximum d'iteration gmres obtenu au cours par les resolution de systemes lineaires au cours d'un pas de tmeps
791 bool _conditionNumber;//computes an estimate of the condition number
793 //simulation monitoring variables
795 bool _initialDataSet;
796 bool _initializedMemory;
797 bool _stationaryMode;//ICoCo V2
798 bool _restartWithNewTimeScheme;
799 bool _restartWithNewFileName;
800 double _timeMax,_time;
801 int _maxNbOfTimeStep,_nbTimeStep;
803 double _precision_Newton;
804 double _erreur_rel;//norme(Un+1-Un)
805 string _fileName;//name of the calculation
807 ofstream * _runLogFile;//for creation of a log file to save the history of the simulation
809 //Heat transfert variables
810 Field _heatPowerField;
811 bool _heatPowerFieldSet;
812 double _heatTransfertCoeff;
814 double _hsatv, _hsatl;//all models appart from DiffusionEquation will need this
817 bool _verbose, _system;
819 string _path;//path to execution directory used for saving results
820 saveFormat _saveFormat;//file saving format : MED, VTK or CSV
822 //MPI related variables
823 PetscMPIInt _mpi_size; /* size of communicator */
824 PetscMPIInt _mpi_rank; /* processor rank */
825 VecScatter _scat; /* For the distribution of a local vector */
826 int _globalNbUnknowns, _localNbUnknowns;
827 int _d_nnz, _o_nnz; /* local and "non local" numbers of non zeros corfficients */
830 #endif /* PROBLEMCOREFLOWS_HXX_ */