Generic system scalar Depends on PR#328

doc/content/source/scalar_kernels/FuelCycleSystemScalarKernel.md

@@ -0,0 +1,29 @@
+# FuelCycleSystemScalarKernel
+
+!syntax description /ScalarKernels/FuelCycleSystemScalarKernel
+
+## Overview
+
+This object implements the time-dependent mass-transport equations for a
+single variable given an arbitrary number of coupled inputs/outputs (coupled [ScalarVariables](/syntax/Variables))
+and/or sources and sinks. Also tracks tritium decay. All parameters are defined as functors,
+which should allow versatility in accepting a variety of input arguments. An example of using this
+kernel for a system is available, following the [Abdou fuel cycle](/examples/fuel_cycle_Abdou) model.
+
+Rather than using [`ParsedODEKernel`](/syntax/ScalarKernels/ParsedODEKernel) and [`ODETimeDerivative`](/syntax/Scalarkernels/ODETimeDerivative) kernels, the scalar kernels block can be simplified.
+
+!listing test/tests/fuel_cycle_Abdou/fuel_cycle_abdou_generic.i link=false block=ScalarKernels
+
+There is still a need to define appropriate [`Postprocessors`](/syntax/Postprocessors) to inform the scalar kernels, as these cannot be assumed for the general case.
+
+!listing test/tests/fuel_cycle_Abdou/fuel_cycle_abdou_generic.i link=false block=Postprocessors
+
+As a reminder, the system of variables should be defined with the [!param](/Variables/family) attribute set to `SCALAR` for each variable.
+
+!listing test/tests/fuel_cycle_Abdou/fuel_cycle_abdou_generic.i link=false block=Variables/T_01_BZ
+
+!syntax parameters /ScalarKernels/FuelCycleSystemScalarKernel
+
+!syntax inputs /ScalarKernels/FuelCycleSystemScalarKernel
+
+!syntax children /ScalarKernels/FuelCycleSystemScalarKernel

include/scalar_kernels/FuelCycleSystemScalarKernel.h

@@ -0,0 +1,49 @@
+/************************************************************/
+/*                DO NOT MODIFY THIS HEADER                 */
+/*   TMAP8: Tritium Migration Analysis Program, Version 8   */
+/*                                                          */
+/*   Copyright 2021 - 2025 Battelle Energy Alliance, LLC    */
+/*                   ALL RIGHTS RESERVED                    */
+/************************************************************/
+
+#pragma once
+
+#include "ODETimeDerivative.h"
+#include "ADScalarTimeDerivative.h"
+#include "FunctorInterface.h"
+#include "MooseTypes.h"
+
+template <bool is_ad>
+class FuelCycleSystemScalarKernelTempl
+  : public std::conditional<is_ad, ADScalarTimeDerivative, ODETimeDerivative>::type,
+    public FunctorInterface
+{
+  using Base = typename std::conditional<is_ad, ADScalarTimeDerivative, ODETimeDerivative>::type;
+
+public:
+  FuelCycleSystemScalarKernelTempl(const InputParameters & parameters);
+  virtual bool isADObject() const override { return is_ad; };
+  static InputParameters validParams();
+
+protected:
+  virtual GenericReal<is_ad> computeQpResidual() override;
+  virtual Real computeQpJacobian();
+  size_t _n_inputs;
+  size_t _n_other_sources;
+  size_t _n_other_sinks;
+  std::vector<const VariableValue *> _input_vals;
+  std::vector<const Moose::Functor<GenericReal<is_ad>> *> _input_fractions;
+  const Moose::Functor<GenericReal<is_ad>> & _decay_constant;
+  const Moose::Functor<GenericReal<is_ad>> & _residence_time;
+  const Moose::Functor<GenericReal<is_ad>> & _leakage_rate;
+  bool _pseudo_steady_state;
+  bool _disable_residence_time;
+  bool _is_implicit;
+  const Moose::Functor<GenericReal<is_ad>> & _TBR;
+  std::vector<const Moose::Functor<GenericReal<is_ad>> *> _other_sources;
+  std::vector<const Moose::Functor<GenericReal<is_ad>> *> _other_sinks;
+  const Moose::Functor<GenericReal<is_ad>> & _burn_rate;
+};
+
+typedef FuelCycleSystemScalarKernelTempl<false> FuelCycleSystemScalarKernel;
+typedef FuelCycleSystemScalarKernelTempl<true> ADFuelCycleSystemScalarKernel;

src/scalar_kernels/FuelCycleSystemScalarKernel.C

