Merge pull request idaholab#28806 from joshuahansel/sol-diff-norm

Added AverageVariableChange post-processor

framework/doc/content/source/postprocessors/AverageVariableChange.md

@@ -0,0 +1,41 @@
+# AverageVariableChange
+
+This post-processor computes a volume-weighted norm of the change of a variable
+or auxiliary variable over a time step or a nonlinear iteration:
+
+!equation
+||\delta u||_{\Omega,p} = \left(\frac{1}{|\Omega|} \int\limits_\Omega |\delta u|^p \, d\Omega\right)^{\frac{1}{p}} \,,
+
+where
+
+- $\delta u$ is the variable change of interest,
+- $\Omega$ is the domain over which the norm is taken,
+- $|\Omega|$ is the volume of the domain, and
+- $p$ is the order of the norm (either 1 or 2).
+
+Denoting the current time index by $n$ and the current nonlinear iteration index
+by $(\ell)$, the variable change may be with respect to either the old time value,
+
+!equation
+\delta u \equiv u_n^{(\ell)} - u_{n-1}
+
+or the previous nonlinear iteration value,
+
+!equation
+\delta u \equiv u_n^{(\ell)} - u_n^{(\ell-1)}
+
+as determined by the value of [!param](/Postprocessors/AverageVariableChange/change_over).
+
+Volume-weighting the norm produces a value that is independent of the size of
+the domain. Without doing this, there is a domain volume bias in the norm.
+Suppose $\delta u$ is constant over the domain with a value $\overline{\delta u}$:
+
+!equation
+||\delta u||_{p} = \left(\int\limits_\Omega |\overline{\delta u}|^p \, d\Omega\right)^{\frac{1}{p}}
+  = |\overline{\delta u}| |\Omega|^{\frac{1}{p}} \,.
+
+!syntax parameters /Postprocessors/AverageVariableChange
+
+!syntax inputs /Postprocessors/AverageVariableChange
+
+!syntax children /Postprocessors/AverageVariableChange

framework/include/postprocessors/AverageVariableChange.h

@@ -0,0 +1,45 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#pragma once
+
+#include "ElementIntegralVariablePostprocessor.h"
+
+/**
+ * Computes the volume-weighted L1 or L2 norm of the change of a variable over a
+ * time step or between nonlinear iterations.
+ */
+class AverageVariableChange : public ElementIntegralVariablePostprocessor
+{
+public:
+  static InputParameters validParams();
+
+  AverageVariableChange(const InputParameters & parameters);
+
+  virtual void initialize() override;
+  virtual void execute() override;
+  virtual void threadJoin(const UserObject & y) override;
+  virtual void finalize() override;
+  virtual Real getValue() const override;
+  virtual Real computeQpIntegral() override;
+
+protected:
+  /// Interval over which to compute change
+  const MooseEnum & _change_over;
+  /// "Old" solution value in the change
+  const VariableValue & _u_change_old;
+
+  /// Selected norm
+  const MooseEnum & _norm;
+  /// Norm exponent
+  const unsigned int _norm_exponent;
+
+  /// Subdomain(s) volume
+  Real _volume;
+};

framework/include/systems/AuxiliarySystem.h

@@ -131,7 +131,8 @@ class AuxiliarySystem : public SystemBase, public PerfGraphInterface
   virtual System & system() override { return _sys; }
   virtual const System & system() const override { return _sys; }
 
-  virtual void setPreviousNewtonSolution();
+  /// Copies the current solution into the previous nonlinear iteration solution
+  virtual void copyCurrentIntoPreviousNL();
 
   void setScalarVariableCoupleableTags(ExecFlagType type);
 

framework/src/postprocessors/AverageVariableChange.C

