Commit 69bbb8b3 authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#987 AssertNumericalValues() provides now more precise informations to deliver...

#987 AssertNumericalValues() provides now more precise informations to deliver without ambiguity the vector for which the problem occurred. SolveLinear() and SolveNonLinear() interface has in consequence been enriched, as AssertNumericalValues() was call to check both the rhs and the solution.
parent 5cc91f5b
......@@ -183,7 +183,8 @@ namespace HappyHeart
*/
template<IsFactorized IsFactorizedT>
void SolveLinear(const NumberingSubset& row_numbering_subset,
const NumberingSubset& col_numbering_subset);
const NumberingSubset& col_numbering_subset,
const char* invoking_file, int invoking_line);
/*!
......@@ -201,7 +202,8 @@ namespace HappyHeart
template<IsFactorized IsFactorizedT>
void SolveLinear(const GlobalMatrix& matrix,
const GlobalVector& rhs,
GlobalVector& out);
GlobalVector& out,
const char* invoking_file, int invoking_line);
/*!
......@@ -217,7 +219,8 @@ namespace HappyHeart
* \copydetails doxygen_hide_varf_solve_numbering_subset_arg
*/
void SolveNonLinear(const NumberingSubset& row_numbering_subset,
const NumberingSubset& col_numbering_subset);
const NumberingSubset& col_numbering_subset,
const char* invoking_file, int invoking_line);
public:
......
......@@ -75,11 +75,13 @@ namespace HappyHeart
template<IsFactorized IsFactorizedT>
inline void VariationalFormulation<DerivedT, SolverIndexT>
::SolveLinear(const NumberingSubset& row_numbering_subset,
const NumberingSubset& col_numbering_subset)
const NumberingSubset& col_numbering_subset,
const char* invoking_file, int invoking_line)
{
return SolveLinear<IsFactorizedT>(GetSystemMatrix(row_numbering_subset, col_numbering_subset),
GetSystemRhs(col_numbering_subset),
GetNonCstSystemSolution(col_numbering_subset));
GetNonCstSystemSolution(col_numbering_subset),
invoking_file, invoking_line);
}
......@@ -92,10 +94,15 @@ namespace HappyHeart
void VariationalFormulation<DerivedT, SolverIndexT>
::SolveLinear(const GlobalMatrix& matrix,
const GlobalVector& rhs,
GlobalVector& out)
GlobalVector& out,
const char* invoking_file, int invoking_line)
{
assert(out.GetNumberingSubset() == rhs.GetNumberingSubset());
#ifndef NDEBUG
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(), "rhs", rhs, invoking_file, invoking_line);
#endif // NDEBUG
DoPrintMessage do_print_message =
(this->GetTimeManager().NtimeModified() % GetDisplayValue()) == 0
? DoPrintMessage::yes
......@@ -130,8 +137,7 @@ namespace HappyHeart
out.UpdateGhosts(__FILE__, __LINE__);
#ifndef NDEBUG
Wrappers::Petsc::AssertNumericValues(rhs, __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(out, __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(), "out", out, invoking_file, invoking_line);
#endif // NDEBUG
}
......@@ -144,7 +150,8 @@ namespace HappyHeart
>
void VariationalFormulation<DerivedT, SolverIndexT>
::SolveNonLinear(const NumberingSubset& row_numbering_subset,
const NumberingSubset& col_numbering_subset)
const NumberingSubset& col_numbering_subset,
const char* invoking_file, int invoking_line)
{
const auto& system_matrix = this->GetSystemMatrix(row_numbering_subset, col_numbering_subset);
const auto& system_rhs = this->GetSystemRhs(col_numbering_subset);
......@@ -153,7 +160,7 @@ namespace HappyHeart
GetNonCstSnes().SolveNonLinear(static_cast<void*>(this),
system_rhs,
system_matrix, system_matrix,
system_solution, __FILE__, __LINE__);
system_solution, invoking_file, invoking_line);
DebugPrintSolutionElasticWithOneUnknown(col_numbering_subset, 0);
}
......
......@@ -70,12 +70,14 @@ namespace HappyHeart
if (GetTimeManager().NtimeModified() == 1)
{
variational_formulation.ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
else
{
variational_formulation.ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
}
......
......@@ -62,12 +62,14 @@ namespace HappyHeart
if (this->GetTimeManager().NtimeModified() == 1)
{
formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
else
{
formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
formulation.template RemoveAverageExtracellularPotentialFromSolution();
......
......@@ -62,12 +62,14 @@ namespace HappyHeart
if (this->GetTimeManager().NtimeModified() == 1)
{
formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
else
{
formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
formulation.template RemoveAverageExtracellularPotentialFromSolution();
......
......@@ -79,7 +79,8 @@ namespace HappyHeart
// Solve the system.
TimeKeep::GetInstance().PrintTimeElapsed("Just before Newton in static case");
formulation.SolveNonLinear(displacement_numbering_subset, displacement_numbering_subset);
formulation.SolveNonLinear(displacement_numbering_subset, displacement_numbering_subset,
__FILE__, __LINE__);
TimeKeep::GetInstance().PrintTimeElapsed("Static Newton done.");
}
......@@ -105,7 +106,8 @@ namespace HappyHeart
const auto& displacement_numbering_subset = *displacement_numbering_subset_ptr;
auto& formulation = GetNonCstVariationalFormulation();
formulation.SolveNonLinear(displacement_numbering_subset, displacement_numbering_subset);
formulation.SolveNonLinear(displacement_numbering_subset, displacement_numbering_subset,
__FILE__, __LINE__);
}
......
......@@ -88,12 +88,12 @@ namespace HappyHeart
if (GetTimeManager().NtimeModified() == 1)
{
formulation.ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset, __FILE__, __LINE__);
}
else
{
formulation.ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset, numbering_subset);
formulation.SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset);
formulation.SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset, __FILE__, __LINE__);
}
}
......
......@@ -86,7 +86,8 @@ namespace HappyHeart
parent::GetNonCstSystemMatrix(numbering_subset, numbering_subset).Copy(GetStiffness(), __FILE__, __LINE__);
parent::template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
parent::template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
parent::template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
parent::DebugPrintSolutionElasticWithOneUnknown(numbering_subset, 10);
}
......
......@@ -236,7 +236,9 @@ namespace HappyHeart
auto& mixed_residual = GetNonCstMixedResidual(); // with velocity and pressure
#ifndef NDEBUG
Wrappers::Petsc::AssertNumericValues(GetSystemSolution(GetNumberingSubset()), __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(),
"system_sol",
GetSystemSolution(GetNumberingSubset()), __FILE__, __LINE__);
#endif // NDEBUG
Wrappers::Petsc::MatMult(GetResidualMatrix(),
......@@ -245,8 +247,11 @@ namespace HappyHeart
__FILE__, __LINE__);
#ifndef NDEBUG
Wrappers::Petsc::AssertNumericValues(GetResidualRhs(), __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(mixed_residual, __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(),
"residual_rhs",
GetResidualRhs(), __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(),
"mixed_residual",mixed_residual, __FILE__, __LINE__);
#endif // NDEBUG
Wrappers::Petsc::AXPY(-1., GetResidualRhs(), mixed_residual,
......@@ -265,7 +270,8 @@ namespace HappyHeart
__FILE__, __LINE__);
#ifndef NDEBUG
Wrappers::Petsc::AssertNumericValues(residual, __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(),
"residual",residual, __FILE__, __LINE__);
#endif // NDEBUG
residual.Scale(factor, __FILE__, __LINE__); // to comply with Freefem convention; probably not that good an idea.
......
......@@ -240,7 +240,8 @@ namespace HappyHeart
explicit_step_formulation.DefineFormulation();
explicit_step_formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset,
numbering_subset);
explicit_step_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
explicit_step_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
explicit_step_formulation.WriteSolution(time_manager, numbering_subset);
explicit_step_formulation.ComputeResidual();
......@@ -270,7 +271,8 @@ namespace HappyHeart
numbering_subset);
const auto& time_manager = this->GetTimeManager();
variational_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
variational_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
variational_formulation.WriteSolution(time_manager, numbering_subset);
......@@ -298,12 +300,17 @@ namespace HappyHeart
implicit_step_fluid_formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset,
numbering_subset);
implicit_step_fluid_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
implicit_step_fluid_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
// implicit_step_fluid_formulation.WriteSolution(time_manager, numbering_subset); #819 Do it carefully: files must not be overwritten (especially with the toggle is_differential).
#ifndef NDEBUG
Wrappers::Petsc::AssertNumericValues(implicit_step_fluid_formulation.GetSystemRhs(numbering_subset), __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(implicit_step_fluid_formulation.GetSystemSolution(numbering_subset), __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(),
"system_rhs",
implicit_step_fluid_formulation.GetSystemRhs(numbering_subset), __FILE__, __LINE__);
Wrappers::Petsc::AssertNumericValues(this->MpiHappyHeart(),
"system_sol",
implicit_step_fluid_formulation.GetSystemSolution(numbering_subset), __FILE__, __LINE__);
#endif // NDEBUG
if (is_differential == differential::no)
......
......@@ -67,7 +67,8 @@ namespace HappyHeart
"/Users/sebastien/Desktop/elastic_rhs.data",
__FILE__, __LINE__);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
......
......@@ -91,7 +91,8 @@ namespace HappyHeart
auto& variational_formulation = GetNonCstVariationalFormulation();
const auto& numbering_subset = variational_formulation.GetNumberingSubset();
variational_formulation.SolveNonLinear(numbering_subset, numbering_subset);
variational_formulation.SolveNonLinear(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
......
......@@ -81,12 +81,14 @@ namespace HappyHeart
if (GetTimeManager().NtimeModified() == 1)
{
variational_formulation.ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
else
{
variational_formulation.ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
}
......
......@@ -81,7 +81,8 @@ namespace HappyHeart
const NumberingSubset& numbering_subset = GetNumberingSubset();
ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
......
......@@ -77,7 +77,8 @@ namespace HappyHeart
// Solve the system.
TimeKeep::GetInstance().PrintTimeElapsed("Just before Newton in static case");
formulation.SolveNonLinear(numbering_subset, numbering_subset);
formulation.SolveNonLinear(numbering_subset, numbering_subset,
__FILE__, __LINE__);
TimeKeep::GetInstance().PrintTimeElapsed("Static Newton done.");
formulation.WriteSolution(this->GetTimeManager(), numbering_subset);
formulation.PrepareDynamicRuns();
......@@ -104,7 +105,8 @@ namespace HappyHeart
const auto& numbering_subset = *numbering_subset_ptr;
auto& formulation = this->GetNonCstVariationalFormulation();
formulation.SolveNonLinear(numbering_subset, numbering_subset);
formulation.SolveNonLinear(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
......
......@@ -297,7 +297,8 @@ namespace HappyHeart
explicit_step_formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset,
numbering_subset);
explicit_step_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
explicit_step_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
explicit_step_formulation.WriteSolution(time_manager, numbering_subset);
const auto& residual = explicit_step_formulation.ComputeResidual();
......@@ -368,7 +369,8 @@ namespace HappyHeart
implicit_step_fluid_formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset,
numbering_subset);
implicit_step_fluid_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
implicit_step_fluid_formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
// implicit_step_fluid_formulation.WriteSolution(time_manager, numbering_subset);
#ifndef NDEBUG
......
......@@ -67,7 +67,8 @@ namespace HappyHeart
"/Users/sebastien/Desktop/elastic_rhs.data",
__FILE__, __LINE__);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
......
......@@ -61,7 +61,8 @@ namespace HappyHeart
mpi, __FILE__, __LINE__);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
variational_formulation.SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
variational_formulation.DebugPrintSolution(numbering_subset);
variational_formulation.DebugPrintNorm(numbering_subset);
variational_formulation.WriteSolution(GetTimeManager(), numbering_subset);
......
......@@ -61,12 +61,14 @@ namespace HappyHeart
if (this->GetTimeManager().NtimeModified() == 1)
{
formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::nonlinear>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
else
{
formulation.template ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset);
formulation.template SolveLinear<IsFactorized::yes>(numbering_subset, numbering_subset,
__FILE__, __LINE__);
}
}
......
......@@ -65,7 +65,8 @@ namespace HappyHeart
// Solve the system.
TimeKeep::GetInstance().PrintTimeElapsed("Just before Newton in static case");
formulation.SolveNonLinear(numbering_subset, numbering_subset);
formulation.SolveNonLinear(numbering_subset, numbering_subset,
__FILE__, __LINE__);
TimeKeep::GetInstance().PrintTimeElapsed("Static Newton done.");
formulation.WriteSolution(GetTimeManager(), numbering_subset);
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment