www.mooseframework.org
MassLumpedTimeDerivative.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
16 
17 // MOOSE includes
18 #include "Assembly.h"
19 #include "MooseVariable.h"
20 
21 #include "libmesh/quadrature.h"
22 
23 template <>
26 {
28  params.addClassDescription(
29  "Lumped formulation of the time derivative $\\frac{\\partial u}{\\partial t}$. Its "
30  "corresponding weak form is $\\dot{u_i}(\\psi_i, 1)$ where $\\dot{u_i}$ denotes the time "
31  "derivative of the solution coefficient associated with node $i$.");
32  return params;
33 }
34 
36  : TimeKernel(parameters), _u_dot_nodal(_var.nodalValueDot())
37 {
38 }
39 
40 Real
42 {
43  return _test[_i][_qp] * _u_dot_nodal[_i];
44 }
45 
46 Real
48 {
49  return _test[_i][_qp] * _du_dot_du[_qp];
50 }
51 
52 void
54 {
55  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
56 
57  for (_i = 0; _i < _test.size(); _i++)
58  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
59  ke(_i, _i) += _JxW[_qp] * _coord[_qp] * computeQpJacobian();
60 }
InputParameters validParams< MassLumpedTimeDerivative >()
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:140
const VariableTestValue & _test
the current test function
Definition: KernelBase.h:152
QBase *& _qrule
active quadrature rule
Definition: KernelBase.h:137
virtual void computeJacobian() override
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
InputParameters validParams< TimeKernel >()
Definition: TimeKernel.C:26
virtual Real computeQpJacobian() override
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
const VariableValue & _du_dot_du
Derivative of u_dot with respect to u.
Definition: Kernel.h:60
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:117
unsigned int _i
current index for the test function
Definition: KernelBase.h:146
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:143
unsigned int number() const
Get variable number coming from libMesh.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
All time kernels should inherit from this class.
Definition: TimeKernel.h:30
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
MassLumpedTimeDerivative(const InputParameters &parameters)
const VariableValue & _u_dot_nodal
virtual Real computeQpResidual() override
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
MooseVariable & _var
Reference to this Kernel&#39;s MooseVariable object.
Definition: KernelBase.h:120
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:131