@@ -0,0 +1,84 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#include "AverageVariableChange.h"
+#include "MooseUtils.h"
+
+registerMooseObject("MooseApp", AverageVariableChange);
+
+InputParameters
+AverageVariableChange::validParams()
+{
+  InputParameters params = ElementIntegralVariablePostprocessor::validParams();
+
+  params.addClassDescription("Computes the volume-weighted L1 or L2 norm of the change of a "
+                             "variable over a time step or between nonlinear iterations.");
+
+  MooseEnum change_over("time_step nonlinear_iteration");
+  params.addRequiredParam<MooseEnum>(
+      "change_over", change_over, "Interval over which to compute the change");
+
+  MooseEnum norm("L1 L2");
+  params.addRequiredParam<MooseEnum>("norm", norm, "Type of norm to compute");
+
+  return params;
+}
+
+AverageVariableChange::AverageVariableChange(const InputParameters & parameters)
+  : ElementIntegralVariablePostprocessor(parameters),
+    _change_over(getParam<MooseEnum>("change_over")),
+    _u_change_old(_change_over == "time_step" ? coupledValueOld("variable")
+                                              : coupledValuePreviousNL("variable")),
+    _norm(getParam<MooseEnum>("norm")),
+    _norm_exponent(_norm == "L1" ? 1 : 2),
+    _volume(0)
+{
+}
+
+void
+AverageVariableChange::initialize()
+{
+  ElementIntegralVariablePostprocessor::initialize();
+  _volume = 0;
+}
+
+void
+AverageVariableChange::execute()
+{
+  ElementIntegralVariablePostprocessor::execute();
+  _volume += _current_elem_volume;
+}
+
+void
+AverageVariableChange::threadJoin(const UserObject & y)
+{
+  ElementIntegralVariablePostprocessor::threadJoin(y);
+  const auto & pps = static_cast<const AverageVariableChange &>(y);
+  _volume += pps._volume;
+}
+
+void
+AverageVariableChange::finalize()
+{
+  gatherSum(_volume);
+  gatherSum(_integral_value);
+}
+
+Real
+AverageVariableChange::getValue() const
+{
+  mooseAssert(!MooseUtils::absoluteFuzzyEqual(_volume, 0.0), "Volume must be nonzero.");
+  return std::pow(_integral_value / _volume, 1.0 / _norm_exponent);
+}
+
+Real
+AverageVariableChange::computeQpIntegral()
+{
+  return std::pow(std::abs(_u[_qp] - _u_change_old[_qp]), _norm_exponent);
+}

framework/src/problems/FEProblemBase.C

@@ -7405,7 +7405,7 @@ FEProblemBase::computePostCheck(NonlinearImplicitSystem & sys,
   if (vectorTagExists(Moose::PREVIOUS_NL_SOLUTION_TAG))
   {
     _current_nl_sys->setPreviousNewtonSolution(old_soln);
-    _aux->setPreviousNewtonSolution();
+    _aux->copyCurrentIntoPreviousNL();
   }
 
   // MOOSE doesn't change the search_direction

framework/src/systems/AuxiliarySystem.C

@@ -792,10 +792,10 @@ AuxiliarySystem::needMaterialOnSide(BoundaryID bnd_id)
 }
 
 void
-AuxiliarySystem::setPreviousNewtonSolution()
+AuxiliarySystem::copyCurrentIntoPreviousNL()
 {
-  // Evaluate aux variables to get the solution vector
-  compute(EXEC_LINEAR);
+  if (solutionPreviousNewton())
+    *solutionPreviousNewton() = *currentSolution();
 }
 
 template <typename AuxKernelType>

test/tests/postprocessors/average_variable_change/base.i

@@ -0,0 +1,72 @@
+aux_execute_on = 'INITIAL NONLINEAR_CONVERGENCE TIMESTEP_END'
+pp_execute_on = 'INITIAL NONLINEAR_CONVERGENCE TIMESTEP_END'
+output_execute_on = 'INITIAL NONLINEAR TIMESTEP_END'
+
+[Mesh]
+  type = GeneratedMesh
+  dim = 2
+  nx = 10
+  ny = 10
+  xmax = 2
+  ymax = 2
+[]
+
+[Variables]
+  [sol]
+  []
+[]
+
+[AuxVariables]
+  [aux]
+  []
+[]
+
+[AuxKernels]
+  [aux_kernel]
+    type = ParsedAux
+    variable = aux
+    functor_names = 'sol'
+    functor_symbols = 'u'
+    expression = 'u'
+    execute_on = ${aux_execute_on}
+  []
+[]
+
+[Kernels]
+  [diff]
+    type = Diffusion
+    variable = sol
+  []
+[]
+
+[BCs]
+  [left]
+    type = DirichletBC
+    variable = sol
+    boundary = left
+    value = 0
+  []
+  [right]
+    type = DirichletBC
+    variable = sol
+    boundary = right
+    value = 1
+  []
+[]
+
+[Problem]
+  previous_nl_solution_required = true
+[]
+
+[Executioner]
+  solve_type = 'PJFNK'
+  petsc_options_iname = '-pc_type'
+  petsc_options_value = 'hypre'
+  nl_rel_tol = 0
+  nl_abs_tol = 1e-8
+[]
+
+[Outputs]
+  csv = true
+  execute_on = ${output_execute_on}
+[]

