From 13032eb7a54327ceec6e8880b39feb76a204aec9 Mon Sep 17 00:00:00 2001 From: landinjm Date: Fri, 3 Jan 2025 21:23:21 -0500 Subject: [PATCH] updating nonlinear printing info --- applications/allenCahn_conserved/customPDE.h | 27 +++---- include/core/matrixFreePDE.h | 2 +- src/core/invM.cc | 14 ---- src/core/solvers/solveIncrement.cc | 74 ++++++++++++-------- 4 files changed, 60 insertions(+), 57 deletions(-) diff --git a/applications/allenCahn_conserved/customPDE.h b/applications/allenCahn_conserved/customPDE.h index 4d0a5bae..0f169a69 100644 --- a/applications/allenCahn_conserved/customPDE.h +++ b/applications/allenCahn_conserved/customPDE.h @@ -177,13 +177,13 @@ customPDE::solveIncrement(bool skip_time_dependent) // parabolic and auxilary equations should also be here if (this->hasNonExplicitEquation) { - bool nonlinear_it_converged = false; - unsigned int nonlinear_it_index = 0; + bool nonlinear_iteration_converged = false; + unsigned int nonlinear_iteration_index = 0; - while (!nonlinear_it_converged) + while (!nonlinear_iteration_converged) { - nonlinear_it_converged = true; // Set to true here and will be set to false if - // any variable isn't converged + nonlinear_iteration_converged = true; // Set to true here and will be set to + // false if any variable isn't converged // Update residualSet for the non-explicitly updated variables // compute_nonexplicit_RHS() @@ -214,8 +214,8 @@ customPDE::solveIncrement(bool skip_time_dependent) this->pcout << buffer; } - nonlinear_it_converged = - this->updateImplicitSolution(fieldIndex, nonlinear_it_index); + nonlinear_iteration_converged = + this->updateImplicitSolution(fieldIndex, nonlinear_iteration_index); // Apply Boundary conditions this->applyBCs(fieldIndex); @@ -223,7 +223,7 @@ customPDE::solveIncrement(bool skip_time_dependent) else if (this->fields[fieldIndex].pdetype == AUXILIARY) { if (this->var_attributes.at(fieldIndex).is_nonlinear || - nonlinear_it_index == 0) + nonlinear_iteration_index == 0) { // If the equation for this field is nonlinear, save the // old solution @@ -281,18 +281,19 @@ customPDE::solveIncrement(bool skip_time_dependent) { this->pcout << "Relative difference between nonlinear " "iterations: " - << diff << " " << nonlinear_it_index << " " - << this->currentIncrement << std::endl; + << diff << " " << nonlinear_iteration_index + << " " << this->currentIncrement + << std::endl; } if (diff > MatrixFreePDE::userInputs .nonlinear_solver_parameters.getToleranceValue( fieldIndex) && - nonlinear_it_index < + nonlinear_iteration_index < MatrixFreePDE::userInputs .nonlinear_solver_parameters.getMaxIterations()) { - nonlinear_it_converged = false; + nonlinear_iteration_converged = false; } } else @@ -318,7 +319,7 @@ customPDE::solveIncrement(bool skip_time_dependent) } } - nonlinear_it_index++; + nonlinear_iteration_index++; } } diff --git a/include/core/matrixFreePDE.h b/include/core/matrixFreePDE.h index 934d4cd2..c94a4d23 100644 --- a/include/core/matrixFreePDE.h +++ b/include/core/matrixFreePDE.h @@ -288,7 +288,7 @@ class MatrixFreePDE : public Subscriptor /*Method to compute an implicit timestep*/ bool - updateImplicitSolution(unsigned int fieldIndex, unsigned int nonlinear_it_index); + updateImplicitSolution(unsigned int fieldIndex, unsigned int nonlinear_iteration_index); /*Method to apply boundary conditions*/ void diff --git a/src/core/invM.cc b/src/core/invM.cc index e5a499d0..2a4cc9d1 100644 --- a/src/core/invM.cc +++ b/src/core/invM.cc @@ -42,20 +42,6 @@ MatrixFreePDE::computeInvM() } } - // Check if invM has been found. If not, print an "error" message - if (!invMscalarFound) - { - pcout << "matrixFreePDE.h: no PARABOLIC scalar field... hence setting " - "parabolicScalarFieldIndex to 0 and marching ahead withn invM " - "computation\n"; - } - else if (!invMvectorFound) - { - pcout << "matrixFreePDE.h: no PARABOLIC vector field... hence setting " - "parabolicVectorFieldIndex to 0 and marching ahead withn invM " - "computation\n"; - } - // Initialize invM and clear its values matrixFreeObject.initialize_dof_vector(invMscalar, parabolicScalarFieldIndex); invMscalar = 0.0; diff --git a/src/core/solvers/solveIncrement.cc b/src/core/solvers/solveIncrement.cc index a8bb6b4f..51bf6b41 100644 --- a/src/core/solvers/solveIncrement.cc +++ b/src/core/solvers/solveIncrement.cc @@ -66,19 +66,14 @@ MatrixFreePDE::solveIncrement(bool skip_time_dependent) // parabolic and auxilary equations should also be here if (hasNonExplicitEquation) { - bool nonlinear_it_converged = false; - unsigned int nonlinear_it_index = 0; + bool nonlinear_iteration_converged = false; + unsigned int nonlinear_iteration_index = 0; - while (!nonlinear_it_converged) + while (!nonlinear_iteration_converged) { - nonlinear_it_converged = true; // Set to true here and will be set to false if - // any variable isn't converged + nonlinear_iteration_converged = true; // Update residualSet for the non-explicitly updated variables - // compute_nonexplicit_RHS() - // Ideally, I'd just do this for the non-explicit variables, but for - // now I'll do all of them this is a little redundant, but hopefully - // not too terrible computeNonexplicitRHS(); for (const auto &[fieldIndex, variable] : var_attributes) @@ -102,15 +97,15 @@ MatrixFreePDE::solveIncrement(bool skip_time_dependent) pcout << buffer; } - nonlinear_it_converged = - updateImplicitSolution(fieldIndex, nonlinear_it_index); + nonlinear_iteration_converged = + updateImplicitSolution(fieldIndex, nonlinear_iteration_index); // Apply Boundary conditions applyBCs(fieldIndex); } else if (fields[fieldIndex].pdetype == AUXILIARY) { - if (variable.is_nonlinear || nonlinear_it_index == 0) + if (variable.is_nonlinear || nonlinear_iteration_index == 0) { // If the equation for this field is nonlinear, save the old // solution @@ -164,19 +159,26 @@ MatrixFreePDE::solveIncrement(bool skip_time_dependent) } if (currentIncrement % userInputs.skip_print_steps == 0) { - pcout << "Relative difference between nonlinear " - "iterations: " - << diff << " " << nonlinear_it_index << " " - << currentIncrement << "\n"; + snprintf(buffer, + sizeof(buffer), + " field '%2s' [nonlinear solve] current " + "increment: %u, nonlinear " + "iteration: " + "%u, dU: %12.6e\n", + fields[fieldIndex].name.c_str(), + currentIncrement, + nonlinear_iteration_index, + diff); + pcout << buffer; } if (diff > userInputs.nonlinear_solver_parameters .getToleranceValue(fieldIndex) && - nonlinear_it_index < + nonlinear_iteration_index < userInputs.nonlinear_solver_parameters .getMaxIterations()) { - nonlinear_it_converged = false; + nonlinear_iteration_converged = false; } } else @@ -201,7 +203,7 @@ MatrixFreePDE::solveIncrement(bool skip_time_dependent) } } - nonlinear_it_index++; + nonlinear_iteration_index++; } } @@ -286,12 +288,12 @@ MatrixFreePDE::updateExplicitSolution(unsigned int fieldIndex) template bool MatrixFreePDE::updateImplicitSolution(unsigned int fieldIndex, - unsigned int nonlinear_it_index) + unsigned int nonlinear_iteration_index) { char buffer[200]; // Assume convergence criterion is met, unless otherwise proven later on. - bool nonlinear_it_converged = true; + bool nonlinear_iteration_converged = true; // Apply Dirichlet BC's. This clears the residual where we want to apply Dirichlet BCs, // otherwise the solver sees a positive residual @@ -469,18 +471,32 @@ MatrixFreePDE::updateImplicitSolution(unsigned int fieldIndex, } if (currentIncrement % userInputs.skip_print_steps == 0) { - pcout << "Relative difference between nonlinear " - "iterations: " - << diff << " " << nonlinear_it_index << " " << currentIncrement - << "\n"; + snprintf(buffer, + sizeof(buffer), + " field '%2s' [nonlinear solve] current increment: %u, nonlinear " + "iteration: " + "%u, dU: %12.6e\n", + fields[fieldIndex].name.c_str(), + currentIncrement, + nonlinear_iteration_index, + diff); + pcout << buffer; } if (diff > userInputs.nonlinear_solver_parameters.getToleranceValue(fieldIndex) && - nonlinear_it_index < + nonlinear_iteration_index < userInputs.nonlinear_solver_parameters.getMaxIterations()) { - nonlinear_it_converged = false; + nonlinear_iteration_converged = false; + } + else if (diff > + userInputs.nonlinear_solver_parameters.getToleranceValue(fieldIndex)) + { + pcout << "\nWarning: nonlinear solver did not converge as " + "per set tolerances. consider increasing the " + "maximum number of iterations or decreasing the " + "solver tolerance.\n"; } } else @@ -492,7 +508,7 @@ MatrixFreePDE::updateImplicitSolution(unsigned int fieldIndex, } else { - if (nonlinear_it_index == 0) + if (nonlinear_iteration_index == 0) { if (fields[fieldIndex].type == SCALAR) { @@ -532,5 +548,5 @@ MatrixFreePDE::updateImplicitSolution(unsigned int fieldIndex, } } - return nonlinear_it_converged; + return nonlinear_iteration_converged; } \ No newline at end of file