@@ -0,0 +1,184 @@
+/************************************************************/
+/*                DO NOT MODIFY THIS HEADER                 */
+/*   TMAP8: Tritium Migration Analysis Program, Version 8   */
+/*                                                          */
+/*   Copyright 2021 - 2025 Battelle Energy Alliance, LLC    */
+/*                   ALL RIGHTS RESERVED                    */
+/************************************************************/
+#include "FuelCycleSystemScalarKernel.h"
+
+// MOOSE includes
+#include "Assembly.h"
+#include "MooseVariableScalar.h"
+#include "FunctorInterface.h"
+#include "TMAP8PhysicalConstants.h"
+#include "ScalarCoupleable.h"
+
+registerMooseObject("TMAP8App", FuelCycleSystemScalarKernel);
+registerMooseObject("TMAP8App", ADFuelCycleSystemScalarKernel);
+
+template <bool is_ad>
+InputParameters
+FuelCycleSystemScalarKernelTempl<is_ad>::validParams()
+{
+  InputParameters params =
+      is_ad ? ADScalarTimeDerivative::validParams() : ODETimeDerivative::validParams();
+  params += FunctorInterface::validParams();
+  params.addClassDescription("Implements a generic system component.");
+  params.addCoupledVar(
+      "inputs", {}, "Variables which feed into this system. Takes a list of scalar variable names");
+  params.addParam<std::vector<MooseFunctorName>>(
+      "input_fractions",
+      std::vector<MooseFunctorName>({}),
+      "Fraction of upstream variable coming into this system. Must be the same length as 'inputs'");
+  params.addParam<MooseFunctorName>("decay_constant",
+                                    PhysicalConstants::tritium_decay_const,
+                                    "The decay constant of tritium (ln(2)/half-life)");
+  params.addParam<MooseFunctorName>(
+      "residence_time", 1, "The residence time for tritium in this system");
+  params.addParam<MooseFunctorName>(
+      "leakage_rate", 0, "The fractional loss rate for tritium in this system");
+  params.addParam<bool>("steady_state",
+                        false,
+                        "Whether to apply a psuedo steady-state approximation (ignore dt term)");
+  params.addParam<bool>(
+      "is_implicit",
+      false,
+      "Whether an explicit (previous value calculation) or implicit (current value) is used");
+  params.addParam<MooseFunctorName>("TBR", 0.0, "Tritium breeding ratio");
+  params.addParam<std::vector<MooseFunctorName>>(
+      "other_sources",
+      std::vector<MooseFunctorName>({}),
+      "Other tritium sources - terms not dependent on any scalar variables");
+  params.addParam<std::vector<MooseFunctorName>>(
+      "other_sinks",
+      std::vector<MooseFunctorName>({}),
+      "Other tritium sinks - terms not dependent on any scalar variables");
+  params.addParam<MooseFunctorName>("burn_rate", 0.0, "Burn rate of tritium within this system");
+  params.addParam<bool>("disable_residence_time",
+                        false,
+                        "Assume an infinite residence time. (no leakage from this system).");
+  return params;
+}
+
+template <bool is_ad>
+FuelCycleSystemScalarKernelTempl<is_ad>::FuelCycleSystemScalarKernelTempl(
+    const InputParameters & parameters)
+  : Base(parameters),
+    FunctorInterface(this),
+    _n_inputs(ScalarCoupleable::isCoupledScalar("inputs")
+                  ? ScalarCoupleable::coupledScalarComponents("inputs")
+                  : 0),
+    _n_other_sources(
+        this->template getParam<std::vector<MooseFunctorName>>("other_sources").size()),
+    _n_other_sinks(this->template getParam<std::vector<MooseFunctorName>>("other_sinks").size()),
+    _input_vals(_n_inputs),
+    _input_fractions(
+        this->template getParam<std::vector<MooseFunctorName>>("input_fractions").size()),
+    _decay_constant(this->template getFunctor<GenericReal<is_ad>>("decay_constant")),
+    _residence_time(this->template getFunctor<GenericReal<is_ad>>("residence_time")),
+    _leakage_rate(this->template getFunctor<GenericReal<is_ad>>("leakage_rate")),
+    _pseudo_steady_state(this->template getParam<bool>("steady_state")),
+    _disable_residence_time(this->template getParam<bool>("disable_residence_time")),
+    _is_implicit(this->template getParam<bool>("is_implicit")),
+    _TBR(this->template getFunctor<GenericReal<is_ad>>("TBR")),
+    _other_sources(this->template getParam<std::vector<MooseFunctorName>>("other_sources").size()),
+    _other_sinks(this->template getParam<std::vector<MooseFunctorName>>("other_sinks").size()),
+    _burn_rate(this->template getFunctor<GenericReal<is_ad>>("burn_rate"))
+{
+  if (_n_inputs != _input_fractions.size())
+  {
+    mooseError("\"input_fractions\" must be defined with the same length as \"inputs\".");
+  }
+  auto & input_functor_names =
+      MooseBase::getParam<std::vector<MooseFunctorName>>("input_fractions");
+  for (size_t i = 0; i < _n_inputs; ++i)
+  {
+    if ((ScalarCoupleable::coupledScalar("inputs", i)) == (Base::_var).number())
+    {
+      mooseError("Primary variable cannot be listed as a coupled variable.");
+    }
+    _input_vals[i] = &(ScalarCoupleable::coupledScalarValue("inputs", i));
+    _input_fractions[i] = &(this->template getFunctor<GenericReal<is_ad>>(input_functor_names[i]));
+  }
+  auto & other_source_names = MooseBase::getParam<std::vector<MooseFunctorName>>("other_sources");
+  for (size_t i = 0; i < _n_other_sources; ++i)
+  {
+    _other_sources[i] = &(this->template getFunctor<GenericReal<is_ad>>(other_source_names[i]));
+  }
+  auto & other_sink_names = MooseBase::getParam<std::vector<MooseFunctorName>>("other_sinks");
+  for (size_t i = 0; i < _n_other_sinks; ++i)
+  {
+    _other_sinks[i] = &(this->template getFunctor<GenericReal<is_ad>>(other_sink_names[i]));
+  }
+}
+
+template <bool is_ad>
+GenericReal<is_ad>
+FuelCycleSystemScalarKernelTempl<is_ad>::computeQpResidual()
+{
+  GenericReal<is_ad> partial_residual = 0;
+  const Moose::ElemArg _qp = Moose::ElemArg();
+  const int _i = 0;
+  const auto _state = _is_implicit ? Moose::currentState() : Moose::oldState();
+  for (unsigned int i = 0; i < _n_inputs; ++i)
+  {
+    partial_residual += -(*(_input_vals[i]))[_i] * (*(_input_fractions[i]))(_qp, _state);
+  }
+  for (unsigned int i = 0; i < _n_other_sources; ++i)
+  {
+    partial_residual += -(*(_other_sources[i]))(_qp, _state);
+  }
+  for (unsigned int i = 0; i < _n_other_sinks; ++i)
+  {
+    partial_residual += (*(_other_sinks[i]))(_qp, _state);
+  }
+  partial_residual += -_TBR(_qp, _state) * _burn_rate(_qp, _state);
+  if (!_disable_residence_time)
+  {
+    partial_residual += _leakage_rate(_qp, _state) / _residence_time(_qp, _state) * Base::_u[_i];
+    partial_residual += Base::_u[_i] / _residence_time(_qp, _state);
+  }
+  partial_residual += _decay_constant(_qp, _state) * Base::_u[_i];
+  if (!_pseudo_steady_state)
+    partial_residual += Base::_u_dot[_i];
+  return partial_residual;
+}
+
+template <bool is_ad>
+Real
+FuelCycleSystemScalarKernelTempl<is_ad>::computeQpJacobian()
+{
+  if constexpr (!is_ad)
+  {
+    Real partial_residual = 0;
+    const Moose::ElemArg _qp = Moose::ElemArg();
+    const int _i = 0;
+    const auto _state = _is_implicit ? Moose::currentState() : Moose::oldState();
+    if (!_disable_residence_time)
+    {
+      partial_residual += _leakage_rate(_qp, _state) / _residence_time(_qp, _state);
+      partial_residual += 1 / _residence_time(_qp, _state);
+    }
+    partial_residual += _decay_constant(_qp, _state);
+    if (!_pseudo_steady_state)
+      partial_residual += Base::_du_dot_du[_i];
+    return partial_residual;
+  }
+  else
+  {
+    mooseError("computeQpJacobian() should not be called in AD mode");
+    return 0;
+  }
+}
+
+template <>
+Real
+FuelCycleSystemScalarKernelTempl<true>::computeQpJacobian()
+{
+  mooseError("Internal error, calling computeQpJacobian in AD class.");
+  return 0.0;
+}
+
+template class FuelCycleSystemScalarKernelTempl<false>;
+template class FuelCycleSystemScalarKernelTempl<true>;

