www.mooseframework.org
CNSFVKernel.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 #include "CNSFVKernel.h"
9 
10 template <>
11 InputParameters
13 {
14  InputParameters params = validParams<DGKernel>();
15 
16  params.addClassDescription("A DGKernel for the CNS equations.");
17 
18  MooseEnum component("mass x-momentum y-momentum z-momentum total-energy");
19 
20  params.addParam<MooseEnum>("component", component, "Choose one of the conservation equations");
21 
22  params.addRequiredCoupledVar("rho", "Conserved variable: rho");
23 
24  params.addRequiredCoupledVar("rhou", "Conserved variable: rhou");
25 
26  params.addCoupledVar("rhov", "Conserved variable: rhov");
27 
28  params.addCoupledVar("rhow", "Conserved variable: rhow");
29 
30  params.addRequiredCoupledVar("rhoe", "Conserved variable: rhoe");
31 
32  params.addRequiredParam<UserObjectName>("flux", "The name of internal side flux object to use");
33 
34  return params;
35 }
36 
37 CNSFVKernel::CNSFVKernel(const InputParameters & parameters)
38  : DGKernel(parameters),
39  _component(getParam<MooseEnum>("component")),
40  _rhoc1(coupledValue("rho")),
41  _rhouc1(coupledValue("rhou")),
42  _rhovc1(isCoupled("rhov") ? coupledValue("rhov") : _zero),
43  _rhowc1(isCoupled("rhow") ? coupledValue("rhow") : _zero),
44  _rhoec1(coupledValue("rhoe")),
45  _rhoc2(coupledNeighborValue("rho")),
46  _rhouc2(coupledNeighborValue("rhou")),
47  _rhovc2(isCoupled("rhov") ? coupledNeighborValue("rhov") : _zero),
48  _rhowc2(isCoupled("rhow") ? coupledNeighborValue("rhow") : _zero),
49  _rhoec2(coupledNeighborValue("rhoe")),
50  _rho1(getMaterialProperty<Real>("rho")),
51  _rhou1(getMaterialProperty<Real>("rhou")),
52  _rhov1(getMaterialProperty<Real>("rhov")),
53  _rhow1(getMaterialProperty<Real>("rhow")),
54  _rhoe1(getMaterialProperty<Real>("rhoe")),
55  _rho2(getNeighborMaterialProperty<Real>("rho")),
56  _rhou2(getNeighborMaterialProperty<Real>("rhou")),
57  _rhov2(getNeighborMaterialProperty<Real>("rhov")),
58  _rhow2(getNeighborMaterialProperty<Real>("rhow")),
59  _rhoe2(getNeighborMaterialProperty<Real>("rhoe")),
60  _flux(getUserObject<InternalSideFluxBase>("flux")),
61  _rho_var(coupled("rho")),
62  _rhou_var(coupled("rhou")),
63  _rhov_var(isCoupled("rhov") ? coupled("rhov") : zero),
64  _rhow_var(isCoupled("rhow") ? coupled("rhow") : zero),
65  _rhoe_var(coupled("rhoe"))
66 {
67  _jmap.insert(std::pair<unsigned int, unsigned int>(_rho_var, 0));
68  _jmap.insert(std::pair<unsigned int, unsigned int>(_rhou_var, 1));
69  _jmap.insert(std::pair<unsigned int, unsigned int>(_rhov_var, 2));
70  _jmap.insert(std::pair<unsigned int, unsigned int>(_rhow_var, 3));
71  _jmap.insert(std::pair<unsigned int, unsigned int>(_rhoe_var, 4));
72 }
73 
75 
76 Real
77 CNSFVKernel::computeQpResidual(Moose::DGResidualType type)
78 {
85 
86  std::vector<Real> uvec1 = {_rho1[_qp], _rhou1[_qp], _rhov1[_qp], _rhow1[_qp], _rhoe1[_qp]};
87 
88  std::vector<Real> uvec2 = {_rho2[_qp], _rhou2[_qp], _rhov2[_qp], _rhow2[_qp], _rhoe2[_qp]};
89 
91  const std::vector<Real> & flux = _flux.getFlux(
92  _current_side, _current_elem->id(), _neighbor_elem->id(), uvec1, uvec2, _normals[_qp], _tid);
93 
94  Real re = 0.;
95  switch (type)
96  {
97  case Moose::Element:
98  re = flux[_component] * _test[_i][_qp];
99  break;
100  case Moose::Neighbor:
101  re = -flux[_component] * _test_neighbor[_i][_qp];
102  break;
103  }
104  return re;
105 }
106 
107 Real
108 CNSFVKernel::computeQpJacobian(Moose::DGJacobianType type)
109 {
111 
112  std::vector<Real> uvec1 = {_rhoc1[_qp], _rhouc1[_qp], _rhovc1[_qp], _rhowc1[_qp], _rhoec1[_qp]};
113 
114  std::vector<Real> uvec2 = {_rhoc2[_qp], _rhouc2[_qp], _rhovc2[_qp], _rhowc2[_qp], _rhoec2[_qp]};
115 
116  const DenseMatrix<Real> & fjac1 = _flux.getJacobian(Moose::Element,
117  _current_side,
118  _current_elem->id(),
119  _neighbor_elem->id(),
120  uvec1,
121  uvec2,
122  _normals[_qp],
123  _tid);
124 
125  const DenseMatrix<Real> & fjac2 = _flux.getJacobian(Moose::Neighbor,
126  _current_side,
127  _current_elem->id(),
128  _neighbor_elem->id(),
129  uvec1,
130  uvec2,
131  _normals[_qp],
132  _tid);
133 
134  Real re = 0.;
135  switch (type)
136  {
137  case Moose::ElementElement:
138  re = fjac1(_component, _component) * _phi[_j][_qp] * _test[_i][_qp];
139  break;
140  case Moose::ElementNeighbor:
141  re = fjac2(_component, _component) * _phi_neighbor[_j][_qp] * _test[_i][_qp];
142  break;
143  case Moose::NeighborElement:
144  re = -fjac1(_component, _component) * _phi[_j][_qp] * _test_neighbor[_i][_qp];
145  break;
146  case Moose::NeighborNeighbor:
147  re = -fjac2(_component, _component) * _phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
148  break;
149  }
150  return re;
151 }
152 
153 Real
154 CNSFVKernel::computeQpOffDiagJacobian(Moose::DGJacobianType type, unsigned int jvar)
155 {
157 
158  std::vector<Real> uvec1 = {_rhoc1[_qp], _rhouc1[_qp], _rhovc1[_qp], _rhowc1[_qp], _rhoec1[_qp]};
159 
160  std::vector<Real> uvec2 = {_rhoc2[_qp], _rhouc2[_qp], _rhovc2[_qp], _rhowc2[_qp], _rhoec2[_qp]};
161 
162  const DenseMatrix<Real> & fjac1 = _flux.getJacobian(Moose::Element,
163  _current_side,
164  _current_elem->id(),
165  _neighbor_elem->id(),
166  uvec1,
167  uvec2,
168  _normals[_qp],
169  _tid);
170 
171  const DenseMatrix<Real> & fjac2 = _flux.getJacobian(Moose::Neighbor,
172  _current_side,
173  _current_elem->id(),
174  _neighbor_elem->id(),
175  uvec1,
176  uvec2,
177  _normals[_qp],
178  _tid);
179 
180  Real re = 0.;
181  switch (type)
182  {
183  case Moose::ElementElement:
184  re = fjac1(_component, _jmap[jvar]) * _phi[_j][_qp] * _test[_i][_qp];
185  break;
186  case Moose::ElementNeighbor:
187  re = fjac2(_component, _jmap[jvar]) * _phi_neighbor[_j][_qp] * _test[_i][_qp];
188  break;
189  case Moose::NeighborElement:
190  re = -fjac1(_component, _jmap[jvar]) * _phi[_j][_qp] * _test_neighbor[_i][_qp];
191  break;
192  case Moose::NeighborNeighbor:
193  re = -fjac2(_component, _jmap[jvar]) * _phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
194  break;
195  }
196  return re;
197 }
const MaterialProperty< Real > & _rho1
Definition: CNSFVKernel.h:71
const MaterialProperty< Real > & _rho2
Definition: CNSFVKernel.h:76
const VariableValue & _rhovc2
Definition: CNSFVKernel.h:66
virtual Real computeQpOffDiagJacobian(Moose::DGJacobianType type, unsigned int jvar)
Definition: CNSFVKernel.C:154
const VariableValue & _rhouc1
Definition: CNSFVKernel.h:60
const VariableValue & _rhoec2
Definition: CNSFVKernel.h:68
unsigned int _rho_var
Definition: CNSFVKernel.h:85
InputParameters validParams< CNSFVKernel >()
Definition: CNSFVKernel.C:12
virtual Real computeQpJacobian(Moose::DGJacobianType type)
Definition: CNSFVKernel.C:108
const MaterialProperty< Real > & _rhoe1
Definition: CNSFVKernel.h:75
const VariableValue & _rhowc2
Definition: CNSFVKernel.h:67
Real component(const SymmTensor &symm_tensor, unsigned int index)
const MaterialProperty< Real > & _rhou1
Definition: CNSFVKernel.h:72
unsigned int _rhov_var
Definition: CNSFVKernel.h:87
virtual const DenseMatrix< Real > & getJacobian(Moose::DGResidualType type, unsigned int iside, dof_id_type ielem, dof_id_type ineig, const std::vector< Real > &uvec1, const std::vector< Real > &uvec2, const RealVectorValue &dwave, THREAD_ID tid) const
Get the Jacobian matrix.
const VariableValue & _rhoc2
Definition: CNSFVKernel.h:64
const VariableValue & _rhouc2
Definition: CNSFVKernel.h:65
virtual ~CNSFVKernel()
Definition: CNSFVKernel.C:74
std::map< unsigned int, unsigned int > _jmap
Definition: CNSFVKernel.h:91
CNSFVKernel(const InputParameters &parameters)
Definition: CNSFVKernel.C:37
const MaterialProperty< Real > & _rhow1
Definition: CNSFVKernel.h:74
const VariableValue & _rhovc1
Definition: CNSFVKernel.h:61
virtual const std::vector< Real > & getFlux(unsigned int iside, dof_id_type ielem, dof_id_type ineig, const std::vector< Real > &uvec1, const std::vector< Real > &uvec2, const RealVectorValue &dwave, THREAD_ID tid) const
Get the flux vector.
MooseEnum _component
choose an equation
Definition: CNSFVKernel.h:53
unsigned int _rhou_var
Definition: CNSFVKernel.h:86
const VariableValue & _rhoc1
Definition: CNSFVKernel.h:59
const MaterialProperty< Real > & _rhow2
Definition: CNSFVKernel.h:79
A base class for computing and caching internal side flux.
const MaterialProperty< Real > & _rhou2
Definition: CNSFVKernel.h:77
const MaterialProperty< Real > & _rhov1
Definition: CNSFVKernel.h:73
const MaterialProperty< Real > & _rhov2
Definition: CNSFVKernel.h:78
const InternalSideFluxBase & _flux
fluid properties object
Definition: CNSFVKernel.h:83
const VariableValue & _rhowc1
Definition: CNSFVKernel.h:62
const VariableValue & _rhoec1
Definition: CNSFVKernel.h:63
virtual Real computeQpResidual(Moose::DGResidualType type)
Definition: CNSFVKernel.C:77
unsigned int _rhow_var
Definition: CNSFVKernel.h:88
unsigned int _rhoe_var
Definition: CNSFVKernel.h:89
const MaterialProperty< Real > & _rhoe2
Definition: CNSFVKernel.h:80