Line data Source code
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 "NavierStokesLHDGKernel.h"
11 : #include "MooseVariableFE.h"
12 : #include "MooseVariableScalar.h"
13 : #include "FEProblemBase.h"
14 :
15 : registerMooseObject("NavierStokesApp", NavierStokesLHDGKernel);
16 :
17 : InputParameters
18 640 : NavierStokesLHDGKernel::validParams()
19 : {
20 640 : auto params = HDGKernel::validParams();
21 640 : params += NavierStokesLHDGAssemblyHelper::validParams();
22 1280 : params.addParam<MooseFunctorName>(
23 1280 : "body_force_x", 0, "Body force for the momentum equation in the x-direction");
24 1280 : params.addParam<MooseFunctorName>(
25 1280 : "body_force_y", 0, "Body force for the momentum equation in the y-direction");
26 1280 : params.addParam<MooseFunctorName>(
27 1280 : "body_force_z", 0, "Body force for the momentum equation in the z-direction");
28 1280 : params.addParam<MooseFunctorName>(
29 : "pressure_mms_forcing_function",
30 1280 : 0,
31 : "A forcing function for the pressure (mass) equation for conducting MMS studies");
32 640 : params.addClassDescription("Implements the steady incompressible Navier-Stokes equations for a "
33 : "hybridized discretization");
34 1280 : params.renameParam("variable", "u", "The x-component of velocity");
35 640 : return params;
36 0 : }
37 :
38 320 : NavierStokesLHDGKernel::NavierStokesLHDGKernel(const InputParameters & parameters)
39 : : HDGKernel(parameters),
40 320 : NavierStokesLHDGAssemblyHelper(this, this, this, this, _fe_problem, _sys, _mesh, _tid),
41 : // body forces
42 320 : _body_force_x(getFunctor<Real>("body_force_x")),
43 640 : _body_force_y(getFunctor<Real>("body_force_y")),
44 640 : _body_force_z(getFunctor<Real>("body_force_z")),
45 640 : _pressure_mms_forcing_function(getFunctor<Real>("pressure_mms_forcing_function")),
46 320 : _qrule_face(_assembly.qRuleFace()),
47 320 : _q_point_face(_assembly.qPointsFace()),
48 320 : _JxW_face(_assembly.JxWFace()),
49 320 : _normals(_assembly.normals()),
50 320 : _current_side(_assembly.side())
51 : {
52 320 : if (_mesh.dimension() > 2)
53 0 : mooseError("This class only supports 2D simulations at this time");
54 :
55 320 : _body_forces.push_back(&_body_force_x);
56 320 : _body_forces.push_back(&_body_force_y);
57 320 : _body_forces.push_back(&_body_force_z);
58 320 : }
59 :
60 : void
61 1730496 : NavierStokesLHDGKernel::computeResidual()
62 : {
63 1730496 : _grad_u_vel_re.resize(_qu_dof_indices.size());
64 1730496 : _grad_v_vel_re.resize(_qv_dof_indices.size());
65 1730496 : _u_vel_re.resize(_u_dof_indices.size());
66 1730496 : _v_vel_re.resize(_v_dof_indices.size());
67 1730496 : _p_re.resize(_p_dof_indices.size());
68 1730496 : _global_lm_re.resize(_global_lm_dof_indices ? _global_lm_dof_indices->size() : 0);
69 :
70 : // qu and u
71 1730496 : vectorVolumeResidual(_qu_sol, _u_sol, _JxW, *_qrule, _grad_u_vel_re);
72 1730496 : scalarVolumeResidual(
73 1730496 : _qu_sol, 0, _body_force_x, _JxW, *_qrule, _current_elem, _q_point, _u_vel_re);
74 :
75 : // qv and v
76 1730496 : vectorVolumeResidual(_qv_sol, _v_sol, _JxW, *_qrule, _grad_v_vel_re);
77 1730496 : scalarVolumeResidual(
78 1730496 : _qv_sol, 1, _body_force_y, _JxW, *_qrule, _current_elem, _q_point, _v_vel_re);
79 :
80 : // p
81 1730496 : pressureVolumeResidual(
82 1730496 : _pressure_mms_forcing_function, _JxW, *_qrule, _current_elem, _q_point, _p_re, _global_lm_re);
83 :
84 1730496 : addResiduals(_assembly, _grad_u_vel_re, _qu_dof_indices, _grad_u_var.scalingFactor());
85 1730496 : addResiduals(_assembly, _grad_v_vel_re, _qv_dof_indices, _grad_v_var.scalingFactor());
86 1730496 : addResiduals(_assembly, _u_vel_re, _u_dof_indices, _u_var.scalingFactor());
87 1730496 : addResiduals(_assembly, _v_vel_re, _v_dof_indices, _v_var.scalingFactor());
88 1730496 : addResiduals(_assembly, _p_re, _p_dof_indices, _pressure_var.scalingFactor());
89 1730496 : if (_global_lm_dof_indices)
90 831152 : addResiduals(
91 831152 : _assembly, _global_lm_re, *_global_lm_dof_indices, _enclosure_lm_var->scalingFactor());
92 1730496 : }
93 :
94 : void
95 492856 : NavierStokesLHDGKernel::computeJacobian()
96 : {
97 492856 : _grad_u_grad_u_jac.resize(_qu_dof_indices.size(), _qu_dof_indices.size());
98 492856 : _grad_v_grad_v_jac.resize(_qv_dof_indices.size(), _qv_dof_indices.size());
99 492856 : _grad_u_u_jac.resize(_qu_dof_indices.size(), _u_dof_indices.size());
100 492856 : _grad_v_v_jac.resize(_qv_dof_indices.size(), _v_dof_indices.size());
101 492856 : _u_grad_u_jac.resize(_u_dof_indices.size(), _qu_dof_indices.size());
102 492856 : _v_grad_v_jac.resize(_v_dof_indices.size(), _qv_dof_indices.size());
103 492856 : _u_u_jac.resize(_u_dof_indices.size(), _u_dof_indices.size());
104 492856 : _u_v_jac.resize(_u_dof_indices.size(), _v_dof_indices.size());
105 492856 : _v_u_jac.resize(_v_dof_indices.size(), _u_dof_indices.size());
106 492856 : _v_v_jac.resize(_v_dof_indices.size(), _v_dof_indices.size());
107 492856 : _u_p_jac.resize(_u_dof_indices.size(), _p_dof_indices.size());
108 492856 : _v_p_jac.resize(_v_dof_indices.size(), _p_dof_indices.size());
109 492856 : _p_u_jac.resize(_p_dof_indices.size(), _u_dof_indices.size());
110 492856 : _p_v_jac.resize(_p_dof_indices.size(), _v_dof_indices.size());
111 492856 : _p_global_lm_jac.resize(_p_dof_indices.size(),
112 492856 : _global_lm_dof_indices ? _global_lm_dof_indices->size() : 0);
113 492856 : _global_lm_p_jac.resize(_global_lm_dof_indices ? _global_lm_dof_indices->size() : 0,
114 492856 : _p_dof_indices.size());
115 :
116 : // qu and u
117 492856 : vectorVolumeJacobian(_JxW, *_qrule, _grad_u_grad_u_jac, _grad_u_u_jac);
118 492856 : scalarVolumeJacobian(0, _JxW, *_qrule, _u_grad_u_jac, _u_u_jac, _u_v_jac, _u_p_jac);
119 :
120 : // qv and v
121 492856 : vectorVolumeJacobian(_JxW, *_qrule, _grad_v_grad_v_jac, _grad_v_v_jac);
122 492856 : scalarVolumeJacobian(1, _JxW, *_qrule, _v_grad_v_jac, _v_u_jac, _v_v_jac, _v_p_jac);
123 :
124 : // p
125 492856 : pressureVolumeJacobian(_JxW, *_qrule, _p_u_jac, _p_v_jac, _p_global_lm_jac, _global_lm_p_jac);
126 :
127 492856 : addJacobian(
128 492856 : _assembly, _grad_u_grad_u_jac, _qu_dof_indices, _qu_dof_indices, _grad_u_var.scalingFactor());
129 492856 : addJacobian(
130 492856 : _assembly, _grad_v_grad_v_jac, _qv_dof_indices, _qv_dof_indices, _grad_v_var.scalingFactor());
131 492856 : addJacobian(
132 492856 : _assembly, _grad_u_u_jac, _qu_dof_indices, _u_dof_indices, _grad_u_var.scalingFactor());
133 492856 : addJacobian(
134 492856 : _assembly, _grad_v_v_jac, _qv_dof_indices, _v_dof_indices, _grad_v_var.scalingFactor());
135 492856 : addJacobian(_assembly, _u_grad_u_jac, _u_dof_indices, _qu_dof_indices, _u_var.scalingFactor());
136 492856 : addJacobian(_assembly, _v_grad_v_jac, _v_dof_indices, _qv_dof_indices, _v_var.scalingFactor());
137 492856 : addJacobian(_assembly, _u_u_jac, _u_dof_indices, _u_dof_indices, _u_var.scalingFactor());
138 492856 : addJacobian(_assembly, _u_v_jac, _u_dof_indices, _v_dof_indices, _u_var.scalingFactor());
139 492856 : addJacobian(_assembly, _v_u_jac, _v_dof_indices, _u_dof_indices, _v_var.scalingFactor());
140 492856 : addJacobian(_assembly, _v_v_jac, _v_dof_indices, _v_dof_indices, _v_var.scalingFactor());
141 492856 : addJacobian(_assembly, _u_p_jac, _u_dof_indices, _p_dof_indices, _u_var.scalingFactor());
142 492856 : addJacobian(_assembly, _v_p_jac, _v_dof_indices, _p_dof_indices, _v_var.scalingFactor());
143 492856 : addJacobian(_assembly, _p_u_jac, _p_dof_indices, _u_dof_indices, _pressure_var.scalingFactor());
144 492856 : addJacobian(_assembly, _p_v_jac, _p_dof_indices, _v_dof_indices, _pressure_var.scalingFactor());
145 492856 : if (_global_lm_dof_indices)
146 : {
147 217408 : addJacobian(_assembly,
148 : _p_global_lm_jac,
149 : _p_dof_indices,
150 : *_global_lm_dof_indices,
151 217408 : _pressure_var.scalingFactor());
152 217408 : addJacobian(_assembly,
153 : _global_lm_p_jac,
154 217408 : *_global_lm_dof_indices,
155 : _p_dof_indices,
156 217408 : _enclosure_lm_var->scalingFactor());
157 : }
158 492856 : }
159 :
160 : void
161 5100416 : NavierStokesLHDGKernel::computeResidualOnSide()
162 : {
163 5100416 : const Elem * const neigh = _current_elem->neighbor_ptr(_current_side);
164 :
165 5100416 : _grad_u_vel_re.resize(_qu_dof_indices.size());
166 5100416 : _u_vel_re.resize(_u_dof_indices.size());
167 5100416 : _lm_u_vel_re.resize(_lm_u_dof_indices.size());
168 5100416 : _grad_v_vel_re.resize(_qv_dof_indices.size());
169 5100416 : _v_vel_re.resize(_v_dof_indices.size());
170 5100416 : _lm_v_vel_re.resize(_lm_v_dof_indices.size());
171 5100416 : _p_re.resize(_p_dof_indices.size());
172 :
173 : // qu, u, lm_u
174 5100416 : vectorFaceResidual(_lm_u_sol, _JxW_face, *_qrule_face, _normals, _grad_u_vel_re);
175 5100416 : scalarFaceResidual(_qu_sol, _u_sol, _lm_u_sol, 0, _JxW_face, *_qrule_face, _normals, _u_vel_re);
176 5100416 : lmFaceResidual(
177 5100416 : _qu_sol, _u_sol, _lm_u_sol, 0, _JxW_face, *_qrule_face, _normals, neigh, _lm_u_vel_re);
178 :
179 : // qv, v, lm_v
180 5100416 : vectorFaceResidual(_lm_v_sol, _JxW_face, *_qrule_face, _normals, _grad_v_vel_re);
181 5100416 : scalarFaceResidual(_qv_sol, _v_sol, _lm_v_sol, 1, _JxW_face, *_qrule_face, _normals, _v_vel_re);
182 5100416 : lmFaceResidual(
183 5100416 : _qv_sol, _v_sol, _lm_v_sol, 1, _JxW_face, *_qrule_face, _normals, neigh, _lm_v_vel_re);
184 :
185 : // p
186 5100416 : pressureFaceResidual(_JxW_face, *_qrule_face, _normals, _p_re);
187 :
188 5100416 : addResiduals(_assembly, _grad_u_vel_re, _qu_dof_indices, _grad_u_var.scalingFactor());
189 5100416 : addResiduals(_assembly, _u_vel_re, _u_dof_indices, _u_var.scalingFactor());
190 5100416 : addResiduals(_assembly, _lm_u_vel_re, _lm_u_dof_indices, _u_face_var.scalingFactor());
191 5100416 : addResiduals(_assembly, _grad_v_vel_re, _qv_dof_indices, _grad_v_var.scalingFactor());
192 5100416 : addResiduals(_assembly, _v_vel_re, _v_dof_indices, _v_var.scalingFactor());
193 5100416 : addResiduals(_assembly, _lm_v_vel_re, _lm_v_dof_indices, _v_face_var.scalingFactor());
194 5100416 : addResiduals(_assembly, _p_re, _p_dof_indices, _pressure_var.scalingFactor());
195 5100416 : }
196 :
197 : void
198 1452352 : NavierStokesLHDGKernel::computeJacobianOnSide()
199 : {
200 1452352 : const Elem * const neigh = _current_elem->neighbor_ptr(_current_side);
201 :
202 1452352 : _grad_u_lm_u_jac.resize(_qu_dof_indices.size(), _lm_u_dof_indices.size());
203 1452352 : _u_grad_u_jac.resize(_u_dof_indices.size(), _qu_dof_indices.size());
204 1452352 : _u_u_jac.resize(_u_dof_indices.size(), _u_dof_indices.size());
205 1452352 : _u_lm_u_jac.resize(_u_dof_indices.size(), _lm_u_dof_indices.size());
206 1452352 : _u_lm_v_jac.resize(_u_dof_indices.size(), _lm_v_dof_indices.size());
207 1452352 : _u_p_jac.resize(_u_dof_indices.size(), _p_dof_indices.size());
208 1452352 : _lm_u_grad_u_jac.resize(_lm_u_dof_indices.size(), _qu_dof_indices.size());
209 1452352 : _lm_u_u_jac.resize(_lm_u_dof_indices.size(), _u_dof_indices.size());
210 1452352 : _lm_u_lm_u_jac.resize(_lm_u_dof_indices.size(), _lm_u_dof_indices.size());
211 1452352 : _lm_u_lm_v_jac.resize(_lm_u_dof_indices.size(), _lm_v_dof_indices.size());
212 1452352 : _lm_u_p_jac.resize(_lm_u_dof_indices.size(), _p_dof_indices.size());
213 1452352 : _grad_v_lm_v_jac.resize(_qv_dof_indices.size(), _lm_v_dof_indices.size());
214 1452352 : _v_grad_v_jac.resize(_v_dof_indices.size(), _qv_dof_indices.size());
215 1452352 : _v_v_jac.resize(_v_dof_indices.size(), _v_dof_indices.size());
216 1452352 : _v_lm_u_jac.resize(_v_dof_indices.size(), _lm_u_dof_indices.size());
217 1452352 : _v_lm_v_jac.resize(_v_dof_indices.size(), _lm_v_dof_indices.size());
218 1452352 : _v_p_jac.resize(_v_dof_indices.size(), _p_dof_indices.size());
219 1452352 : _lm_v_grad_v_jac.resize(_lm_v_dof_indices.size(), _qv_dof_indices.size());
220 1452352 : _lm_v_v_jac.resize(_lm_v_dof_indices.size(), _v_dof_indices.size());
221 1452352 : _lm_v_lm_u_jac.resize(_lm_v_dof_indices.size(), _lm_u_dof_indices.size());
222 1452352 : _lm_v_lm_v_jac.resize(_lm_v_dof_indices.size(), _lm_v_dof_indices.size());
223 1452352 : _lm_v_p_jac.resize(_lm_v_dof_indices.size(), _p_dof_indices.size());
224 1452352 : _p_lm_u_jac.resize(_p_dof_indices.size(), _lm_u_dof_indices.size());
225 1452352 : _p_lm_v_jac.resize(_p_dof_indices.size(), _lm_v_dof_indices.size());
226 :
227 : // qu, u, lm_u
228 1452352 : vectorFaceJacobian(_JxW_face, *_qrule_face, _normals, _grad_u_lm_u_jac);
229 1452352 : scalarFaceJacobian(0,
230 : _JxW_face,
231 1452352 : *_qrule_face,
232 : _normals,
233 : _u_grad_u_jac,
234 : _u_u_jac,
235 : _u_lm_u_jac,
236 : _u_p_jac,
237 : _u_lm_u_jac,
238 : _u_lm_v_jac);
239 1452352 : lmFaceJacobian(0,
240 : _JxW_face,
241 1452352 : *_qrule_face,
242 : _normals,
243 : neigh,
244 : _lm_u_grad_u_jac,
245 : _lm_u_u_jac,
246 : _lm_u_lm_u_jac,
247 : _lm_u_p_jac,
248 : _lm_u_lm_u_jac,
249 : _lm_u_lm_v_jac);
250 :
251 : // qv, v, lm_v
252 1452352 : vectorFaceJacobian(_JxW_face, *_qrule_face, _normals, _grad_v_lm_v_jac);
253 1452352 : scalarFaceJacobian(1,
254 : _JxW_face,
255 1452352 : *_qrule_face,
256 : _normals,
257 : _v_grad_v_jac,
258 : _v_v_jac,
259 : _v_lm_v_jac,
260 : _v_p_jac,
261 : _v_lm_u_jac,
262 : _v_lm_v_jac);
263 1452352 : lmFaceJacobian(1,
264 : _JxW_face,
265 1452352 : *_qrule_face,
266 : _normals,
267 : neigh,
268 : _lm_v_grad_v_jac,
269 : _lm_v_v_jac,
270 : _lm_v_lm_v_jac,
271 : _lm_v_p_jac,
272 : _lm_v_lm_u_jac,
273 : _lm_v_lm_v_jac);
274 :
275 : // p
276 1452352 : pressureFaceJacobian(_JxW_face, *_qrule_face, _normals, _p_lm_u_jac, _p_lm_v_jac);
277 :
278 1452352 : addJacobian(
279 1452352 : _assembly, _grad_u_lm_u_jac, _qu_dof_indices, _lm_u_dof_indices, _grad_u_var.scalingFactor());
280 1452352 : addJacobian(_assembly, _u_grad_u_jac, _u_dof_indices, _qu_dof_indices, _u_var.scalingFactor());
281 1452352 : addJacobian(_assembly, _u_u_jac, _u_dof_indices, _u_dof_indices, _u_var.scalingFactor());
282 1452352 : addJacobian(_assembly, _u_lm_u_jac, _u_dof_indices, _lm_u_dof_indices, _u_var.scalingFactor());
283 1452352 : addJacobian(_assembly, _u_lm_v_jac, _u_dof_indices, _lm_v_dof_indices, _u_var.scalingFactor());
284 1452352 : addJacobian(_assembly, _u_p_jac, _u_dof_indices, _p_dof_indices, _u_var.scalingFactor());
285 1452352 : addJacobian(
286 1452352 : _assembly, _lm_u_grad_u_jac, _lm_u_dof_indices, _qu_dof_indices, _u_face_var.scalingFactor());
287 1452352 : addJacobian(
288 1452352 : _assembly, _lm_u_u_jac, _lm_u_dof_indices, _u_dof_indices, _u_face_var.scalingFactor());
289 1452352 : addJacobian(
290 1452352 : _assembly, _lm_u_lm_u_jac, _lm_u_dof_indices, _lm_u_dof_indices, _u_face_var.scalingFactor());
291 1452352 : addJacobian(
292 1452352 : _assembly, _lm_u_lm_v_jac, _lm_u_dof_indices, _lm_v_dof_indices, _u_face_var.scalingFactor());
293 1452352 : addJacobian(
294 1452352 : _assembly, _lm_u_p_jac, _lm_u_dof_indices, _p_dof_indices, _u_face_var.scalingFactor());
295 1452352 : addJacobian(
296 1452352 : _assembly, _grad_v_lm_v_jac, _qv_dof_indices, _lm_v_dof_indices, _grad_v_var.scalingFactor());
297 1452352 : addJacobian(_assembly, _v_grad_v_jac, _v_dof_indices, _qv_dof_indices, _v_var.scalingFactor());
298 1452352 : addJacobian(_assembly, _v_v_jac, _v_dof_indices, _v_dof_indices, _v_var.scalingFactor());
299 1452352 : addJacobian(_assembly, _v_lm_u_jac, _v_dof_indices, _lm_u_dof_indices, _v_var.scalingFactor());
300 1452352 : addJacobian(_assembly, _v_lm_v_jac, _v_dof_indices, _lm_v_dof_indices, _v_var.scalingFactor());
301 1452352 : addJacobian(_assembly, _v_p_jac, _v_dof_indices, _p_dof_indices, _v_var.scalingFactor());
302 1452352 : addJacobian(
303 1452352 : _assembly, _lm_v_grad_v_jac, _lm_v_dof_indices, _qv_dof_indices, _v_face_var.scalingFactor());
304 1452352 : addJacobian(
305 1452352 : _assembly, _lm_v_v_jac, _lm_v_dof_indices, _v_dof_indices, _v_face_var.scalingFactor());
306 1452352 : addJacobian(
307 1452352 : _assembly, _lm_v_lm_u_jac, _lm_v_dof_indices, _lm_u_dof_indices, _v_face_var.scalingFactor());
308 1452352 : addJacobian(
309 1452352 : _assembly, _lm_v_lm_v_jac, _lm_v_dof_indices, _lm_v_dof_indices, _v_face_var.scalingFactor());
310 1452352 : addJacobian(
311 1452352 : _assembly, _lm_v_p_jac, _lm_v_dof_indices, _p_dof_indices, _v_face_var.scalingFactor());
312 1452352 : addJacobian(
313 1452352 : _assembly, _p_lm_u_jac, _p_dof_indices, _lm_u_dof_indices, _pressure_var.scalingFactor());
314 1452352 : addJacobian(
315 1452352 : _assembly, _p_lm_v_jac, _p_dof_indices, _lm_v_dof_indices, _pressure_var.scalingFactor());
316 1452352 : }
317 :
318 : void
319 320 : NavierStokesLHDGKernel::initialSetup()
320 : {
321 : // This check must occur after FEProblemBase::init()
322 320 : checkCoupling();
323 320 : }
324 :
325 : void
326 1171 : NavierStokesLHDGKernel::jacobianSetup()
327 : {
328 1171 : _cached_elem = nullptr;
329 1171 : }
330 :
331 : void
332 3449992 : NavierStokesLHDGKernel::computeOffDiagJacobian(const unsigned int)
333 : {
334 3449992 : if (_cached_elem != _current_elem)
335 : {
336 492856 : computeJacobian();
337 492856 : _cached_elem = _current_elem;
338 : }
339 3449992 : }
340 :
341 : std::set<std::string>
342 320 : NavierStokesLHDGKernel::additionalROVariables()
343 : {
344 320 : std::set<std::string> covered_vars = {_grad_u_var.name(),
345 320 : _u_face_var.name(),
346 320 : _v_var.name(),
347 320 : _grad_v_var.name(),
348 320 : _v_face_var.name(),
349 2240 : _pressure_var.name()};
350 320 : if (_enclosure_lm_var)
351 160 : covered_vars.insert(_enclosure_lm_var->name());
352 320 : return covered_vars;
353 320 : }
|
