19 params.
addClassDescription(
"Conservative form of $\\nabla \\cdot \\vec{v} u$ which in its weak " 20 "form is given by: $(-\\nabla \\psi_i, \\vec{v} u)$.");
22 MooseEnum upwinding_type(
"none full",
"none");
25 "Type of upwinding used. None: Typically results in overshoots and " 26 "undershoots, but numerical diffusion is minimized. Full: Overshoots " 27 "and undershoots are avoided, but numerical diffusion is large");
33 _velocity(coupledVectorValue(
"velocity")),
35 _u_nodal(_var.dofValues()),
95 const unsigned int num_nodes =
_test.size();
108 for (
_i = 0;
_i < num_nodes; ++
_i)
116 Real total_mass_out = 0.0;
121 for (
unsigned int n = 0; n < num_nodes; ++n)
147 for (
unsigned int n = 0; n < num_nodes; ++n)
154 _local_re(n) *= total_mass_out / total_in;
165 Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
167 var->sys().solution().add_vector(
_local_re, var->dofIndices());
179 for (
unsigned int i = 0; i < rows; i++)
182 Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
184 var->sys().solution().add_vector(diag, var->dofIndices());
Advection of the variable by the velocity provided by the user.
virtual void computeResidual() override
Compute this Kernel's contribution to the residual.
virtual void computeResidual() override
Compute this Kernel's contribution to the residual.
static InputParameters validParams()
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
std::vector< MooseVariableFEBase * > _save_in
const MooseArray< Real > & _JxW
The current quadrature point weight value.
unsigned int number() const
Get variable number coming from libMesh.
enum ConservativeAdvection::UpwindingType _upwinding
virtual Real computeQpResidual() override
Compute this Kernel's contribution to the residual at the current quadrature point.
ConservativeAdvection(const InputParameters ¶meters)
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
virtual Real computeQpJacobian() override
Compute this Kernel's contribution to the Jacobian at the current quadrature point.
registerMooseObject("MooseApp", ConservativeAdvection)
unsigned int size() const
The number of elements that can currently be stored in the array.
const VariableTestValue & _test
the current test function
virtual void computeJacobian() override
Compute this Kernel's contribution to the diagonal Jacobian entries.
std::vector< MooseVariableFEBase * > _diag_save_in
void fullUpwind(JacRes res_or_jac)
Calculates the fully-upwind Residual and Jacobian (depending on res_or_jac)
std::vector< Real > _dtotal_mass_out
In the full-upwind scheme d(total_mass_out)/d(variable_at_node_i)
const QBase *const & _qrule
active quadrature rule
bool _has_save_in
The aux variables to save the residual contributions to.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Real negSpeedQp() const
Returns - _grad_test * velocity.
unsigned int _i
current index for the test function
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const VariableValue & _u_nodal
Nodal value of u, used for full upwinding.
const VectorVariableValue & _velocity
advection velocity
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Assembly & _assembly
Reference to this Kernel's assembly object.
UpwindingType
Type of upwinding.
virtual void computeJacobian() override
Compute this Kernel's contribution to the diagonal Jacobian entries.
unsigned int _j
current index for the shape function
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< bool > _upwind_node
In the full-upwind scheme, whether a node is an upwind node.
const VariableTestGradient & _grad_test
gradient of the test function
DenseVector< Number > _local_re
Holds local residual entries as they are accumulated by this Kernel.
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual according to active tags.
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
const VariablePhiValue & _phi
the current shape functions
static InputParameters validParams()
const VariableValue & _u
Holds the solution at current quadrature points.
unsigned int _qp
The current quadrature point index.
JacRes
enum to make the code clearer