diff --git a/framework/doc/content/source/interfacekernels/ADMatInterfaceReaction.md b/framework/doc/content/source/interfacekernels/ADMatInterfaceReaction.md new file mode 100644 index 000000000000..7f5d8860e418 --- /dev/null +++ b/framework/doc/content/source/interfacekernels/ADMatInterfaceReaction.md @@ -0,0 +1,41 @@ +# ADMatInterfaceReaction + +## Description + +Specie M transports between two domains (domain 1 and domain 2), at the interface consider the following reaction is taking place: + +\begin{equation} +M(1)\xrightleftharpoons[k_b]{k_f}M(2) +\end{equation} + +With the first order reaction rate assuming a quasi-steady-state + +\begin{equation} +\textrm{Reaction Rate} = \frac {\partial C_1} {\partial t} = k_f C_1 - k_b C_2 \approx 0 +\end{equation} + +where $C_1$ is the specie concentration in domain 1, $C_2$ is the specie concentration in domain 2, $k_f$ is the forward reaction coefficient, and $k_b$ is the backward reaction coefficient. `ADMatInterfaceReaction` object is used to impose this condition. + +[InterfaceDiffusion.md] is also used in this case to control flux at the interface. + +\begin{equation} +D_1 \frac {\partial C_1} {\partial n} = D_2 \frac {\partial C_2} {\partial n} +\end{equation} + +However, the flux is not [well-defined](https://en.wikipedia.org/wiki/Well-defined_expression) across the interface. The `ADMatInterfaceReaction` interfacekernel applies a condition to constrain the potential discontinuity across the interface. + +Both kernels at the interface work together to give full mathematical and physical meaning of the problem. Together, the implicit equations represented by `ADMatInterfaceReaction` and [InterfaceDiffusion.md] combine to provide the following relationship at the interface. + +\begin{equation} +D_1 \frac {\partial C_1} {\partial n} = D_2 \frac {\partial C_2} {\partial n} + k_f C_1 - k_b C_2 +\end{equation} + +## Example Input Syntax + +!listing test/tests/interfacekernels/1d_interface/ADMatreaction_1D_steady.i start=[interface_reaction] end=[] include-end=true + +!listing test/tests/interfacekernels/1d_interface/ADMatreaction_1D_transient.i start=[interface_reaction] end=[] include-end=true + +!syntax parameters /InterfaceKernels/InterfaceReaction + +!syntax children /InterfaceKernels/ADMatInterfaceReaction diff --git a/framework/include/interfacekernels/ADMatInterfaceReaction.h b/framework/include/interfacekernels/ADMatInterfaceReaction.h new file mode 100644 index 000000000000..edf376188a2b --- /dev/null +++ b/framework/include/interfacekernels/ADMatInterfaceReaction.h @@ -0,0 +1,33 @@ +//* 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 "ADInterfaceKernel.h" + +/** + * Implements a reaction to establish ReactionRate=k_f*u-k_b*v + * at interface. + */ +class ADMatInterfaceReaction : public ADInterfaceKernel +{ +public: + static InputParameters validParams(); + + ADMatInterfaceReaction(const InputParameters & parameters); + +protected: + virtual ADReal computeQpResidual(Moose::DGResidualType type) override; + + /// Forward reaction rate coefficient + const ADMaterialProperty & _kf; + + /// Backward reaction rate coefficient + const ADMaterialProperty & _kb; +}; diff --git a/framework/src/interfacekernels/ADMatInterfaceReaction.C b/framework/src/interfacekernels/ADMatInterfaceReaction.C new file mode 100644 index 000000000000..8f45cb26002c --- /dev/null +++ b/framework/src/interfacekernels/ADMatInterfaceReaction.C @@ -0,0 +1,54 @@ +//* 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 "ADMatInterfaceReaction.h" + +registerMooseObject("MooseApp", ADMatInterfaceReaction); + +InputParameters +ADMatInterfaceReaction::validParams() + +{ + InputParameters params = ADInterfaceKernel::validParams(); + params.addParam("forward_rate", "kf", "Forward reaction rate coefficient."); + params.addParam( + "backward_rate", "kb", "Backward reaction rate coefficient."); + params.addClassDescription("Implements a reaction to establish ReactionRate=k_f*u-k_b*v " + "at interface."); + return params; +} + +ADMatInterfaceReaction::ADMatInterfaceReaction(const InputParameters & parameters) + : ADInterfaceKernel(parameters), + _kf(getADMaterialProperty("forward_rate")), + _kb(getNeighborADMaterialProperty("backward_rate")) +{ +} + +ADReal +ADMatInterfaceReaction::computeQpResidual(Moose::DGResidualType type) +{ + ADReal r = 0; + switch (type) + { + // Move all the terms to the LHS to get residual, for primary domain + // Residual = kf*u - kb*v = kf*u - kb*v + // Weak form for primary domain is: (test, kf*u - kb*v) + case Moose::Element: + r = _test[_i][_qp] * (_kf[_qp] * _u[_qp] - _kb[_qp] * _neighbor_value[_qp]); + break; + + // Similarly, weak form for secondary domain is: -(test, kf*u - kb*v), + // flip the sign because the direction is opposite. + case Moose::Neighbor: + r = -_test_neighbor[_i][_qp] * (_kf[_qp] * _u[_qp] - _kb[_qp] * _neighbor_value[_qp]); + break; + } + return r; +} diff --git a/test/tests/interfacekernels/1d_interface/ADMatreaction_1D_steady.i b/test/tests/interfacekernels/1d_interface/ADMatreaction_1D_steady.i new file mode 100644 index 000000000000..cef35499d532 --- /dev/null +++ b/test/tests/interfacekernels/1d_interface/ADMatreaction_1D_steady.i @@ -0,0 +1,170 @@ +# Steady-state test for the ADMatInterfaceReaction kernel. +# +# Specie M transport from domain 1 (0<=x<=1) to domain 2 (1