www.mooseframework.org
NavierStokesMaterial.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 // Navier-Stokes includes
9 #include "NavierStokesMaterial.h"
10 #include "NS.h"
11 
12 // FluidProperties includes
14 
15 // MOOSE includes
16 #include "Assembly.h"
17 #include "MooseMesh.h"
18 
19 #include "libmesh/quadrature.h"
20 
21 template <>
22 InputParameters
24 {
25  InputParameters params = validParams<Material>();
26 
27  params.addClassDescription("This is the base class all materials should use if you are trying to "
28  "use the Navier-Stokes Kernels.");
29  params.addRequiredCoupledVar(NS::velocity_x, "x-velocity");
30  params.addCoupledVar(NS::velocity_y, "y-velocity"); // only required in >= 2D
31  params.addCoupledVar(NS::velocity_z, "z-velocity"); // only required in 3D
32 
33  params.addRequiredCoupledVar(NS::temperature, "temperature");
34  params.addRequiredCoupledVar(NS::enthalpy, "total enthalpy");
35 
36  params.addRequiredCoupledVar(NS::density, "density");
37  params.addRequiredCoupledVar(NS::momentum_x, "x-momentum");
38  params.addCoupledVar(NS::momentum_y, "y-momentum"); // only required in >= 2D
39  params.addCoupledVar(NS::momentum_z, "z-momentum"); // only required in 3D
40  params.addRequiredCoupledVar(NS::total_energy, "energy");
41  params.addRequiredParam<UserObjectName>("fluid_properties",
42  "The name of the user object for fluid properties");
43 
44  return params;
45 }
46 
47 NavierStokesMaterial::NavierStokesMaterial(const InputParameters & parameters)
48  : Material(parameters),
49  _mesh_dimension(_mesh.dimension()),
50  _grad_u(coupledGradient(NS::velocity_x)),
51  _grad_v(_mesh_dimension >= 2 ? coupledGradient(NS::velocity_y) : _grad_zero),
52  _grad_w(_mesh_dimension == 3 ? coupledGradient(NS::velocity_z) : _grad_zero),
53 
54  _viscous_stress_tensor(declareProperty<RealTensorValue>("viscous_stress_tensor")),
55  _thermal_conductivity(declareProperty<Real>("thermal_conductivity")),
56 
57  // Declared here but _not_ calculated here
58  // (See e.g. derived class, bighorn/include/materials/FluidTC1.h)
59  _dynamic_viscosity(declareProperty<Real>("dynamic_viscosity")),
60 
61  // The momentum components of the inviscid flux Jacobians.
62  _calA(declareProperty<std::vector<RealTensorValue>>("calA")),
63 
64  // "Velocity column" matrices
65  _calC(declareProperty<std::vector<RealTensorValue>>("calC")),
66 
67  // Energy equation inviscid flux matrices, "cal E_{kl}" in the notes.
68  _calE(declareProperty<std::vector<std::vector<RealTensorValue>>>("calE")),
69  _vel_grads({&_grad_u, &_grad_v, &_grad_w}),
70 
71  // Coupled solution values needed for computing SUPG stabilization terms
72  _u_vel(coupledValue(NS::velocity_x)),
73  _v_vel(_mesh.dimension() >= 2 ? coupledValue(NS::velocity_y) : _zero),
74  _w_vel(_mesh.dimension() == 3 ? coupledValue(NS::velocity_z) : _zero),
75 
76  _temperature(coupledValue(NS::temperature)),
77  _enthalpy(coupledValue(NS::enthalpy)),
78 
79  // Coupled solution values
80  _rho(coupledValue(NS::density)),
81  _rho_u(coupledValue(NS::momentum_x)),
82  _rho_v(_mesh.dimension() >= 2 ? coupledValue(NS::momentum_y) : _zero),
83  _rho_w(_mesh.dimension() == 3 ? coupledValue(NS::momentum_z) : _zero),
84  _rho_E(coupledValue(NS::total_energy)),
85 
86  // Time derivative values
87  _drho_dt(coupledDot(NS::density)),
88  _drhou_dt(coupledDot(NS::momentum_x)),
89  _drhov_dt(_mesh.dimension() >= 2 ? coupledDot(NS::momentum_y) : _zero),
90  _drhow_dt(_mesh.dimension() == 3 ? coupledDot(NS::momentum_z) : _zero),
91  _drhoE_dt(coupledDot(NS::total_energy)),
92 
93  // Gradients
94  _grad_rho(coupledGradient(NS::density)),
95  _grad_rho_u(coupledGradient(NS::momentum_x)),
96  _grad_rho_v(_mesh.dimension() >= 2 ? coupledGradient(NS::momentum_y) : _grad_zero),
97  _grad_rho_w(_mesh.dimension() == 3 ? coupledGradient(NS::momentum_z) : _grad_zero),
98  _grad_rho_E(coupledGradient(NS::total_energy)),
99 
100  // Material properties for stabilization
101  _hsupg(declareProperty<Real>("hsupg")),
102  _tauc(declareProperty<Real>("tauc")),
103  _taum(declareProperty<Real>("taum")),
104  _taue(declareProperty<Real>("taue")),
105  _strong_residuals(declareProperty<std::vector<Real>>("strong_residuals")),
106  _fp(getUserObject<IdealGasFluidProperties>("fluid_properties"))
107 {
108 }
109 
113 void
115 {
116  for (unsigned int qp = 0; qp < _qrule->n_points(); ++qp)
117  {
118  /******* Viscous Stress Tensor *******/
119  // Technically... this _is_ the transpose (since we are loading these by rows)
120  // But it doesn't matter....
121  RealTensorValue grad_outer_u(_grad_u[qp], _grad_v[qp], _grad_w[qp]);
122 
123  grad_outer_u += grad_outer_u.transpose();
124 
125  Real div_vel = 0.0;
126  for (unsigned int i = 0; i < 3; ++i)
127  div_vel += (*_vel_grads[i])[qp](i);
128 
129  // Add diagonal terms
130  for (unsigned int i = 0; i < 3; ++i)
131  grad_outer_u(i, i) -= 2.0 / 3.0 * div_vel;
132 
133  grad_outer_u *= _dynamic_viscosity[qp];
134 
135  _viscous_stress_tensor[qp] = grad_outer_u;
136 
137  // Tabulated values of thermal conductivity vs. Temperature for air (k increases slightly with
138  // T):
139  // T (K) k (W/m-K)
140  // 273 0.0243
141  // 373 0.0314
142  // 473 0.0386
143  // 573 0.0454
144  // 673 0.0515
145 
146  // Pr = (mu * cp) / k ==> k = (mu * cp) / Pr = (mu * gamma * cv) / Pr.
147  // TODO: We are using a fixed value of the Prandtl number which is
148  // valid for air, it may or may not depend on temperature? Since
149  // this is a property of the fluid, it could possibly be moved to
150  // the FluidProperties module...
151  const Real Pr = 0.71;
152  _thermal_conductivity[qp] = (_dynamic_viscosity[qp] * _fp.cp()) / Pr;
153 
154  // Compute stabilization parameters:
155 
156  // .) Compute SUPG element length scale.
157  computeHSUPG(qp);
158  // Moose::out << "_hsupg[" << qp << "]=" << _hsupg[qp] << std::endl;
159 
160  // .) Compute SUPG parameter values. (Must call this after computeHSUPG())
161  computeTau(qp);
162  // Moose::out << "_tauc[" << qp << "]=" << _tauc[qp] << ", ";
163  // Moose::out << "_taum[" << qp << "]=" << _taum[qp] << ", ";
164  // Moose::out << "_taue[" << qp << "]=" << _taue[qp] << std::endl;
165 
166  // .) Compute strong residual values.
168  // Moose::out << "_strong_residuals[" << qp << "]=";
169  // for (unsigned i=0; i<_strong_residuals[qp].size(); ++i)
170  // Moose::out << _strong_residuals[qp][i] << " ";
171  // Moose::out << std::endl;
172  }
173 }
174 
175 void
177 {
178  // // Grab reference to linear Lagrange finite element object pointer,
179  // // currently this is always a linear Lagrange element, so this might need to
180  // // be generalized if we start working with higher-order elements...
181  // FEBase*& fe(_assembly.getFE(FEType(), _current_elem->dim()));
182  //
183  // // Grab references to FE object's mapping data from the _subproblem's FE object
184  // const std::vector<Real> & dxidx(fe->get_dxidx());
185  // const std::vector<Real> & dxidy(fe->get_dxidy());
186  // const std::vector<Real> & dxidz(fe->get_dxidz());
187  // const std::vector<Real> & detadx(fe->get_detadx());
188  // const std::vector<Real> & detady(fe->get_detady());
189  // const std::vector<Real> & detadz(fe->get_detadz());
190  // const std::vector<Real> & dzetadx(fe->get_dzetadx()); // Empty in 2D
191  // const std::vector<Real> & dzetady(fe->get_dzetady()); // Empty in 2D
192  // const std::vector<Real> & dzetadz(fe->get_dzetadz()); // Empty in 2D
193  //
194  // // Bounds checking on element data
195  // mooseAssert(qp < dxidx.size(), "Insufficient data in dxidx array!");
196  // mooseAssert(qp < dxidy.size(), "Insufficient data in dxidy array!");
197  // mooseAssert(qp < dxidz.size(), "Insufficient data in dxidz array!");
198  //
199  // mooseAssert(qp < detadx.size(), "Insufficient data in detadx array!");
200  // mooseAssert(qp < detady.size(), "Insufficient data in detady array!");
201  // mooseAssert(qp < detadz.size(), "Insufficient data in detadz array!");
202  //
203  // if (_mesh_dimension == 3)
204  // {
205  // mooseAssert(qp < dzetadx.size(), "Insufficient data in dzetadx array!");
206  // mooseAssert(qp < dzetady.size(), "Insufficient data in dzetady array!");
207  // mooseAssert(qp < dzetadz.size(), "Insufficient data in dzetadz array!");
208  // }
209  //
210  // // The velocity vector at this quadrature point.
211  // RealVectorValue U(_u_vel[qp],_v_vel[qp],_w_vel[qp]);
212  //
213  // // Pull out element inverse map values at the current qp into a little dense matrix
214  // Real dxi_dx[3][3] = {{0.,0.,0.}, {0.,0.,0.}, {0.,0.,0.}};
215  //
216  // dxi_dx[0][0] = dxidx[qp]; dxi_dx[0][1] = dxidy[qp];
217  // dxi_dx[1][0] = detadx[qp]; dxi_dx[1][1] = detady[qp];
218  //
219  // // OK to access third entries on 2D elements if LIBMESH_DIM==3, though they
220  // // may be zero...
221  // if (LIBMESH_DIM == 3)
222  // {
223  // /**/ /**/ dxi_dx[0][2] = dxidz[qp];
224  // /**/ /**/ dxi_dx[1][2] = detadz[qp];
225  // }
226  //
227  // // The last row of entries available only for 3D elements.
228  // if (_mesh_dimension == 3)
229  // {
230  // dxi_dx[2][0] = dzetadx[qp]; dxi_dx[2][1] = dzetady[qp]; dxi_dx[2][2] = dzetadz[qp];
231  // }
232  //
233  // // Construct the g_ij = d(xi_k)/d(x_j) * d(xi_k)/d(x_i) matrix
234  // // from Ben and Bova's paper by summing over k...
235  // Real g[3][3] = {{0.,0.,0.}, {0.,0.,0.}, {0.,0.,0.}};
236  // for (unsigned int i = 0; i < 3; ++i)
237  // for (unsigned int j = 0; j < 3; ++j)
238  // for (unsigned int k = 0; k < 3; ++k)
239  // g[i][j] += dxi_dx[k][j] * dxi_dx[k][i];
240  //
241  // // Compute the denominator of the h_supg term: U * (g) * U
242  // Real denom = 0.;
243  // for (unsigned int i = 0; i < 3; ++i)
244  // for (unsigned int j = 0; j < 3; ++j)
245  // denom += U(j) * g[i][j] * U(i);
246  //
247  // // Compute h_supg. Some notes:
248  // // .) The 2 coefficient in this term should be a 1 if we are using tets/triangles.
249  // // .) The denominator will be identically zero only if the velocity
250  // // is identically zero, in which case we can't divide by it.
251  // if (denom != 0.0)
252  // _hsupg[qp] = 2.* sqrt( U.norm_sq() / denom );
253  // else
254  // _hsupg[qp] = 0.;
255 
256  // Simple (and fast) implementation: Just use hmin for the element!
257  _hsupg[qp] = _current_elem->hmin();
258 }
259 
260 void
262 {
263  Real velmag =
264  std::sqrt(_u_vel[qp] * _u_vel[qp] + _v_vel[qp] * _v_vel[qp] + _w_vel[qp] * _w_vel[qp]);
265 
266  // Moose::out << "velmag=" << velmag << std::endl;
267 
268  // Make sure temperature >= 0 before trying to take sqrt
269  // if (_temperature[qp] < 0.)
270  // {
271  // Moose::err << "Negative temperature "
272  // << _temperature[qp]
273  // << " found at quadrature point "
274  // << qp
275  // << ", element "
276  // << _current_elem->id()
277  // << std::endl;
278  // mooseError("Can't continue, would be nice to throw an exception here?");
279  // }
280 
281  // The speed of sound for an ideal gas, sqrt(gamma * R * T). Not needed unless
282  // we want to use a form of Tau that requires it.
283  // Real soundspeed = _fp.c(_specific_volume[_qp], _internal_energy[_qp]);
284 
285  // If velmag == 0, then _hsupg should be zero as well. Then tau
286  // will have only the time-derivative contribution (or zero, if we
287  // are not including dt terms in our taus!) Note that using the
288  // time derivative contribution in this way assumes we are solving
289  // unsteady, and guarantees *some* stabilization is added even when
290  // u -> 0 in certain regions of the flow.
291  if (velmag == 0.)
292  {
293  // 1.) Tau without dt terms
294  // _tauc[qp] = 0.;
295  // _taum[qp] = 0.;
296  // _taue[qp] = 0.;
297 
298  // 2.) Tau *with* dt terms
299  _tauc[qp] = _taum[qp] = _taue[qp] = 0.5 * _dt;
300  }
301  else
302  {
303  // The element length parameter, squared
304  Real h2 = _hsupg[qp] * _hsupg[qp];
305 
306  // The viscosity-based term
307  Real visc_term = _dynamic_viscosity[qp] / _rho[qp] / h2;
308 
309  // The thermal conductivity-based term, cp = gamma * cv
310  Real k_term = _thermal_conductivity[qp] / _rho[qp] / _fp.cp() / h2;
311 
312  // 1a.) Standard compressible flow tau. Does not account for low Mach number
313  // limit.
314  // _tauc[qp] = _hsupg[qp] / (velmag + soundspeed);
315 
316  // 1b.) Inspired by Hauke, the sum of the compressible and incompressible tauc.
317  // _tauc[qp] =
318  // _hsupg[qp] / (velmag + soundspeed) +
319  // _hsupg[qp] / (velmag);
320 
321  // 1c.) From Wong 2001. This tau is O(M^2) for small M. At small M,
322  // tauc dominates the inverse square sums and basically makes
323  // taum=taue=tauc. However, all my flows occur at low Mach numbers,
324  // so there would basically never be any stabilization...
325  // _tauc[qp] = (_hsupg[qp] * velmag) / (velmag*velmag + soundspeed*soundspeed);
326 
327  // For use with option "1",
328  // (tau_c)^{-2}
329  // Real taucm2 = 1./_tauc[qp]/_tauc[qp];
330  // _taum[qp] = 1. / std::sqrt(taucm2 + visc_term*visc_term);
331  // _taue[qp] = 1. / std::sqrt(taucm2 + k_term*k_term);
332 
333  // 2.) Tau with timestep dependence (guarantees stabilization even
334  // in zero-velocity limit) incorporated via the "r-switch" method,
335  // with r=2.
336  Real sqrt_term = 4. / _dt / _dt + velmag * velmag / h2;
337 
338  // For use with option "2", i.e. the option that uses dt in the definition of tau
339  _tauc[qp] = 1. / std::sqrt(sqrt_term);
340  _taum[qp] = 1. / std::sqrt(sqrt_term + visc_term * visc_term);
341  _taue[qp] = 1. / std::sqrt(sqrt_term + k_term * k_term);
342  }
343 
344  // Debugging
345  // Moose::out << "_tauc[" << qp << "]=" << _tauc[qp] << std::endl;
346  // Moose::out << "_hsupg[" << qp << "]=" << _hsupg[qp] << std::endl;
347  // Moose::out << "velmag[" << qp << "]=" << velmag << std::endl;
348 }
349 
350 void
352 {
353  // Create storage at this qp for the strong residuals of all the equations.
354  // In 2D, the value for the z-velocity equation will just be zero.
355  _strong_residuals[qp].resize(5);
356 
357  // The timestep is stored in the Problem object, which can be accessed through
358  // the parent pointer of the SubProblem. Don't need this if we are not
359  // approximating time derivatives ourselves.
360  // Real dt = _subproblem.parent()->dt();
361  // Moose::out << "dt=" << dt << std::endl;
362 
363  // Vector object for the velocity
364  RealVectorValue vel(_u_vel[qp], _v_vel[qp], _w_vel[qp]);
365 
366  // A VectorValue object containing all zeros. Makes it easier to
367  // construct type tensor objects
368  RealVectorValue zero(0., 0., 0.);
369 
370  // Velocity vector magnitude squared
371  Real velmag2 = vel.norm_sq();
372 
373  // Debugging: How large are the time derivative parts of the strong residuals?
374  // Moose::out << "drho_dt=" << _drho_dt
375  // << ", drhou_dt=" << _drhou_dt
376  // << ", drhov_dt=" << _drhov_dt
377  // << ", drhow_dt=" << _drhow_dt
378  // << ", drhoE_dt=" << _drhoE_dt
379  // << std::endl;
380 
381  // Momentum divergence
382  Real divU = _grad_rho_u[qp](0) + _grad_rho_v[qp](1) + _grad_rho_w[qp](2);
383 
384  // Enough space to hold three space dimensions of velocity components at each qp,
385  // regardless of what dimension we are actually running in.
386  _calC[qp].resize(3);
387 
388  // Explicitly zero the calC
389  for (unsigned int i = 0; i < 3; ++i)
390  _calC[qp][i].zero();
391 
392  // x-column matrix
393  _calC[qp][0](0, 0) = _u_vel[qp];
394  _calC[qp][0](1, 0) = _v_vel[qp];
395  _calC[qp][0](2, 0) = _w_vel[qp];
396 
397  // y-column matrix
398  _calC[qp][1](0, 1) = _u_vel[qp];
399  _calC[qp][1](1, 1) = _v_vel[qp];
400  _calC[qp][1](2, 1) = _w_vel[qp];
401 
402  // z-column matrix (this assumes LIBMESH_DIM==3!)
403  _calC[qp][2](0, 2) = _u_vel[qp];
404  _calC[qp][2](1, 2) = _v_vel[qp];
405  _calC[qp][2](2, 2) = _w_vel[qp];
406 
407  // The matrix S can be computed from any of the calC via calC_1*calC_1^T
408  RealTensorValue calS = _calC[qp][0] * _calC[qp][0].transpose();
409 
410  // Enough space to hold five (=n_sd + 2) 3*3 calA matrices at this qp, regarless of dimension
411  _calA[qp].resize(5);
412 
413  // 0.) _calA_0 = diag( (gam - 1)/2*|u|^2 ) - S
414  _calA[qp][0].zero(); // zero this calA entry
415  _calA[qp][0](0, 0) = _calA[qp][0](1, 1) = _calA[qp][0](2, 2) =
416  0.5 * (_fp.gamma() - 1.0) * velmag2; // set diag. entries
417  _calA[qp][0] -= calS;
418 
419  for (unsigned int m = 1; m <= 3; ++m)
420  {
421  // Use m_local when indexing into matrices and vectors
422  unsigned int m_local = m - 1;
423 
424  // For m=1,2,3, calA_m = C_m + C_m^T + diag( (1.-gam)*u_m )
425  _calA[qp][m].zero(); // zero this calA entry
426  _calA[qp][m](0, 0) = _calA[qp][m](1, 1) = _calA[qp][m](2, 2) =
427  (1. - _fp.gamma()) * vel(m_local); // set diag. entries
428  _calA[qp][m] += _calC[qp][m_local]; // Note: use m_local for indexing into _calC!
429  _calA[qp][m] += _calC[qp][m_local].transpose(); // Note: use m_local for indexing into _calC!
430  }
431 
432  // 4.) calA_4 = diag(gam - 1)
433  _calA[qp][4].zero(); // zero this calA entry
434  _calA[qp][4](0, 0) = _calA[qp][4](1, 1) = _calA[qp][4](2, 2) = (_fp.gamma() - 1.0);
435 
436  // Enough space to hold the 3*5 "cal E" matrices which comprise the inviscid flux term
437  // of the energy equation. See notes for additional details
438  _calE[qp].resize(3); // Three rows, 5 entries in each row
439 
440  for (unsigned int k = 0; k < 3; ++k)
441  {
442  // Make enough room to store all 5 E matrices for this k
443  _calE[qp][k].resize(5);
444 
445  // Store and reuse the velocity column transpose matrix for the
446  // current value of k.
447  RealTensorValue Ck_T = _calC[qp][k].transpose();
448 
449  // E_{k0} (density gradient term)
450  _calE[qp][k][0].zero();
451  _calE[qp][k][0] = (0.5 * (_fp.gamma() - 1.0) * velmag2 - _enthalpy[qp]) * Ck_T;
452 
453  for (unsigned int m = 1; m <= 3; ++m)
454  {
455  // Use m_local when indexing into matrices and vectors
456  unsigned int m_local = m - 1;
457 
458  // E_{km} (momentum gradient terms)
459  _calE[qp][k][m].zero();
460  _calE[qp][k][m](k, m_local) = _enthalpy[qp]; // H * D_{km}
461  _calE[qp][k][m] += (1. - _fp.gamma()) * vel(m_local) * Ck_T; // (1-gam) * u_m * C_k^T
462  }
463 
464  // E_{k4} (energy gradient term)
465  _calE[qp][k][4].zero();
466  _calE[qp][k][4] = _fp.gamma() * Ck_T;
467  }
468 
469  // Compute the sum over ell of: A_ell grad(U_ell), store in DenseVector or Gradient object?
470  // The gradient object might be more useful, since we are multiplying by VariableGradient
471  // (which is a MooseArray of RealGradients) objects?
472  RealVectorValue mom_resid = _calA[qp][0] * _grad_rho[qp] + _calA[qp][1] * _grad_rho_u[qp] +
473  _calA[qp][2] * _grad_rho_v[qp] + _calA[qp][3] * _grad_rho_w[qp] +
474  _calA[qp][4] * _grad_rho_E[qp];
475 
476  // No matrices/vectors for the energy residual strong form... just write it out like
477  // the mass equation residual. See "Momentum SUPG terms prop. to energy residual"
478  // section of the notes.
479  Real energy_resid =
480  (0.5 * (_fp.gamma() - 1.0) * velmag2 - _enthalpy[qp]) * (vel * _grad_rho[qp]) +
481  _enthalpy[qp] * divU +
482  (1. - _fp.gamma()) * (vel(0) * (vel * _grad_rho_u[qp]) + vel(1) * (vel * _grad_rho_v[qp]) +
483  vel(2) * (vel * _grad_rho_w[qp])) +
484  _fp.gamma() * (vel * _grad_rho_E[qp]);
485 
486  // Now for the actual residual values...
487 
488  // The density strong-residual
489  _strong_residuals[qp][0] = _drho_dt[qp] + divU;
490 
491  // The x-momentum strong-residual, viscous terms neglected.
492  // TODO: If we want to add viscous contributions back in, should this kernel
493  // not inherit from NSViscousFluxBase so it can get tau values? This would
494  // also involve shape function second derivative values.
495  _strong_residuals[qp][1] = _drhou_dt[qp] + mom_resid(0);
496 
497  // The y-momentum strong residual, viscous terms neglected.
498  _strong_residuals[qp][2] = _drhov_dt[qp] + mom_resid(1);
499 
500  // The z-momentum strong residual, viscous terms neglected.
501  if (_mesh_dimension == 3)
502  _strong_residuals[qp][3] = _drhow_dt[qp] + mom_resid(2);
503  else
504  _strong_residuals[qp][3] = 0.;
505 
506  // The energy equation strong residual
507  _strong_residuals[qp][4] = _drhoE_dt[qp] + energy_resid;
508 }
Definition: NS.h:13
const std::string momentum_x
Definition: NS.h:16
MaterialProperty< Real > & _dynamic_viscosity
MaterialProperty< Real > & _hsupg
const unsigned int _mesh_dimension
void computeTau(unsigned int qp)
const VariableValue & _drhow_dt
MaterialProperty< std::vector< RealTensorValue > > & _calC
const std::string velocity_z
Definition: NS.h:23
void computeHSUPG(unsigned int qp)
const VariableGradient & _grad_rho_v
const std::string density
Definition: NS.h:15
const std::string velocity_x
Definition: NS.h:21
const VariableValue & _u_vel
const std::string temperature
Definition: NS.h:25
MaterialProperty< Real > & _tauc
const VariableGradient & _grad_rho_w
const std::string enthalpy
Definition: NS.h:26
MaterialProperty< Real > & _taue
void computeStrongResiduals(unsigned int qp)
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _calE
std::vector< const VariableGradient * > _vel_grads
const VariableValue & _enthalpy
MaterialProperty< std::vector< RealTensorValue > > & _calA
virtual Real cp(Real v=0., Real u=0.) const override
Specific heat.
const VariableValue & _v_vel
MaterialProperty< std::vector< Real > > & _strong_residuals
const VariableGradient & _grad_rho_E
const std::string velocity_y
Definition: NS.h:22
const VariableGradient & _grad_rho_u
NavierStokesMaterial(const InputParameters &parameters)
const VariableValue & _drhou_dt
const std::string momentum_y
Definition: NS.h:17
const VariableGradient & _grad_rho
const VariableValue & _rho_w
const VariableValue & _w_vel
const VariableGradient & _grad_u
MaterialProperty< RealTensorValue > & _viscous_stress_tensor
virtual Real gamma(Real v=0., Real u=0.) const override
Compute the ratio of specific heats.
const VariableValue & _drhoE_dt
MaterialProperty< Real > & _thermal_conductivity
const std::string total_energy
Definition: NS.h:19
const std::string momentum_z
Definition: NS.h:18
const VariableGradient & _grad_v
virtual void computeProperties()
Must be called after the child class computes dynamic_viscocity.
const VariableValue & _drhov_dt
const IdealGasFluidProperties & _fp
const VariableValue & _drho_dt
const VariableValue & _rho_u
const VariableValue & _rho_E
const VariableValue & _temperature
const VariableValue & _rho_v
const VariableGradient & _grad_w
const VariableValue & _rho
MaterialProperty< Real > & _taum