From fef1637a6f9140cfacea566eccd29d25407ad81c Mon Sep 17 00:00:00 2001 From: "Joshua E. Hansel" Date: Mon, 14 Oct 2024 15:10:35 -0500 Subject: [PATCH] more cleanup --- framework/doc/content/bib/moose.bib | 12 ++ .../convergence/IterationCountConvergence.md | 22 +++ .../convergence/PostprocessorConvergence.md | 20 +++ .../ReferenceResidualConvergence.md | 11 +- .../source/convergence/ResidualConvergence.md | 70 +++++---- .../doc/content/syntax/Convergence/index.md | 136 ++++++++++++++++-- framework/include/convergence/Convergence.h | 8 +- .../ReferenceResidualConvergence.h | 1 - .../interfaces/ReferenceResidualInterface.h | 13 +- framework/include/problems/FEProblemBase.h | 4 + .../problems/ReferenceResidualProblem.h | 6 +- framework/src/convergence/Convergence.C | 3 +- .../ReferenceResidualConvergence.C | 60 -------- .../src/convergence/ResidualConvergence.C | 15 +- .../interfaces/ReferenceResidualInterface.C | 58 +++++++- .../src/problems/ReferenceResidualProblem.C | 2 +- framework/src/utils/PetscSupport.C | 2 +- .../AugmentedLagrangianContactConvergence.h | 6 +- .../AugmentedLagrangianContactProblem.h | 1 + .../AugmentedLagrangianContactConvergence.C | 19 ++- .../gold/reference_residual_out.e | Bin .../gold/residual_convergence_csv.csv | 0 .../gold/residual_default_csv.csv | 0 .../reference_residual.i | 0 .../residual_default.i | 0 .../tests | 1 + .../diffusion_convergence_default.i | 44 ------ .../gold/diffusion_convergence_out.e | Bin 41140 -> 0 bytes .../gold/residual_convergence_out.csv | 2 + ...n_convergence.i => residual_convergence.i} | 14 +- .../convergence/residual_convergence/tests | 23 +-- 31 files changed, 336 insertions(+), 217 deletions(-) create mode 100644 framework/doc/content/source/convergence/IterationCountConvergence.md create mode 100644 framework/doc/content/source/convergence/PostprocessorConvergence.md rename test/tests/convergence/{reference_residual => reference_residual_convergence}/gold/reference_residual_out.e (100%) rename test/tests/convergence/{reference_residual => reference_residual_convergence}/gold/residual_convergence_csv.csv (100%) rename test/tests/convergence/{reference_residual => reference_residual_convergence}/gold/residual_default_csv.csv (100%) rename test/tests/convergence/{reference_residual => reference_residual_convergence}/reference_residual.i (100%) rename test/tests/convergence/{reference_residual => reference_residual_convergence}/residual_default.i (100%) rename test/tests/convergence/{reference_residual => reference_residual_convergence}/tests (99%) delete mode 100644 test/tests/convergence/residual_convergence/diffusion_convergence_default.i delete mode 100644 test/tests/convergence/residual_convergence/gold/diffusion_convergence_out.e create mode 100644 test/tests/convergence/residual_convergence/gold/residual_convergence_out.csv rename test/tests/convergence/residual_convergence/{diffusion_convergence.i => residual_convergence.i} (67%) diff --git a/framework/doc/content/bib/moose.bib b/framework/doc/content/bib/moose.bib index 40a28ddf9884..40183a629b63 100644 --- a/framework/doc/content/bib/moose.bib +++ b/framework/doc/content/bib/moose.bib @@ -609,3 +609,15 @@ @article{liu2015finite year={2015}, publisher={Taylor \& Francis} } + +@article{rao2018stopping, + title = {A stopping criterion for the iterative solution of partial differential equations}, + journal = {Journal of Computational Physics}, + volume = {352}, + pages = {265-284}, + year = {2018}, + issn = {0021-9991}, + doi = {https://doi.org/10.1016/j.jcp.2017.09.033}, + url = {https://www.sciencedirect.com/science/article/pii/S0021999117306939}, + author = {Kaustubh Rao and Paul Malan and J. Blair Perot} +} diff --git a/framework/doc/content/source/convergence/IterationCountConvergence.md b/framework/doc/content/source/convergence/IterationCountConvergence.md new file mode 100644 index 000000000000..fba414e3d3b9 --- /dev/null +++ b/framework/doc/content/source/convergence/IterationCountConvergence.md @@ -0,0 +1,22 @@ +# IterationCountConvergence + +This [Convergence](Convergence/index.md) specifies: + +- $\ell_\text{max}$, the maximum number of iterations, + via [!param](/Convergence/IterationCountConvergence/max_iterations), and +- $\ell_\text{min}$, the minimum number of iterations, + via [!param](/Convergence/IterationCountConvergence/min_iterations). + +If [!param](/Convergence/IterationCountConvergence/converge_at_max_iterations) +is set to `true`, then the solve will converge when $\ell = \ell_\text{max}$ +instead of diverge. + +Other `Convergence` objects may inherit from this class and override +`checkConvergenceInner(iter)` instead of the usual `checkConvergence(iter)`, +to inherit the iteration bounds. An example is [/PostprocessorConvergence.md]. + +!syntax parameters /Convergence/IterationCountConvergence + +!syntax inputs /Convergence/IterationCountConvergence + +!syntax children /Convergence/IterationCountConvergence diff --git a/framework/doc/content/source/convergence/PostprocessorConvergence.md b/framework/doc/content/source/convergence/PostprocessorConvergence.md new file mode 100644 index 000000000000..76bd9f5248d8 --- /dev/null +++ b/framework/doc/content/source/convergence/PostprocessorConvergence.md @@ -0,0 +1,20 @@ +# PostprocessorConvergence + +This [Convergence](Convergence/index.md) derives from [/IterationCountConvergence.md] +and compares a [Postprocessor](Postprocessors/index.md) value $y$ to a tolerance $\tau$: + +!equation +y \leq \tau \,. + +For this to work as expected, the `execute_on` parameter of the post-processor +must include values that trigger execution before the desired check. For example, for nonlinear +convergence, the value `NONLINEAR_CONVERGENCE` should be used. + +Typically the post-processor used should attempt to approximate the error in a system, +such as [/AverageVariableChange.md]. + +!syntax parameters /Convergence/PostprocessorConvergence + +!syntax inputs /Convergence/PostprocessorConvergence + +!syntax children /Convergence/PostprocessorConvergence diff --git a/framework/doc/content/source/convergence/ReferenceResidualConvergence.md b/framework/doc/content/source/convergence/ReferenceResidualConvergence.md index adbc739416ce..c6d75d1e7ea1 100644 --- a/framework/doc/content/source/convergence/ReferenceResidualConvergence.md +++ b/framework/doc/content/source/convergence/ReferenceResidualConvergence.md @@ -1,13 +1,8 @@ # ReferenceResidualConvergence -This object uses a reference residual for computing relative convergence criteria. -The `ReferenceResidualConvergence` object facilitates the assesment of convergence properties of MOOSE applications and allows extensions to provide the user more control over the application behaviour. This convergence object replaces [ReferenceResidualProblem.md], and as such the capabilities are similar to `ReferenceResidualProblem`. - -## Description - -See [ReferenceResidualProblem.md]. Similar input as the one used in the `[Problem]` block can now be prescribed in the `[Convergence]` block as below - -!listing test/tests/convergence/reference_residual/reference_residual.i block=Convergence +This [Convergence](Convergence/index.md) is a [/ResidualConvergence.md] with a +customized reference residual for its relative convergence checks. See +[ReferenceResidualProblem.md] for more information. !syntax parameters /Convergence/ReferenceResidualConvergence diff --git a/framework/doc/content/source/convergence/ResidualConvergence.md b/framework/doc/content/source/convergence/ResidualConvergence.md index c10b67306c20..2eb5e518648e 100644 --- a/framework/doc/content/source/convergence/ResidualConvergence.md +++ b/framework/doc/content/source/convergence/ResidualConvergence.md @@ -1,44 +1,64 @@ # ResidualConvergence -By default, MOOSE checks convergence using relative and absolute criteria. Once the residual drops -below either an absolute tolerance or the residual divided by the initial residual for the current -timestep drops below a relative tolerance, the solution is considered converged. +This [Convergence](Convergence/index.md) uses a combination of criteria to +determine convergence, focused on the residual vector. -### +Solver convergence criteria parameters+ +Consider the system of algebraic equations: -Parameters for setting absolute convergence, relative convergence etc. are usually set in the `Executioner` block. However, now they can also prescribed in the `Convergence` block. +!equation +\mathbf{r}(\mathbf{u}) = \mathbf{0} \,. -We shall provide a few guidelines for setting convergence parameters. In the following presume we have a partial differential equation, given as -$\mathcal{P}(\mathbf u)=f$, and subsequent residual $\mathcal R(\mathbf u)=\mathcal P(\mathbf u)-f=0$. Convergence with respect to a tolerance $\tau$, implies $|\mathcal R (\mathbf u)|<\tau$, while divergence is encountered when the residual cannot decrease below the value of $\tau$. +This class reports convergence of the solution to this system if +any of the following conditions are true: -If we consider an iterative process to determine the solution $\mathbf u$, we have a set of intermediary solutions $\mathbf u_i$ required to verify the equation $\mathcal R(\mathbf u)\approx 0$, but do not meet the required tolerance $\tau$. +!equation +\|\mathbf{r}\|_2 < \tau_\text{abs} \,, -#### 1. Absolute and relative tolerances +!equation +\frac{\|\mathbf{r}\|_2}{\|\mathbf{r}_0\|_2} > \tau_\text{rel} \,, -- +Absolute Tolerance (`abs_tol`)+: This is a parameter that determines the threshold at which the residual norm of the solution is deemed small enough in absolute terms. For a system of equation this translates into $|\mathbf u_{i+1}-\mathbf u_i|<\tau$. - -- +Relative Tolerance (`rel_tol`)+: This parameter sets the threshold for the residual norm in relation to the norm of the right-hand side of the equation. For a system of equation this translates into $|\mathbf u_{i+1}-\mathbf u_i|/|\mathbf u_{i+1}|<\tau$. +!equation +\frac{\|\mathbf{\delta u}\|_2}{\|\mathbf{u}\|_2} < \tau_{\delta u,\text{rel}} \,. +This class reports divergence if any of the following conditions are true: -Considering that nonlinear systems are ultimately solved via linearization the user should append `l_` for linear systems, or `nl_` for nonlinear ones, on a per case base. +!equation +\|\mathbf{r}\|_2 = \text{NaN} \,, -#### 2. Choosing appropiate `abs_tol` and `rel_tol` +!equation +n_\text{evals} \geq n_\text{evals,max} \,, -- To start, apply the same values for both `abs_tol` and `rel_tol` in both preconditioned and non-preconditioned systems. This allows for a straightforward comparison under similar stopping criteria. -- +Preconditioning+ aims to improve the numerical properties of an iterative solver (specifically the condition number). This often results in the preconditioned system converging more quickly or requiring fewer iterations. However, preconditioning modifies the scale and the properties of the problem, which can affect the scale of the residuals. -- The user must make sure that `abs_tol` and `rel_tol` are set in a manner that mirrors the +scaling+ of the problem. If preconditioning significantly changes the scale (which is often the case), the user needs to adjust these tolerances to accommodate the changes. -- The user should +experiment+ with different settings of `abs_tol` and `rel_tol` to observe their effect on convergence and the quality of the solution. To make a fair comparison, one should ensure that the solver’s behavior (in terms of convergence and accuracy) is similar under both settings. +!equation +\|\mathbf{r}\|_2 > \tau_\text{div,abs} \,, -#### Other stopping criteria +!equation +\frac{\|\mathbf{r}\|_2}{\|\mathbf{r}_0\|_2} > \tau_\text{div,rel} \,, -Presuming the user gathered enough knowledge about the solver behaviour and the nature of the problem, addtional stopping requirements can be added, such as `nl_max_its` which instructs after how many iterations to abort the solver, or `nl_abs_step_tol`, which indicates what tolerance to be accepted at each solver step. +!equation +n_\text{ping} > n_\text{ping,max} \,, -## Example input syntax +where -!listing test/tests/convergence/residual_convergence/diffusion_convergence.i block=Convergence - -!listing test/tests/convergence/residual_convergence/diffusion_convergence.i block=Executioner -where the [!param](/Executioner/nonlinear_convergence) indicates the convergence type to be `ResidualConvergence` and additional parameters. Curently convergence specific parameters can be still specified in the Executioner block. +- $\|\cdot\|_2$ is the discrete $L_2$ norm, +- $\mathbf{r}_0$ is the initial (guess) residual vector, +- $\mathbf{\delta u} = \mathbf{u}^{(\ell)} - \mathbf{u}^{(\ell-1)}$ is the solution step vector, +- "NaN" is a not-a-number value, +- $\tau_\text{abs}$ is the absolute residual tolerance, provided by + [!param](/Convergence/ResidualConvergence/nl_abs_tol). +- $\tau_\text{rel}$ is the relative residual tolerance, provided by + [!param](/Convergence/ResidualConvergence/nl_rel_tol). +- $\tau_{\delta u,\text{abs}}$ is the relative step tolerance., provided by + [!param](/Convergence/ResidualConvergence/nl_rel_step_tol). +- $\tau_\text{div,abs}$ is the absolute residual divergence tolerance, provided by + [!param](/Convergence/ResidualConvergence/nl_abs_div_tol). +- $\tau_\text{div,rel}$ is the relative residual divergence tolerance, provided by + [!param](/Convergence/ResidualConvergence/nl_div_tol). +- $n_\text{evals}$ is the number of residual evaluations. +- $n_\text{evals,max}$ is the maximum number of residual evaluations, provided by + [!param](/Convergence/ResidualConvergence/nl_max_funcs). +- $n_\text{ping}$ is the number of ping-pong iterations. +- $n_\text{ping,max}$ is the maximum number of ping-pong iterations, provided by + [!param](/Convergence/ResidualConvergence/n_max_nonlinear_pingpong). !syntax parameters /Convergence/ResidualConvergence diff --git a/framework/doc/content/syntax/Convergence/index.md b/framework/doc/content/syntax/Convergence/index.md index 04182773c688..805cab77c556 100644 --- a/framework/doc/content/syntax/Convergence/index.md +++ b/framework/doc/content/syntax/Convergence/index.md @@ -1,24 +1,130 @@ -# Convergence system +# Convergence System -The Convergence system provides the infrastructure for creating MOOSE objects that interact with the -solvers to give the user more control over the application behaviour. +The Convergence system allows users to customize the stopping criteria for the +iteration in various solves: -## Description +- Nonlinear system solves +- Linear system solves (not yet implemented) +- Steady-state detection in [Transient.md] (not yet implemented) +- Fixed point solves (not yet implemented) -By default, MOOSE checks convergence using relative and absolute criteria. Once the residual drops -below either an absolute tolerance, or the residual divided by the initial residual for the current -time step drops below a relative tolerance, the solution is considered converged. This works well for -many problems, but there are some scenarios that are problematic, where the user may desire -interaction with the solver at runtime for better control or analysis. +Instead of supplying convergence-related parameters directly to the executioner, +the user creates `Convergence` objects whose names are then supplied to the +executioner, e.g., -Currently this object supports algebraic two types of convergence -[/ResidualConvergence.md] and -[ReferenceResidualConvergence.md]. +``` +[Convergence] + [my_convergence1] + type = MyCustomConvergenceClass + # some convergence parameters, like tolerances + [] +[] -### +Methods supported+ +[Executioner] + type = Steady + nonlinear_convergence = my_convergence1 +[] +``` -- [/ResidualConvergence.md] tracks the residual decay across iterations. This MOOSE object interfaces directly with the PETSc callback function at each iteration and allows the user to prescribe additional requirements and tests which can be performed at every iteration step. +Currently only the nonlinear solve convergence is supported, but others are planned +for the near future. -- [ReferenceResidualConvergence.md] uses a user-specified reference vector for convergence checks, instead of the initial residual. This is beneficial when solution variables have different scaling. Convergence is achieved when the $L_2$ norm of each solution variable’s residual is less than either the relative tolerance times the $L_2$ norm of the corresponding reference variable or the absolute tolerance. As this method computes relative convergence differently, the required nonlinear relative tolerance to achieve the same error can vary from the default approach in MOOSE. Users must ensure appropriate tolerances are used. +## Available Classes +The `Convergence` classes provided by MOOSE are the following: + +- [/ResidualConvergence.md]: The default convergence criteria, which includes + several convergence criteria combined together. +- [/ReferenceResidualConvergence.md]: Like `ResidualConvergence`, but + uses a custom norm to define the relative residual convergence criteria. + [/ReferenceResidualProblem.md] uses this as its convergence criteria. +- [/IterationCountConvergence.md]: Specifies minimum and maximum numbers of iterations. +- [/PostprocessorConvergence.md]: Compares the value of a + [Postprocessor](Postprocessors/index.md) to a tolerance. + +## Convergence Criteria Design Considerations + +Here we provide some considerations to make in designing convergence criteria +and choosing appropriate parameter values. +Consider a system of algebraic system of equations + +!equation +\mathbf{r}(\mathbf{u}) = \mathbf{0} \,, + +where $\mathbf{u}$ is the unknown solution vector, and $\mathbf{r}$ is the residual +function. To solve this system using an iterative method, we must decide on +criteria to stop the iteration. +In general, iteration for a solve should halt when the approximate solution $\tilde{\mathbf{u}}$ +has reached a satisfactory level of error $\mathbf{e} \equiv \tilde{\mathbf{u}} - \mathbf{u}$, +using a condition such as + +!equation id=error_criteria +\|\mathbf{e}\| \leq \tau_u \,, + +where $\|\cdot\|$ denotes some norm, and $\tau_u$ denotes some tolerance. +Unfortunately, since we do not know $\mathbf{u}$, the error $\mathbf{e}$ is +also unknown and thus may not be computed directly. Thus some approximation of +the condition [!eqref](error_criteria) must be made. This may entail some +approximation of the error $\mathbf{e}$ or some criteria which implies the +desired criteria. For example, a very common approach is to use a residual +criteria such as + +!equation +\|\mathbf{r}\| \leq \tau_{r,\text{abs}} \,. + +While it is true that $\|\mathbf{r}\| = 0$ implies $\|\mathbf{e}\| = 0$, a +zero-tolerance is impractical, and the translation between the tolerance +$\tau_u$ to the tolerance is $\tau_r$ is difficult. The "acceptable" absolute +residual tolerance is tricky to determine and is highly dependent on the +equations being solved. To attempt to circumvent this issue, relative residual +criteria have been used, dividing the residual norm by another value in an +attempt to normalize it. A common approach that has been used is to use the +initial residual vector $\mathbf{r}_0$ to normalize: + +!equation +\frac{\|\mathbf{r}\|}{\|\mathbf{r}_0\|} \leq \tau_{r,\text{rel}} \,, + +where $\tau_{r,\text{rel}}$ is the relative residual tolerance. The disadvantage +with this particular choice is that this is highly dependent on how good the +initial guess is: if the initial guess is very good, it will be nearly impossible +to converge to the tolerance, and if the initial guess is very bad, it will be +too easy to converge to the tolerance, resulting in an erroneous solution. + +Some other considerations are the following: + +- Consider round-off error: if error ever reaches values around round-off error, + the solve should definitely be considered converged, as iterating further + provides no benefit. +- Consider the other sources of error in the model that produced the system of + algebraic equations that you're solving. For example, if solving a system of + partial differential equations, consider the model error and the discretization + error; it is not beneficial to require algebraic error less than the other + sources of error. +- Since each convergence criteria typically has some weak point where they break + down, it is usually advisable to use a combination of criteria. + +For more information on convergence criteria, see [!cite](rao2018stopping) for +example. + +!alert tip title=Create your own convergence action +The `Convergence` system provides a lot of flexibility by providing several +pieces that can be combined together to create a desired set of convergence +criteria. Since this may involve a large number of objects (including objects +from other systems), it may be beneficial to create an [Action](/Action.md) +to create more compact and convenient syntax for your application. + +## Implementing a New Convergence Class + +`Convergence` objects are responsible for overriding the virtual method + +``` +MooseConvergenceStatus checkConvergence(unsigned int iter) +``` + +The returned type `MooseConvergenceStatus` is one of the following values: + +- `CONVERGED`: The system has converged. +- `DIVERGED`: The system has diverged. +- `ITERATING`: The system has neither converged nor diverged and thus will + continue to iterate. diff --git a/framework/include/convergence/Convergence.h b/framework/include/convergence/Convergence.h index c936469f7540..075cfd16d34c 100644 --- a/framework/include/convergence/Convergence.h +++ b/framework/include/convergence/Convergence.h @@ -14,6 +14,8 @@ #include "PostprocessorInterface.h" #include "PerfGraphInterface.h" +class FEProblemBase; + /** * Base class for convergence criteria. */ @@ -47,5 +49,9 @@ class Convergence : public MooseObject, virtual MooseConvergenceStatus checkConvergence(unsigned int iter) = 0; protected: - PerfID _perf_check_convergence; + /// FE problem + FEProblemBase & _fe_problem_base; + + /// Performance ID for \c checkConvergence + PerfID _perfid_check_convergence; }; diff --git a/framework/include/convergence/ReferenceResidualConvergence.h b/framework/include/convergence/ReferenceResidualConvergence.h index 9b8c6fd79120..a82d1be11df8 100644 --- a/framework/include/convergence/ReferenceResidualConvergence.h +++ b/framework/include/convergence/ReferenceResidualConvergence.h @@ -28,7 +28,6 @@ class ReferenceResidualConvergence : public ResidualConvergence, public Referenc { public: static InputParameters validParams(); - static InputParameters validCommonReferenceResidualProblemParams(); ReferenceResidualConvergence(const InputParameters & parameters); diff --git a/framework/include/interfaces/ReferenceResidualInterface.h b/framework/include/interfaces/ReferenceResidualInterface.h index 01bdc461dc7a..71761ca698fe 100644 --- a/framework/include/interfaces/ReferenceResidualInterface.h +++ b/framework/include/interfaces/ReferenceResidualInterface.h @@ -16,20 +16,15 @@ class InputParameters; class MooseObject; /** - * The ReferenceResidualInterface class is designed to provide an interface for - * the ReferenceConvergence class, and allow access to parameters shared by other classes. + * Interface class shared between ReferenceResidualProblem and ReferenceResidualConvergence. */ - -InputParameters validParams(); - class ReferenceResidualInterface { public: - ReferenceResidualInterface(const MooseObject * moose_object); - virtual ~ReferenceResidualInterface(); - static InputParameters validParams(); + ReferenceResidualInterface(const MooseObject * moose_object); + /** * Add a set of variables that need to be grouped together. For use in * actions that create variables. This is templated for backwards compatibility to allow passing @@ -40,8 +35,6 @@ class ReferenceResidualInterface void addGroupVariables(const std::set & group_vars); protected: - const InputParameters & _fi_params; - /// Name of variables that are grouped together to check convergence std::vector> _group_variables; diff --git a/framework/include/problems/FEProblemBase.h b/framework/include/problems/FEProblemBase.h index 1b2f3421431e..e97b087f04eb 100644 --- a/framework/include/problems/FEProblemBase.h +++ b/framework/include/problems/FEProblemBase.h @@ -2205,7 +2205,11 @@ class FEProblemBase : public SubProblem, public Restartable _set_nonlinear_convergence_name = true; } + /** + * Gets the nonlinear convergence object name if there is one + */ ConvergenceName getNonlinearConvergenceName() const { return _nonlinear_convergence_name; } + /** * Setter for whether we're computing the scaling jacobian */ diff --git a/framework/include/problems/ReferenceResidualProblem.h b/framework/include/problems/ReferenceResidualProblem.h index bd179f74a0dc..4de0cd51921e 100644 --- a/framework/include/problems/ReferenceResidualProblem.h +++ b/framework/include/problems/ReferenceResidualProblem.h @@ -13,8 +13,8 @@ #include "ReferenceResidualInterface.h" /** - * FEProblemBase derived class to enable convergence checking relative to a user-specified - * postprocessor + * Problem that checks for convergence relative to a user-supplied reference quantity + * rather than the initial residual. */ class ReferenceResidualProblem : public FEProblem, public ReferenceResidualInterface { @@ -24,6 +24,4 @@ class ReferenceResidualProblem : public FEProblem, public ReferenceResidualInter ReferenceResidualProblem(const InputParameters & params); virtual void addDefaultConvergence() override; - -protected: }; diff --git a/framework/src/convergence/Convergence.C b/framework/src/convergence/Convergence.C index 20fd669e3187..8c1f39c771a4 100644 --- a/framework/src/convergence/Convergence.C +++ b/framework/src/convergence/Convergence.C @@ -27,6 +27,7 @@ Convergence::Convergence(const InputParameters & parameters) SetupInterface(this), PostprocessorInterface(this), PerfGraphInterface(this), - _perf_check_convergence(registerTimedSection("checkConvergence", 5, "Checking Convergence")) + _fe_problem_base(*parameters.getCheckedPointerParam("_fe_problem_base")), + _perfid_check_convergence(registerTimedSection("checkConvergence", 5, "Checking Convergence")) { } diff --git a/framework/src/convergence/ReferenceResidualConvergence.C b/framework/src/convergence/ReferenceResidualConvergence.C index 8b700749718c..6731db3ab493 100644 --- a/framework/src/convergence/ReferenceResidualConvergence.C +++ b/framework/src/convergence/ReferenceResidualConvergence.C @@ -36,7 +36,6 @@ InputParameters ReferenceResidualConvergence::validParams() { InputParameters params = ResidualConvergence::validParams(); - params += validCommonReferenceResidualProblemParams(); params += ReferenceResidualInterface::validParams(); params.addClassDescription( @@ -46,65 +45,6 @@ ReferenceResidualConvergence::validParams() return params; } -InputParameters -ReferenceResidualConvergence::validCommonReferenceResidualProblemParams() -{ - InputParameters params = emptyInputParameters(); - params.addParam>( - "solution_variables", "Set of solution variables to be checked for relative convergence"); - params.addParam>( - "reference_residual_variables", - "Set of variables that provide reference residuals for relative convergence check"); - params.addParam("reference_vector", "The tag name of the reference residual vector."); - params.addParam("acceptable_multiplier", - 1.0, - "Multiplier applied to relative tolerance for acceptable limit"); - params.addParam("acceptable_iterations", - 0, - "Iterations after which convergence to acceptable limits is accepted"); - /*params.addParam>>( - "group_variables", - "Name of variables that are grouped together to check convergence. (Multiple groups can be " - "provided, separated by semicolon)");*/ - params.addParam>( - "converge_on", - {}, - "If supplied, use only these variables in the individual variable convergence check"); - MooseEnum Lnorm("global_L2 local_L2 global_Linf local_Linf", "global_L2"); - params.addParam( - "normalization_type", - Lnorm, - "The normalization type used to compare the reference and actual residuals."); - Lnorm.addDocumentation("global_L2", - "Compare the L2 norm of the residual vector to the L2 norm of the " - "absolute reference vector to determine relative convergence"); - Lnorm.addDocumentation( - "local_L2", - "Compute the L2 norm of the residual vector divided componentwise by the absolute reference " - "vector to the L2 norm of the absolute reference vector to determine relative convergence"); - Lnorm.addDocumentation( - "global_Linf", - "Compare the L-infinity norm of the residual vector to the L-infinity norm of the " - "absolute reference vector to determine relative convergence"); - Lnorm.addDocumentation("local_Linf", - "Compute the L-infinity norm of the residual vector divided componentwise " - "by the absolute reference " - "vector to the L-infinity norm of the absolute reference vector to " - "determine relative convergence"); - - MooseEnum zero_ref_res("zero_tolerance relative_tolerance", "relative_tolerance"); - params.addParam("zero_reference_residual_treatment", - zero_ref_res, - "Determine behavior if a reference residual value of zero is present " - "for a particular variable."); - zero_ref_res.addDocumentation("zero_tolerance", - "Solve is treated as converged if the residual is zero"); - zero_ref_res.addDocumentation( - "relative_tolerance", - "Solve is treated as converged if the residual is below the relative tolerance"); - return params; -} - ReferenceResidualConvergence::ReferenceResidualConvergence(const InputParameters & parameters) : ResidualConvergence(parameters), ReferenceResidualInterface(this), diff --git a/framework/src/convergence/ResidualConvergence.C b/framework/src/convergence/ResidualConvergence.C index c00223d5ce1f..039b1bca4b35 100644 --- a/framework/src/convergence/ResidualConvergence.C +++ b/framework/src/convergence/ResidualConvergence.C @@ -43,9 +43,11 @@ ResidualConvergence::ResidualConvergence(const InputParameters & parameters) _fe_problem(*getCheckedPointerParam("_fe_problem_base")) { EquationSystems & es = _fe_problem.es(); + _l_tol = getSharedESParam("l_tol", "linear solver tolerance", es); _l_abs_tol = getSharedESParam("l_abs_tol", "linear solver absolute tolerance", es); _l_max_its = getSharedESParam("l_max_its", "linear solver maximum iterations", es); + _nl_max_its = getSharedESParam("nl_max_its", "nonlinear solver maximum iterations", es); _nl_max_funcs = getSharedESParam( @@ -55,6 +57,10 @@ ResidualConvergence::ResidualConvergence(const InputParameters & parameters) _nl_rel_tol = getSharedESParam("nl_rel_tol", "nonlinear solver relative residual tolerance", es); _divtol = getSharedESParam("nl_div_tol", "nonlinear solver divergence tolerance", es); + + // These two are not used at all in the check below; it appears a relative step tolerance is + // obtained from PETSc, which may come from 'nl_rel_step_tol' in an indirect fashion. + // 'nl_abs_step_tol' is used only in libMesh's own nonlinear solver. _nl_abs_step_tol = getSharedESParam("nl_abs_step_tol", "nonlinear solver absolute step tolerance", es); _nl_rel_step_tol = @@ -66,9 +72,7 @@ ResidualConvergence::ResidualConvergence(const InputParameters & parameters) _nl_forced_its = getParam("nl_forced_its"); } else - { _nl_forced_its = _fe_problem.getNonlinearForcedIterations(); - } if (isParamSetByUser("nl_abs_div_tol")) _fe_problem.setNonlinearAbsoluteDivergenceTolerance(getParam("nl_abs_div_tol")); @@ -81,9 +85,7 @@ ResidualConvergence::ResidualConvergence(const InputParameters & parameters) _n_max_nl_pingpong = getParam("n_max_nonlinear_pingpong"); } else - { _n_max_nl_pingpong = _fe_problem.getMaxNLPingPong(); - } } bool @@ -106,7 +108,7 @@ ResidualConvergence::checkRelativeConvergence(const PetscInt /*it*/, Convergence::MooseConvergenceStatus ResidualConvergence::checkConvergence(unsigned int iter) { - TIME_SECTION(_perf_check_convergence); + TIME_SECTION(_perfid_check_convergence); NonlinearSystemBase & system = _fe_problem.currentNonlinearSystem(); MooseConvergenceStatus status = MooseConvergenceStatus::ITERATING; @@ -116,12 +118,15 @@ ResidualConvergence::checkConvergence(unsigned int iter) SNES snes = system.getSNES(); + // ||u|| ierr = SNESGetSolutionNorm(snes, &xnorm_petsc); CHKERRABORT(_fe_problem.comm().get(), ierr); + // ||r|| ierr = SNESGetFunctionNorm(snes, &fnorm_petsc); CHKERRABORT(_fe_problem.comm().get(), ierr); + // ||du|| ierr = SNESGetUpdateNorm(snes, &snorm_petsc); CHKERRABORT(_fe_problem.comm().get(), ierr); diff --git a/framework/src/interfaces/ReferenceResidualInterface.C b/framework/src/interfaces/ReferenceResidualInterface.C index ff85553a3af7..b3ec9f739e75 100644 --- a/framework/src/interfaces/ReferenceResidualInterface.C +++ b/framework/src/interfaces/ReferenceResidualInterface.C @@ -14,22 +14,70 @@ InputParameters ReferenceResidualInterface::validParams() { InputParameters params = emptyInputParameters(); + + params.addParam>( + "solution_variables", "Set of solution variables to be checked for relative convergence"); + params.addParam>( + "reference_residual_variables", + "Set of variables that provide reference residuals for relative convergence check"); + params.addParam("reference_vector", "The tag name of the reference residual vector."); + params.addParam("acceptable_multiplier", + 1.0, + "Multiplier applied to relative tolerance for acceptable limit"); + params.addParam("acceptable_iterations", + 0, + "Iterations after which convergence to acceptable limits is accepted"); params.addParam>>( "group_variables", "Name of variables that are grouped together to check convergence. (Multiple groups can be " "provided, separated by semicolon)"); + params.addParam>( + "converge_on", + {}, + "If supplied, use only these variables in the individual variable convergence check"); + MooseEnum Lnorm("global_L2 local_L2 global_Linf local_Linf", "global_L2"); + params.addParam( + "normalization_type", + Lnorm, + "The normalization type used to compare the reference and actual residuals."); + Lnorm.addDocumentation("global_L2", + "Compare the L2 norm of the residual vector to the L2 norm of the " + "absolute reference vector to determine relative convergence"); + Lnorm.addDocumentation( + "local_L2", + "Compute the L2 norm of the residual vector divided componentwise by the absolute reference " + "vector to the L2 norm of the absolute reference vector to determine relative convergence"); + Lnorm.addDocumentation( + "global_Linf", + "Compare the L-infinity norm of the residual vector to the L-infinity norm of the " + "absolute reference vector to determine relative convergence"); + Lnorm.addDocumentation("local_Linf", + "Compute the L-infinity norm of the residual vector divided componentwise " + "by the absolute reference " + "vector to the L-infinity norm of the absolute reference vector to " + "determine relative convergence"); + + MooseEnum zero_ref_res("zero_tolerance relative_tolerance", "relative_tolerance"); + params.addParam("zero_reference_residual_treatment", + zero_ref_res, + "Determine behavior if a reference residual value of zero is present " + "for a particular variable."); + zero_ref_res.addDocumentation("zero_tolerance", + "Solve is treated as converged if the residual is zero"); + zero_ref_res.addDocumentation( + "relative_tolerance", + "Solve is treated as converged if the residual is below the relative tolerance"); + return params; } ReferenceResidualInterface::ReferenceResidualInterface(const MooseObject * moose_object) - : _fi_params(moose_object->parameters()), _use_group_variables(false) + : _use_group_variables(false) { - if (_fi_params.isParamValid("group_variables")) + if (moose_object->isParamValid("group_variables")) { _group_variables = - _fi_params.get>>("group_variables"); + moose_object->getParam>>("group_variables"); _use_group_variables = true; } } - -ReferenceResidualInterface::~ReferenceResidualInterface() {} diff --git a/framework/src/problems/ReferenceResidualProblem.C b/framework/src/problems/ReferenceResidualProblem.C index 4ca5e530d368..0521592046d2 100644 --- a/framework/src/problems/ReferenceResidualProblem.C +++ b/framework/src/problems/ReferenceResidualProblem.C @@ -16,8 +16,8 @@ InputParameters ReferenceResidualProblem::validParams() { InputParameters params = FEProblem::validParams(); - params += ReferenceResidualConvergence::validCommonReferenceResidualProblemParams(); params += ReferenceResidualInterface::validParams(); + params.addClassDescription("Problem that checks for convergence relative to " "a user-supplied reference quantity rather than " "the initial residual"); diff --git a/framework/src/utils/PetscSupport.C b/framework/src/utils/PetscSupport.C index c6ad00731907..d58e6c3f1ec1 100644 --- a/framework/src/utils/PetscSupport.C +++ b/framework/src/utils/PetscSupport.C @@ -284,7 +284,7 @@ petscSetupOutput(CommandLine * cmd_line) } PetscErrorCode -petscNonlinearConverged(SNES snes, +petscNonlinearConverged(SNES /*snes*/, PetscInt it, PetscReal /*xnorm*/, PetscReal /*snorm*/, diff --git a/modules/contact/include/convergence/AugmentedLagrangianContactConvergence.h b/modules/contact/include/convergence/AugmentedLagrangianContactConvergence.h index bda904db7148..1ef811568d30 100644 --- a/modules/contact/include/convergence/AugmentedLagrangianContactConvergence.h +++ b/modules/contact/include/convergence/AugmentedLagrangianContactConvergence.h @@ -9,10 +9,6 @@ #pragma once -#include "ReferenceResidualProblem.h" -#include "FEProblem.h" -#include "NodeFaceConstraint.h" -#include "MechanicalContactConstraint.h" #include "AugmentedLagrangianContactProblemInterface.h" #include "ReferenceResidualConvergence.h" #include "ResidualConvergence.h" @@ -33,9 +29,9 @@ class AugmentedLagrangianContactConvergence : public T, checkConvergence(unsigned int iter) override; protected: - FEProblemBase & _fe_problem; using AugmentedLagrangianContactProblemInterface::_lagrangian_iteration_number; using AugmentedLagrangianContactProblemInterface::_maximum_number_lagrangian_iterations; + using Convergence::_fe_problem_base; }; typedef AugmentedLagrangianContactConvergence diff --git a/modules/contact/include/problems/AugmentedLagrangianContactProblem.h b/modules/contact/include/problems/AugmentedLagrangianContactProblem.h index 1cbd72f7aee9..bddf7a4502b3 100644 --- a/modules/contact/include/problems/AugmentedLagrangianContactProblem.h +++ b/modules/contact/include/problems/AugmentedLagrangianContactProblem.h @@ -40,6 +40,7 @@ class AugmentedLagrangianContactProblemTempl : public T, using AugmentedLagrangianContactProblemInterface::_lagrangian_iteration_number; using AugmentedLagrangianContactProblemInterface::_maximum_number_lagrangian_iterations; }; + typedef AugmentedLagrangianContactProblemTempl AugmentedLagrangianContactProblem; typedef AugmentedLagrangianContactProblemTempl AugmentedLagrangianContactFEProblem; diff --git a/modules/contact/src/convergence/AugmentedLagrangianContactConvergence.C b/modules/contact/src/convergence/AugmentedLagrangianContactConvergence.C index cdecb578862f..2fe6444c67b2 100644 --- a/modules/contact/src/convergence/AugmentedLagrangianContactConvergence.C +++ b/modules/contact/src/convergence/AugmentedLagrangianContactConvergence.C @@ -40,16 +40,14 @@ AugmentedLagrangianContactConvergence::validParams() { InputParameters params = T::validParams(); params += AugmentedLagrangianContactProblemInterface::validParams(); - params.addClassDescription("Manages convergence for augmented Lagrange contact"); + params.addClassDescription("Convergence for augmented Lagrangian contact"); return params; } template AugmentedLagrangianContactConvergence::AugmentedLagrangianContactConvergence( const InputParameters & params) - : T(params), - AugmentedLagrangianContactProblemInterface(params), - _fe_problem(*params.getCheckedPointerParam("_fe_problem_base")) + : T(params), AugmentedLagrangianContactProblemInterface(params) { } @@ -59,7 +57,7 @@ AugmentedLagrangianContactConvergence::checkConvergence(unsigned int iter) { auto reason = T::checkConvergence(iter); - auto aug_contact = dynamic_cast(&_fe_problem); + auto aug_contact = dynamic_cast(&_fe_problem_base); _lagrangian_iteration_number = aug_contact->getLagrangianIterationNumber(); bool repeat_augmented_lagrange_step = false; @@ -68,15 +66,15 @@ AugmentedLagrangianContactConvergence::checkConvergence(unsigned int iter) { if (_lagrangian_iteration_number < _maximum_number_lagrangian_iterations) { - auto & nonlinear_sys = _fe_problem.currentNonlinearSystem(); + auto & nonlinear_sys = _fe_problem_base.currentNonlinearSystem(); nonlinear_sys.update(); // Get the penetration locator from the displaced mesh if it exist, otherwise get // it from the undisplaced mesh. - const auto displaced_problem = _fe_problem.getDisplacedProblem(); - const auto & penetration_locators = - (displaced_problem ? displaced_problem->geomSearchData() : _fe_problem.geomSearchData()) - ._penetration_locators; + const auto displaced_problem = _fe_problem_base.getDisplacedProblem(); + const auto & penetration_locators = (displaced_problem ? displaced_problem->geomSearchData() + : _fe_problem_base.geomSearchData()) + ._penetration_locators; // loop over contact pairs (penetration locators) const ConstraintWarehouse & constraints = nonlinear_sys.getConstraintWarehouse(); @@ -135,7 +133,6 @@ AugmentedLagrangianContactConvergence::checkConvergence(unsigned int iter) if (repeat_augmented_lagrange_step) { _lagrangian_iteration_number++; - // setLagrangianIterationNumber(_lagrangian_iteration_number); Moose::out << "Augmented Lagrangian contact repeat " << _lagrangian_iteration_number << '\n'; diff --git a/test/tests/convergence/reference_residual/gold/reference_residual_out.e b/test/tests/convergence/reference_residual_convergence/gold/reference_residual_out.e similarity index 100% rename from test/tests/convergence/reference_residual/gold/reference_residual_out.e rename to test/tests/convergence/reference_residual_convergence/gold/reference_residual_out.e diff --git a/test/tests/convergence/reference_residual/gold/residual_convergence_csv.csv b/test/tests/convergence/reference_residual_convergence/gold/residual_convergence_csv.csv similarity index 100% rename from test/tests/convergence/reference_residual/gold/residual_convergence_csv.csv rename to test/tests/convergence/reference_residual_convergence/gold/residual_convergence_csv.csv diff --git a/test/tests/convergence/reference_residual/gold/residual_default_csv.csv b/test/tests/convergence/reference_residual_convergence/gold/residual_default_csv.csv similarity index 100% rename from test/tests/convergence/reference_residual/gold/residual_default_csv.csv rename to test/tests/convergence/reference_residual_convergence/gold/residual_default_csv.csv diff --git a/test/tests/convergence/reference_residual/reference_residual.i b/test/tests/convergence/reference_residual_convergence/reference_residual.i similarity index 100% rename from test/tests/convergence/reference_residual/reference_residual.i rename to test/tests/convergence/reference_residual_convergence/reference_residual.i diff --git a/test/tests/convergence/reference_residual/residual_default.i b/test/tests/convergence/reference_residual_convergence/residual_default.i similarity index 100% rename from test/tests/convergence/reference_residual/residual_default.i rename to test/tests/convergence/reference_residual_convergence/residual_default.i diff --git a/test/tests/convergence/reference_residual/tests b/test/tests/convergence/reference_residual_convergence/tests similarity index 99% rename from test/tests/convergence/reference_residual/tests rename to test/tests/convergence/reference_residual_convergence/tests index fde4445244db..2172ba7f01a7 100644 --- a/test/tests/convergence/reference_residual/tests +++ b/test/tests/convergence/reference_residual_convergence/tests @@ -1,6 +1,7 @@ [Tests] design = 'ReferenceResidualConvergence.md' issues = '#25931' + [reference_residual] type = 'Exodiff' input = 'reference_residual.i' diff --git a/test/tests/convergence/residual_convergence/diffusion_convergence_default.i b/test/tests/convergence/residual_convergence/diffusion_convergence_default.i deleted file mode 100644 index fd237781689f..000000000000 --- a/test/tests/convergence/residual_convergence/diffusion_convergence_default.i +++ /dev/null @@ -1,44 +0,0 @@ -[Mesh] - type = GeneratedMesh - dim = 2 - nx = 10 - ny = 10 -[] - -[Variables] - [u] - [] -[] - -[Kernels] - [diff] - type = Diffusion - variable = u - [] -[] - -[BCs] - [left] - type = DirichletBC - variable = u - boundary = left - value = 0 - [] - [right] - type = DirichletBC - variable = u - boundary = right - value = 1 - [] -[] - -[Executioner] - type = Steady - solve_type = 'PJFNK' - petsc_options_iname = '-pc_type' - petsc_options_value = 'hypre' -[] - -[Outputs] - exodus = true -[] diff --git a/test/tests/convergence/residual_convergence/gold/diffusion_convergence_out.e b/test/tests/convergence/residual_convergence/gold/diffusion_convergence_out.e deleted file mode 100644 index 936ae63210ff95afe64cac109c29e99a638fc5e8..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 41140 zcmeHQd#oJQdB5gie2sY%VhAJy2J8U7Hon(j8+QX@4|qblYRp{i2Ss((~ftCm)gP*ka*A}VVBNU7A8HW31qR%!cB|7gC}bFuvMnjH@f6vF zl3lIe;|sWvt07Ow2-htmT(^vHow5Zdq>OPQ;?46b0WM3tW_r{-AL%jjBxnujD$vzV z!?PFMF!BN)mm^n#R!pt<;cvmeDnq`uQVX2e^Q(3v@ch_56NF{N&pBcceJRCgA=I@23AVnexVJP&c}c?EFN81^T z-6)PWl>#41S6#p4ZYq0~N|WKT6H+^9#$)O^#`~8*tA%ZznGX~o^&@Gw63i-jasFob zqi(US6s~zM<4rs2NrD^QgP&EnzYzuv+K#@L@lXB=nMXtzE}!My@>%`-K3tFbQP2CR zlXkrKdU-LrAO0%bpL`$hUFv=M`DbzOXw>^Y+KP6`_thX=h5M85&_;9`Y28-aGhVHB+^IFC zLGT{FpANicZ@3-TjMF5Rj{gC%RC12oD6kM;Nh}?=hFG@{Yc=Dpg(kzvaL^5ir*FTF zcxzbx+AN^z4vXsU^guvX{G8w;Be#pNAaD6@DW+``i1}m`q~ATtvfGZ&HdLPDJu9C) zHnm5TUc%N_?t~2Pi~Lz5?(clRjy(R)_kWAb$jIG0@v}!kG86tk5ojd#?qd#$D=O=4Fk`0-d;wC+(8BjC+l{)V<%~y-Ht3 zD{-#@X_v%(6GDx<)V&|@UTF`l_`Ozx!YA$n@_POC@b$eq9fo$pQ2I`Zw|1H49D+xmDuP$l-{| zIL;`0&Wn7+y4Mcf5=u5*{@jFTZU$)woW@akvEfbl#H6F2BrW56`Of5j$!GNcynkf! zjq`LbYhMT{e@n@=|C4dyCbj>r`d|C+$tNeZ|1)*J_TR4vH>v#()$`i_=xw?mPip@c z3QzlQ$or)vDA1zApcV?dx)W*uJihhV1Knd8zuv zAdUku#B*uJKBGhrvxP(JLGE#FZmvZutEa(5_{C_DI0GE3Hzm(H{IqV_^ zI7X$%rZ`^5ad?iObKH9sh~vc^x8)cs$LKhA$#F)GA6^gQcpt~`^!OR)7jFT59<&y8 zE9f@R7eE|)y930*WX|napgTc#f!2ZU2HgX?7qlL<0rX|ieW3e69M9Yc8Ut+tZ3c1J zk7Jq~tDFFB1w90M7_<$v9kc`V2#8~w1<+%lNzhKv6lfP{H)s!N8Z-mi3)%;o1w9UW z0<<4=0CW&^2y__qB4585`gU1%H8_MvS>+sVy%H+7;O+Gfa8584i8<_1FTUj6C&Vo}!FA^_T!rMxJ^wuTe&xdNAKn zMxMGo45DldL_KtVBu_ndC>eR`@raUULPk$Wyl`m5e-fJECOdsoR{Ak*97)m5e;~ zIHqLesmF08BTqe^QZn+??Szt1C+hYz=o!URx06anp1M7&WaO#ab4o^@x;?LCQ(|#MxCfz8RROSx>b~nJawxo8F}h9uVmz@o2O*t zsmCcLBmV;EG>9_t)T0KXj6C(IgD4|U-Fy&b)QP$UpoZe9+lxv@p1OrfMxMGwN=BZ# z#Y#q=dNh@cJoQ*mGV;{pjFOS39*aswp1Pe?vK8Zf0(2#4CFm;9D$posHRx*4 z8qhVMYeCn6t_R%!x)F2}=w{I8K(~NC4_XVl6?7Zu3!vLUcYwYK`Vz1vCzt0Br?51bP^>4YVDk^Z29S3!uk9 zlc1fTDbOy^ZqOdkG-w927qkyF3wj*%1ZY3#0O%ly`9j;iwsCFCv|nwfv`yNmwnJ@i z+O8b@R|L^!XeTQ`mw}c_58gKjG3}TxjB~~(Jr)|+@c7d42x{T3z zC4j5}dJze>@Zn3@txIti6d&uX z-J$EmZrLiHwGJFQG&gNA0rI;P1uq!T+8xyEj$gL+V_Bdz6;_*d*N-C$@5-U5HMXf9 z1d(fVA>PKI5nCHQYrsW$V;zHUmrrte(@bL7}xl%8jWsQ$***wNN%h1j3pFW%o$coU=?QVu4tj;=$<-zUnA}WtI zU2DcITGm6>gA+T(x9`}xO+?>P_|uDSsTn)Pnkz7|c(mRSz&(KD;{yRM^}tY>tGoo$8@zLRWj4AWVb-SoW|o38DYn2p5>$aV&_tb#RFE+27gB9lk+UgLn{ zJ5`reW~#H0Vf{D^YIf{Yqa>2z3)bwx<5T<965Ujc4xwdvK6P7gQ@|)HSR4GnF1r<{ zS&P}+Si%v|i#?}iSHhrfi?!8UxXd;6;uQsJ|J2dxIkx@=a%q5AG}l%)PmraC<+zQK z;LHJQ!5OT@F9uCiQz_0o4xJhxZQrfJFcw_<#WODFXIlwgu*?cYm`yXn{5Ccc1f|n0 zn6|3s8Y)1n8$^BB2GX2)9CBqSidb79YM8c0DpKT>A>*UxbPKFIjO{Qz#NoPRs!~|1kAJE5|w~?;Q+wsQ1KKN+{<85=Y%sazAUS6UMX06ywEGn z*W7s5?$lW86jB$b7zMRv>|zBrc2;0>Q7Oh21y!x)i9=ERJK@oJz6XcA*O8A1cK01r0Q{7@4ilcn}9iM81E;bUTsRvme4YN#8 zydsYsMlOaNvRuJ3>z@&)i4MIs?y`->=6nJ@51=R$5*Qt9KETk+%qs=UY-o$xWSF#6 zv|zyoymHg2^*0y&FKfNAO~b14H0m&YQ7?-{ibK+p5tXN2&Ri%-z~y?GvCwuGL%p1t zpsD7|ur1fiokI`H^>V7IB!x7a<$C#Jua`|nMNo}8msGuc5eyORZl8@Xq^drd5r@Fb zridGf^BRVsA`6G5n<9yLk)7&J;?vGqC(jL#mNp5pR|fgBNPHij|mIH>i)iN2n@aqHl{DG7)ZPi(I})SLuH=W#RK zK-%#p37DZB{{*0$HOrXx7Cz%+cwiA*-F+St5dwGakK&c{+~@Fg+q=9C8uQV2;rQ!1Um%ICSrS&@7b!h`H*LTTU~boIz$X zSDJpQXQWD|vOzklAZJ^4!7@lxi(j zgEJL!YF&>%+C>)?gG-yopfu*TqJg54hAjb=en_2;py$gePAL|r%i*LwCoIj|VwxSP zU32QivSXudM)kPpb4e$&UN>N^U3%w;LpevE)uIy>y+J(%4JBYhF1$z>jAV0bm{m>5 zrie);RrXL8#&)A2XW_b4G{|)#!8kS-ahM(JPTIM+Ws5$R)7hnYw{&_SBLrj{Q03wa zz)~9a?NsC`50gPnu^n@>=>V)%gV2lT>+Pf|STZh_%p&Rr%Sr-@W+IyJVZGI^Gf(9t zU`ngbym06M#DSY#9iR&XOvx_hImN3Y&Nx{u$EzZ?w`~^Vx&SxJ%}aD0g&i7p9&4gw z))jL}voY;FX0?SCri{b%C8}w*i6EoLsEYo!UTv;msby|yPp8N_0PTiTI*p}!Z8*kZ z11GwcQXWkgSr(0Xv4W!rn~gFSW+%bO7p#db+b6a^^zfGL1BDxal}cHhYuF7XuM5_h zM!GKP(j~gz!fXsWW-MS;=QgW9%2YP%fGfcC?l(0YxFsjzprv8EYVV9oVTVTba^LQ^ zh(p5;c;D`K&LC%Zzf~$5q_c>_l!Ceyk>(w61Ub9=oioTO1$9D04po`7ZQ4Ur=M$GE5gx4DO zwqgs`{8{vUmcnWM4?x*nOlippMa@ea1(?$CQh=uAp#sdx&*pZh0JE~62(nSu1Cp~S zI_`D_U9=#vTAX1v*H8eWqu6SGiKuJ+m@k47{Rn#@vQIg!J~)~X*!v?u%`2%^#{^rY z2pSYU$H!n|)s5%fFh^bnVC;mkT)|&-qKM2Cd1(`!_8PXk$l)dK{*YZy4aYb%&3cjVM)eA6 zK88FkQB11aeXE2sl4olMQX^A%#Ffs7Y9OToawBo5B!za;_8}dd2ZsKa|)XAAwQp*b(yZ2J&i9@+_A!Ai` z9>7m}qk-IK(9f8Y|Iz69F%=u_OH7j4QL}GqZ z+^`W~6EIpN{EU#Us7>o-#*=PJtV*6yyF@k5 z_Mu#sut#n8DjyC5%!$=ytxt72Q)r25K2Hl)Z^n^Fb+ipoapvXW(08db<7q8>whS@I zradXoFV&j(C=n`WY<|aJT6I2%+8rVDyetRPx?tK-X$~khf7@Qx)+_Twk$~M}Xe||U z0CERCwoqm$XYiDauo9H} zN>k769Og9Fuon=(i_k-7qDUO{)N>)|hwfFhCOy#@7>Wa}A4b+y^g5UdEtsWG9>Hol2V5dIB(g%R8ZS z0-ozsD{++ruwQREPkC5udbKi_j;ShuXTKy`5axzm!m%e|;9=3BF1Ju-$(2zwBXhH9 zS5mE7Nmf#=I?E_w*gnLuC-U%^L|U}iJQ@-(vk)s-C!Xf9e}g&>VI)>4wG?YH11{=o zy$^do^amkZ7#)mxkcYnRtLorfykV9$Ri}whe%6AbgKy!8ndYt_o6c6~7BCDu#%fD` zVmFCs9*6P+yd33BJ>S3_kaL(d*QxaACiO&VV-#l0(c^H9$m?vwCBkf8IQZG6HF0xj zm|a@u1k5h26<{k{4Z>ht3Dcmp8n;5Eo{xkv*y)+6WBUg^K}-P3_F%_^Bb5`*eBD>N zJ!?E9%7mfWNI*63gAN87Zq-CO6<}AG!x<#{#66Z3p2N}-ma+JLb&53>ES|JGxJv-0 zA0cL*^;#OXYPt?RF8~k=RxwZ3frg|p2Nb;k*>|wd%;Ur4olRQvrU)jLTY+sEOPev* za1rHK`}rQuT*JkB%27jYT=U}V?65UFJ8TWl4O_!=!`1L2m~0*P2Hk3WcI5tiv4t)7yIeM*>NDVcnS{n@)%b(JgPO}g*FTN6hTyT#qQ|-b~0Gcjn z6ia}_=CTeTePOO)DX4WXZ_{A(&8q-*9s``uo8fW+_CK-KbRn5ql+azlYlU;fVYHwh zPV>zH%R!~?WbC@a8_WS+aa?r9!2twZKYMWY=*(wO9IiWk*q%PPC-b~GbZ-}f zax*Ii^O%h|=(l%8nvPHNwoSD@>qpjrCv96tS$Z{gQLtS&gp?AR5G@Mrv>n*upKgRfTE9;mW zcEb%T7+rB1$wLY+Bo5nv*xcHCkhgg%cLCx+#X6ArodDqgXt4K#6V}~;BDLFy*!_`O zp=SxX^aQ3Yw({Zx5G0A%UQx#6g(D0!N05_%dM$M;k#kg|QbP;|&12cG<5w}ohiSM# zr?a_XOTgYuC9P1K1zN{p8?Yt_8j8fAa)6s6wlDGU(WozByntADrpL z=(sYUh7&9MIMFzdnctx2?B5(3p7pW+1chOGpv^Um-O_Mty@fb)?I9!6Jjjb>ry~`* zF7f~-_Z4)(<l}58Wr=Bu?!?g7c{JpK zG&ov{pC0RVdF&9Mba8tZi70H-I)#t*ong+SVWnRWJdcLhw2fuvz0oHmnB=ONOq=l` z**VKV%#+GuSt!m;Pu0;p$U{C>{j=Y4Xcarl%yJAY$G|5Z1Czh#eC6aD=g*&i^`CxU z`sew5SK%h#qI~V2OaH*W$)Ajmk8l0r?+Sn7KmU00UHt#bZ~mF|fA96l_wd+?O@Au= z+isiu8SbCB{okZN{tuJyS1Oh2wpWz@Uy~o=`Il$bOaIWnPJZ_d|g@@$dU2?l(?I+?C&#xQ~2Q;vW4M zi97d8g}c5J_j~;Qv2fDP9{3xn?;Wp8+_`aqJM!V*O5DPKNZb?C68G_MN!&+&u5te# Dnuzu~ diff --git a/test/tests/convergence/residual_convergence/gold/residual_convergence_out.csv b/test/tests/convergence/residual_convergence/gold/residual_convergence_out.csv new file mode 100644 index 000000000000..5f09f2d350a7 --- /dev/null +++ b/test/tests/convergence/residual_convergence/gold/residual_convergence_out.csv @@ -0,0 +1,2 @@ +time,num_nl_its +2,1 diff --git a/test/tests/convergence/residual_convergence/diffusion_convergence.i b/test/tests/convergence/residual_convergence/residual_convergence.i similarity index 67% rename from test/tests/convergence/residual_convergence/diffusion_convergence.i rename to test/tests/convergence/residual_convergence/residual_convergence.i index de545a3d588a..4ee07e0e7b4a 100644 --- a/test/tests/convergence/residual_convergence/diffusion_convergence.i +++ b/test/tests/convergence/residual_convergence/residual_convergence.i @@ -35,7 +35,16 @@ [Convergence] [res_conv] type = ResidualConvergence - nl_abs_tol=1e-6 + # With default tolerances, the solve takes 2 iterations, but with this + # value, the solve takes only 1: + nl_abs_tol = 1e-5 + [] +[] + +[Postprocessors] + [num_nl_its] + type = NumNonlinearIterations + execute_on = 'TIMESTEP_END' [] [] @@ -48,5 +57,6 @@ [] [Outputs] - exodus = true + csv = true + execute_on = 'FINAL' [] diff --git a/test/tests/convergence/residual_convergence/tests b/test/tests/convergence/residual_convergence/tests index 35c81255b8be..40c85b8fcea7 100644 --- a/test/tests/convergence/residual_convergence/tests +++ b/test/tests/convergence/residual_convergence/tests @@ -1,24 +1,11 @@ [Tests] design = 'ResidualConvergence.md' issues = '#25931' + [nonlinear] - type = 'Exodiff' - input = 'diffusion_convergence.i' - exodiff = 'diffusion_convergence_out.e' - requirement = 'The system shall query the convergence status of the nonlinear solvers.' - [] - [default] - type = 'Exodiff' - input = 'diffusion_convergence_default.i' - cli_args = "Outputs/file_base=diffusion_convergence_out" - exodiff = 'diffusion_convergence_out.e' - requirement = 'The system shall query the convergence status using the default convergence object.' - [] - [block_test] - type = 'Exodiff' - input = 'diffusion_convergence_default.i' - cli_args = "Outputs/file_base=diffusion_convergence_out Convergence/res_conv/type=ResidualConvergence Convergence/res_conv/nl_abs_tol=1e-6 Executioner/nonlinear_convergence=res_conv" - exodiff = 'diffusion_convergence_out.e' - requirement = 'The system shall query the convergence status of the nonlinear solvers.' + type = CSVDiff + input = 'residual_convergence.i' + csvdiff = 'residual_convergence_out.csv' + requirement = 'The system shall be able to determine convergence of the nonlinear solve using residual-based criteria through an object.' [] []