23 params.
addClassDescription(
"Conservative form of $\\nabla \\cdot \\vec{v} u$ which in its weak " 24 "form is given by: $(-\\nabla \\psi_i, \\vec{v} u)$. Velocity can be " 25 "given as 1) a variable, for which the gradient is automatically " 26 "taken, 2) a vector variable, or a 3) vector material.");
27 params.
addParam<MaterialPropertyName>(
28 "velocity_scalar_coef",
30 "Name of material property multiplied against the velocity to scale advection strength.");
31 MooseEnum upwinding_type(
"none full",
"none");
34 "Type of upwinding used. None: Typically results in overshoots and " 35 "undershoots, but numerical diffusion is minimized. Full: Overshoots " 36 "and undershoots are avoided, but numerical diffusion is large");
37 params.
addParam<MaterialPropertyName>(
"advected_quantity",
38 "An optional material property to be advected. If not " 39 "supplied, then the variable will be used.");
41 "velocity_as_variable_gradient",
42 "Gradient of this coupled variable is used to define the advection velocity. " 43 "Can be supplied instead of velocity material or velocity variable.");
54 params.
addParam<MaterialPropertyName>(
"velocity_material",
"Velocity vector given as a material");
63 params.
addCoupledVar(
"velocity_variable",
"Velocity vector given as a variable");
64 params.
addParam<MaterialPropertyName>(
"velocity",
"Velocity vector given as a material");
65 params.
deprecateParam(
"velocity",
"velocity_material",
"12/31/2025");
72 _scalar(this->template getGenericMaterialProperty<
Real, is_ad>(
"velocity_scalar_coef")),
73 _coupled_variable_present(isParamValid(
"velocity_as_variable_gradient")),
74 _coupled_variable_var(_coupled_variable_present ? coupled(
"velocity_as_variable_gradient") : 0),
76 _coupled_variable_present
77 ? &this->template coupledGenericGradient<is_ad>(
"velocity_as_variable_gradient")
78 : (this->isParamValid(
"velocity_variable")
79 ? &this->template coupledGenericVectorValue<is_ad>(
"velocity_variable")
80 : (this->isParamValid(
"velocity_material")
85 _user_supplied_adv_quant(isParamValid(
"advected_quantity")),
87 _user_supplied_adv_quant
88 ? this->template getGenericMaterialProperty<
Real, is_ad>(
"advected_quantity").
get()
92 _u_nodal(_var.template genericDofValues<is_ad>()),
98 "Use a different kernel (i.e., diffusion) if the gradient used as the velocity is " 99 "the same as the member variable");
103 "Upwinding is not compatible with an advected quantity that is not the primary variable.");
109 "velocity_as_variable_gradient",
110 "One and only one of the following input variables must be specified: velocity_variable, " 111 "velocity_material, or velocity_as_variable_gradient.");
115 "_local_ke not computed for global AD indexing. Save-in is deprecated anyway. Use " 116 "the tagging system instead.");
119 template <
bool is_ad>
123 return -_grad_test[_i][_qp] * (*_velocity)[_qp] * _scalar[_qp];
126 template <
bool is_ad>
132 return negSpeedQp() * _adv_quant[_qp];
141 if (!_user_supplied_adv_quant)
142 return negSpeedQp() * _phi[_j][_qp];
150 mooseError(
"Internal error, should never get here when using AD");
160 if (_coupled_variable_present && _coupled_variable_var == jvar)
161 return -_grad_test[_i][_qp] * _grad_phi[_j][_qp] * _adv_quant[_qp] * _scalar[_qp];
170 mooseError(
"Internal error, should never get here when using AD");
174 template <
bool is_ad>
180 case UpwindingType::none:
183 case UpwindingType::full:
184 fullUpwind(JacRes::CALCULATE_RESIDUAL);
189 template <
bool is_ad>
195 case UpwindingType::none:
198 case UpwindingType::full:
199 fullUpwind(JacRes::CALCULATE_JACOBIAN);
204 template <
bool is_ad>
209 const unsigned int num_nodes = _test.size();
213 _my_local_re.resize(_var.dofIndices().size());
215 if (!is_ad && (res_or_jac == JacRes::CALCULATE_JACOBIAN))
216 prepareMatrixTag(this->_assembly, _var.number(), _var.number());
221 _upwind_node.resize(num_nodes);
222 for (_i = 0; _i < num_nodes; ++_i)
224 for (_qp = 0; _qp < this->_qrule->n_points(); _qp++)
225 _my_local_re(_i) += this->_JxW[_qp] * this->_coord[_qp] * negSpeedQp();
232 if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
233 _dtotal_mass_out.assign(num_nodes, 0.0);
239 if constexpr (!is_ad)
240 if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
242 if (_test.size() == _phi.size())
249 _local_ke(n, n) += _my_local_re(n);
251 _dtotal_mass_out[n] += _local_ke(n, n);
253 _my_local_re(n) *= _u_nodal[n];
254 total_mass_out += _my_local_re(n);
257 total_in -= _my_local_re(n);
263 if (!_upwind_node[n])
265 if constexpr (!is_ad)
266 if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
267 for (_j = 0; _j < _phi.size(); _j++)
268 _local_ke(n, _j) += _my_local_re(n) * _dtotal_mass_out[_j] / total_in;
269 _my_local_re(n) *= total_mass_out / total_in;
273 if (res_or_jac == JacRes::CALCULATE_RESIDUAL)
275 this->addResiduals(this->_assembly, _my_local_re, _var.dofIndices(), _var.scalingFactor());
277 if (this->_has_save_in)
279 Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
280 for (
const auto & var : this->_save_in)
285 if (res_or_jac == JacRes::CALCULATE_JACOBIAN)
288 this->addJacobian(this->_assembly, _my_local_re, _var.dofIndices(), _var.scalingFactor());
290 accumulateTaggedLocalMatrix();
Moose::GenericType< Real, is_ad > GenericReal
virtual void computeResidual() override
Compute this Kernel's contribution to the residual.
virtual Real computeQpOffDiagJacobian(unsigned int jvar) override
For coupling standard variables.
virtual void computeResidual() override
Compute this Kernel's contribution to the residual.
const unsigned int _coupled_variable_var
Coupled variable variable number.
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
void paramError(const std::string ¶m, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
unsigned int number() const
Get variable number coming from libMesh.
static InputParameters validParams()
T * get(const std::unique_ptr< T > &u)
The MooseUtils::get() specializations are used to support making forwards-compatible code changes fro...
const MooseArray< GenericRealVectorValue< is_ad > > * _velocity
advection velocity
UpwindingType
Type of upwinding.
const bool _coupled_variable_present
Flag to determine if coupled variable is present.
enum ConservativeAdvectionTempl::UpwindingType _upwinding
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
virtual GenericReal< is_ad > negSpeedQp() const
Returns - _grad_test * velocity.
registerMooseObject("MooseApp", ConservativeAdvection)
static InputParameters generalParams()
virtual GenericReal< is_ad > computeQpResidual() override
Compute this Kernel's contribution to the residual at the current quadrature point.
ConservativeAdvectionTempl(const InputParameters ¶meters)
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Advection of the variable by the velocity provided by the user.
void fullUpwind(JacRes res_or_jac)
Calculates the fully-upwind Residual and Jacobian (depending on res_or_jac)
virtual Real computeQpJacobian() override
Compute this Kernel's contribution to the Jacobian at the current quadrature point.
virtual void computeJacobian() override
Compute this Kernel's contribution to the diagonal Jacobian entries.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
JacRes
enum to make the code clearer
IntRange< T > make_range(T beg, T end)
static InputParameters validParams()
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
virtual void computeJacobian() override
Compute this Kernel's contribution to the diagonal Jacobian entries.