https://mooseframework.inl.gov
ArrayHFEMDiffusion.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 "ArrayHFEMDiffusion.h"
11 
13 
16 {
18  params.addClassDescription("Imposes the constraints on internal sides with HFEM.");
19  return params;
20 }
21 
23  : ArrayDGLowerDKernel(parameters)
24 {
25 }
26 
27 void
29 {
30  switch (type)
31  {
32  case Moose::Element:
33  r -= _lambda[_qp] * _test[_i][_qp];
34  break;
35 
36  case Moose::Neighbor:
37  r += _lambda[_qp] * _test_neighbor[_i][_qp];
38  break;
39  }
40 }
41 
42 void
44 {
45  r += (_u_neighbor[_qp] - _u[_qp]) * _test_lambda[_i][_qp];
46 }
47 
50 {
51  RealEigenVector r = RealEigenVector::Zero(_count);
52  switch (type)
53  {
55  return RealEigenVector::Constant(_count, -_test_lambda[_i][_qp] * _phi[_j][_qp]);
56 
58  return RealEigenVector::Constant(_count, _test_lambda[_i][_qp] * _phi_neighbor[_j][_qp]);
59 
61  return RealEigenVector::Constant(_count, -_phi_lambda[_j][_qp] * _test[_i][_qp]);
62 
64  return RealEigenVector::Constant(_count, _phi_lambda[_j][_qp] * _test_neighbor[_i][_qp]);
65 
66  default:
67  break;
68  }
69 
70  return r;
71 }
const ArrayVariableValue & _u
Holds the current solution at the current quadrature point on the face.
Definition: ArrayDGKernel.h:98
The DGKernel class is responsible for calculating the residuals for various physics on internal sides...
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
Definition: DGKernelBase.h:136
static InputParameters validParams()
Factory constructor initializes all internal references needed for residual computation.
DGResidualType
Definition: MooseTypes.h:743
const ArrayVariableTestValue & _test_lambda
test functions
const ArrayVariableTestValue & _test
test functions
unsigned int _qp
Definition: DGKernelBase.h:134
const unsigned int _count
Number of components of the array variable.
unsigned int _j
Definition: DGKernelBase.h:136
virtual void computeQpResidual(Moose::DGResidualType type, RealEigenVector &residual) override
This is the virtual that derived classes should override for computing the residual on neighboring el...
virtual void computeLowerDQpResidual(RealEigenVector &residual) override
Method for computing the Lower part of residual at quadrature points, to be filled in residual...
const std::string & type() const
Get the type of this class.
Definition: MooseBase.h:51
registerMooseObject("MooseApp", ArrayHFEMDiffusion)
virtual RealEigenVector computeLowerDQpJacobian(Moose::ConstraintJacobianType type) override
Computes one of the five pieces of Jacobian involving lower-d at quadrature points.
ConstraintJacobianType
Definition: MooseTypes.h:796
ArrayHFEMDiffusion(const InputParameters &parameters)
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...
const ArrayVariablePhiValue & _phi_neighbor
Side shape function.
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
const ArrayVariablePhiValue & _phi_lambda
Shape functions.
const ArrayVariableValue & _lambda
Holds the current solution at the current quadrature point on the face.
const ArrayVariableTestValue & _test_neighbor
Side test function.
const ArrayVariableValue & _u_neighbor
Holds the current solution at the current quadrature point.
static InputParameters validParams()
const ArrayVariablePhiValue & _phi
Shape functions.