//============================================================================
#include "ProblemCoreFlows.hxx"
+#include "SparseMatrixPetsc.hxx"
+
#include <limits.h>
#include <unistd.h>
_freqSave = 1;
_initialDataSet=false;
_initializedMemory=false;
- _spaceScheme=upwind;
+ _restartWithNewTimeScheme=false;
+ _restartWithNewFileName=false;
_timeScheme=Explicit;
_wellBalancedCorrection=false;
+ _FECalculation=false;
_maxPetscIts=50;
_MaxIterLinearSolver=0;//During several newton iterations, stores the max petssc interations
_maxNewtonIts=50;
int _PetscIts=0;//the number of iterations of the linear solver
_ksptype = (char*)&KSPGMRES;
_pctype = (char*)&PCLU;
+ _nonLinearSolver = Newton_SOLVERLAB;
_heatPowerFieldSet=false;
_heatTransfertCoeff=0;
_rodTemperatureFieldSet=false;
_saveFormat=VTK;
}
+TimeScheme ProblemCoreFlows::getTimeScheme()
+{
+ return _timeScheme;
+}
+
+void ProblemCoreFlows::setTimeScheme(TimeScheme timeScheme)
+{
+ if( _nbTimeStep>0 && timeScheme!=_timeScheme)//This is a change of time scheme during a simulation
+ _restartWithNewTimeScheme=true;
+ _timeScheme = timeScheme;
+}
+
bool ProblemCoreFlows::isStationary() const
{
return _isStationary;
{
_precision=precision;
}
-void ProblemCoreFlows::setNumericalScheme(SpaceScheme spaceScheme, TimeScheme timeScheme){
- _timeScheme = timeScheme;
- _spaceScheme = spaceScheme;
-}
-
void ProblemCoreFlows::setInitialField(const Field &VV)
{
void ProblemCoreFlows::setInitialFieldConstant(string fileName, const vector<double> Vconstant)
{
Mesh M(fileName);
- Field VV("Primitive", CELLS, M, Vconstant.size());
+ Field VV("SOLVERLAB results", CELLS, M, Vconstant.size());
for (int j = 0; j < M.getNumberOfCells(); j++) {
for (int i=0; i< VV.getNumberOfComponents(); i++)
}
void ProblemCoreFlows:: setInitialFieldConstant(const Mesh& M, const Vector Vconstant)
{
- Field VV("Primitive", CELLS, M, Vconstant.getNumberOfRows());
+ Field VV("SOLVERLAB results", CELLS, M, Vconstant.getNumberOfRows());
for (int j = 0; j < M.getNumberOfCells(); j++) {
for (int i=0; i< VV.getNumberOfComponents(); i++)
}
void ProblemCoreFlows:: setInitialFieldConstant(const Mesh& M, const vector<double> Vconstant)
{
- Field VV("Primitive", CELLS, M, Vconstant.size());
+ Field VV("SOLVERLAB results", CELLS, M, Vconstant.size());
for (int j = 0; j < M.getNumberOfCells(); j++) {
for (int i=0; i< VV.getNumberOfComponents(); i++)
_runLogFile->close();
throw CdmathException( "ProblemCoreFlows::setStepFunctionInitialField: Vectors VV_Left and VV_Right have different sizes");
}
- Field VV("Primitive", CELLS, M, VV_Left.getNumberOfRows());
+ Field VV("SOLVERLAB results", CELLS, M, VV_Left.getNumberOfRows());
double component_value;
{
return _precision;
}
-SpaceScheme ProblemCoreFlows::getSpaceScheme()
-{
- return _spaceScheme;
-}
-TimeScheme ProblemCoreFlows::getTimeScheme()
-{
- return _timeScheme;
-}
Mesh ProblemCoreFlows::getMesh()
{
return _mesh;
}
}
+void ProblemCoreFlows::setNewtonSolver(double precision, int iterations, nonLinearSolver solverName)
+{
+ _maxNewtonIts=iterations;
+ _precision_Newton=precision;
+ _nonLinearSolver=solverName;
+}
+
+
// Description:
// Cette methode lance une execution du ProblemCoreFlows
// Elle peut etre utilisee si le probleme n'est couple a aucun autre.
// (s'il n'a besoin d'aucun champ d'entree).
// Precondition: initialize
-// Seule la methode terminate peut etre appelee apres
+// Seule la methode terminate peut etre appelée apres
bool ProblemCoreFlows::run()
{
if(!_initializedMemory)
if (!ok) // The resolution failed, try with a new time interval.
{
- abortTimeStep();
if(_dt>_precision){
- cout << "Failed solving time step "<<_nbTimeStep<<", time = " << _time <<" _dt= "<<_dt<<", cfl= "<<_cfl<<", trying again with cfl/2"<< endl;
- *_runLogFile << "Failed solving time step "<<_nbTimeStep<<", time = " << _time <<" _dt= "<<_dt<<", cfl= "<<_cfl<<", trying again with cfl/2"<< endl;
- _dt*=0.5;
- _cfl*=0.5;
+ cout<<"ComputeTimeStep returned _dt="<<_dt<<endl;
+ cout << "Failed solving time step "<<_nbTimeStep<<", time = " << _time <<" _dt= "<<_dt<<", cfl= "<<_cfl<<", trying again with dt/2"<< endl;
+ *_runLogFile << "Failed solving time step "<<_nbTimeStep<<", time = " << _time <<" _dt= "<<_dt<<", cfl= "<<_cfl<<", trying again with dt/2"<< endl;
+ double dt=_dt/2;//We chose to divide the time step by 2
+ abortTimeStep();//Cancel the initTimeStep
+ _dt=dt;//new value of time step is previous time step divided by 2 (we do not call computeTimeStep
+ //_cfl*=0.5;//If we change the cfl, we must compute the new time step with computeTimeStep
+ //_dt=computeTimeStep(stop);
}
else{
cout << "Failed solving time step "<<_nbTimeStep<<", _time = " << _time<<" _dt= "<<_dt<<", cfl= "<<_cfl <<", stopping calculation"<< endl;
cout << endl;
}
void ProblemCoreFlows::setFileName(string fileName){
+ if( _nbTimeStep>0 && fileName!=_fileName)//This is a change of file name during a simulation
+ _restartWithNewFileName=true;
_fileName = fileName;
}
void ProblemCoreFlows::setFreqSave(int freqSave){
_freqSave = freqSave;
}
+
bool ProblemCoreFlows::solveTimeStep(){
_NEWTON_its=0;
bool converged=false, ok=true;
*/
delete _runLogFile;
}
+
+double
+ProblemCoreFlows::getConditionNumber(bool isSingular, double tol) const
+{
+ SparseMatrixPetsc A = SparseMatrixPetsc(_A);
+ return A.getConditionNumber( isSingular, tol);
+}
+std::vector< double >
+ProblemCoreFlows::getEigenvalues(int nev, EPSWhich which, double tol) const
+{
+ SparseMatrixPetsc A = SparseMatrixPetsc(_A);
+ return A.getEigenvalues( nev, which, tol);
+}
+std::vector< Vector >
+ProblemCoreFlows::getEigenvectors(int nev, EPSWhich which, double tol) const
+{
+ SparseMatrixPetsc A = SparseMatrixPetsc(_A);
+ return A.getEigenvectors( nev, which, tol);
+}
+Field
+ProblemCoreFlows::getEigenvectorsField(int nev, EPSWhich which, double tol) const
+{
+ SparseMatrixPetsc A = SparseMatrixPetsc(_A);
+ MEDCoupling::DataArrayDouble * d = A.getEigenvectorsDataArrayDouble( nev, which, tol);
+ Field my_eigenfield;
+
+ if(_FECalculation)
+ my_eigenfield = Field("Eigenvectors field", NODES, _mesh, nev);
+ else
+ my_eigenfield = Field("Eigenvectors field", CELLS, _mesh, nev);
+
+ my_eigenfield.setFieldByDataArrayDouble(d);
+
+ return my_eigenfield;
+}
+
+std::vector< double >
+ProblemCoreFlows::getSingularValues( int nsv, SVDWhich which, double tol) const
+{
+ SparseMatrixPetsc A = SparseMatrixPetsc(_A);
+ return A.getSingularValues( nsv, which, tol);
+}
+std::vector< Vector >
+ProblemCoreFlows::getSingularVectors(int nsv, SVDWhich which, double tol) const
+{
+ SparseMatrixPetsc A = SparseMatrixPetsc(_A);
+ return A.getSingularVectors( nsv, which, tol);
+}
+
+Field
+ProblemCoreFlows::getUnknownField() const
+{
+ if(!_initializedMemory)
+ {
+ _runLogFile->close();
+ throw CdmathException("ProblemCoreFlows::getUnknownField() call initialize() first");
+ }
+ return _VV;
+}
