https://mooseframework.inl.gov
MassMatrixDGKernel.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "MassMatrixDGKernel.h"
11 #include "MassMatrix.h"
12 
13 registerMooseObject("NavierStokesApp", MassMatrixDGKernel);
14 
17 {
19  params.addClassDescription(
20  "Computes a finite element mass matrix on internal faces meant for use in "
21  "preconditioning schemes which require one");
23  params.addParam<Real>("density", 1, "The density");
24  return params;
25 }
26 
28  : DGKernel(parameters), _density(getParam<Real>("density"))
29 {
30  if (!isParamValid("matrix_tags") && !isParamValid("extra_matrix_tags"))
31  mooseError("One of 'matrix_tags' or 'extra_matrix_tags' must be provided");
32 }
33 
34 Real
36 {
37  mooseAssert(false, "should never be called");
38  return 0;
39 }
40 
41 Real
43 {
44  Real jac = 0;
45 
46  switch (type)
47  {
49  jac = _test[_i][_qp] * _density * _phi[_j][_qp];
50  break;
51 
52  default:
53  jac = 0;
54  break;
55  }
56 
57  return jac;
58 }
This class can be used to compute a mass matrix for facet unknowns, e.g.
static void setMassMatrixParams(InputParameters &params)
MassMatrixDGKernel(const InputParameters &parameters)
static InputParameters validParams()
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
unsigned int _i
DGResidualType
ElementElement
static InputParameters validParams()
unsigned int _qp
unsigned int _j
const std::string & type() const
const VariableTestValue & _test
DGJacobianType
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
registerMooseObject("NavierStokesApp", MassMatrixDGKernel)
void mooseError(Args &&... args) const
void addClassDescription(const std::string &doc_string)
const VariablePhiValue & _phi
bool isParamValid(const std::string &name) const
virtual Real computeQpJacobian(const Moose::DGJacobianType type) override
virtual Real computeQpResidual(Moose::DGResidualType) override