test/tests/postprocessors/average_variable_change/gold/steady_out.csv

@@ -0,0 +1,6 @@
+time,aux_avg,aux_diff_l1,sol_avg,sol_diff_l1,sol_diff_l2
+0,0,0,0,0,0
+0.001,0.05,0.05,0.05,0.05,0.18257418583506
+0.002,0.5000000662596,0.4500000662596,0.5000000662596,0.4500000662596,0.51961533896988
+0.003,0.50000000000575,1.1618620291334e-07,0.50000000000575,1.1618620291334e-07,1.4653850459149e-07
+1,0.50000000000575,1.1618620291334e-07,0.50000000000575,1.1618620291334e-07,1.4653850459149e-07

test/tests/postprocessors/average_variable_change/gold/transient_out.csv

@@ -0,0 +1,10 @@
+time,aux_nlit,aux_ts,sol_nlit,sol_ts
+0,0,0,0,0
+0.0001,0.05,0.05,0.05,0.05
+0.0002,0.11021258306293,0.16021258306293,0.11021258306293,0.16021258306293
+0.0003,1.5959019350813e-07,0.16021256040952,1.5959019350813e-07,0.16021256040952
+0.1,0,0.16021256040952,0,0.16021256040952
+0.1001,0,0,0,0
+0.1002,0.075808624717138,0.075808624717138,0.075808624717138,0.075808624717138
+0.1003,5.3715533178459e-08,0.0758085898486,5.3715533178459e-08,0.0758085898486
+0.2,0,0.0758085898486,0,0.0758085898486

test/tests/postprocessors/average_variable_change/steady.i

@@ -0,0 +1,40 @@
+!include base.i
+
+[Postprocessors]
+  [sol_avg]
+    type = ElementAverageValue
+    variable = sol
+    execute_on = ${pp_execute_on}
+  []
+  [sol_diff_l1]
+    type = AverageVariableChange
+    variable = sol
+    change_over = nonlinear_iteration
+    norm = l1
+    execute_on = ${pp_execute_on}
+  []
+  [sol_diff_l2]
+    type = AverageVariableChange
+    variable = sol
+    change_over = nonlinear_iteration
+    norm = l2
+    execute_on = ${pp_execute_on}
+  []
+
+  [aux_avg]
+    type = ElementAverageValue
+    variable = aux
+    execute_on = ${pp_execute_on}
+  []
+  [aux_diff_l1]
+    type = AverageVariableChange
+    variable = aux
+    change_over = nonlinear_iteration
+    norm = l1
+    execute_on = ${pp_execute_on}
+  []
+[]
+
+[Executioner]
+  type = Steady
+[]

test/tests/postprocessors/average_variable_change/tests

@@ -0,0 +1,22 @@
+[Tests]
+  design = 'AverageVariableChange.md'
+  issues = '#28801'
+
+  [test]
+    requirement = 'The system shall be able to compute volume-weighted norms of the difference of a variable between nonlinear iterations and time steps'
+    [steady]
+      type = CSVDiff
+      input = 'steady.i'
+      csvdiff = 'steady_out.csv'
+      abs_zero = 1e-5
+      detail = 'for a steady executioner.'
+    []
+    [transient]
+      type = CSVDiff
+      input = 'transient.i'
+      csvdiff = 'transient_out.csv'
+      abs_zero = 1e-5
+      detail = 'for a transient executioner.'
+    []
+  []
+[]