test/tests/fuel_cycle_Abdou/fuel_cycle.i

@@ -79,7 +79,6 @@
 # on the ParsedODEKernels due to moving the terms to the
 # left hand side.
 
-
 [ScalarKernels]
   [I1t]
     type = ODETimeDerivative
@@ -296,7 +295,7 @@
     value = 0.95
   []
   [TBR] #According to the PhD Thesis referenced in the paper,
-        # this is the required Tritium Breeding Ratio (TBR)
+    # this is the required Tritium Breeding Ratio (TBR)
     type = ConstantPostprocessor
     execute_on = 'TIMESTEP_BEGIN INITIAL LINEAR NONLINEAR'
     value = 1.9247
@@ -339,7 +338,7 @@
   [epsilon8] #FCU
     type = ConstantPostprocessor
     execute_on = 'TIMESTEP_BEGIN INITIAL LINEAR NONLINEAR'
-    value =  1e-4
+    value = 1e-4
   []
   [epsilon9] #ISS
     type = ConstantPostprocessor
@@ -419,7 +418,7 @@
   [tdecay]
     type = ConstantPostprocessor
     execute_on = 'TIMESTEP_BEGIN INITIAL LINEAR NONLINEAR'
-    value = ${fparse log(2)/388800000.0}
+    value = '${fparse log(2)/388800000.0}'
   []
   [tritium_burn_fraction]
     type = ConstantPostprocessor