www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | List of all members
NavierStokesMaterial Class Reference

This is the base class all materials should use if you are trying to use the Navier-Stokes Kernels. More...

#include <NavierStokesMaterial.h>

Inheritance diagram for NavierStokesMaterial:
[legend]

Public Member Functions

 NavierStokesMaterial (const InputParameters &parameters)
 

Protected Member Functions

virtual void computeProperties ()
 Must be called after the child class computes dynamic_viscocity. More...
 

Protected Attributes

const unsigned int _mesh_dimension
 
const VariableGradient & _grad_u
 
const VariableGradient & _grad_v
 
const VariableGradient & _grad_w
 
MaterialProperty< RealTensorValue > & _viscous_stress_tensor
 
MaterialProperty< Real > & _thermal_conductivity
 
MaterialProperty< Real > & _dynamic_viscosity
 
MaterialProperty< std::vector< RealTensorValue > > & _calA
 
MaterialProperty< std::vector< RealTensorValue > > & _calC
 
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _calE
 
std::vector< const VariableGradient * > _vel_grads
 
const VariableValue & _u_vel
 
const VariableValue & _v_vel
 
const VariableValue & _w_vel
 
const VariableValue & _temperature
 
const VariableValue & _enthalpy
 
const VariableValue & _rho
 
const VariableValue & _rho_u
 
const VariableValue & _rho_v
 
const VariableValue & _rho_w
 
const VariableValue & _rho_E
 
const VariableValue & _drho_dt
 
const VariableValue & _drhou_dt
 
const VariableValue & _drhov_dt
 
const VariableValue & _drhow_dt
 
const VariableValue & _drhoE_dt
 
const VariableGradient & _grad_rho
 
const VariableGradient & _grad_rho_u
 
const VariableGradient & _grad_rho_v
 
const VariableGradient & _grad_rho_w
 
const VariableGradient & _grad_rho_E
 
MaterialProperty< Real > & _hsupg
 
MaterialProperty< Real > & _tauc
 
MaterialProperty< Real > & _taum
 
MaterialProperty< Real > & _taue
 
MaterialProperty< std::vector< Real > > & _strong_residuals
 
const IdealGasFluidProperties_fp
 

Private Member Functions

void computeHSUPG (unsigned int qp)
 
void computeTau (unsigned int qp)
 
void computeStrongResiduals (unsigned int qp)
 

Detailed Description

This is the base class all materials should use if you are trying to use the Navier-Stokes Kernels.

Note that the derived class just needs to compute dynamic_viscocity then call this class's computeProperties() function.

Also make sure that the derived class's validParams function just adds to this class's validParams.

Finally, note that this Material isn't registered with the MaterialFactory. The reason is that by itself this material doesn't work! You must derive from this material and compute dynamic_viscocity!

Definition at line 37 of file NavierStokesMaterial.h.

Constructor & Destructor Documentation

◆ NavierStokesMaterial()

NavierStokesMaterial::NavierStokesMaterial ( const InputParameters &  parameters)

Definition at line 49 of file NavierStokesMaterial.C.

50  : Material(parameters),
51  _mesh_dimension(_mesh.dimension()),
52  _grad_u(coupledGradient(NS::velocity_x)),
53  _grad_v(_mesh_dimension >= 2 ? coupledGradient(NS::velocity_y) : _grad_zero),
54  _grad_w(_mesh_dimension == 3 ? coupledGradient(NS::velocity_z) : _grad_zero),
55 
56  _viscous_stress_tensor(declareProperty<RealTensorValue>("viscous_stress_tensor")),
57  _thermal_conductivity(declareProperty<Real>("thermal_conductivity")),
58 
59  // Declared here but _not_ calculated here
60  // (See e.g. derived class, bighorn/include/materials/FluidTC1.h)
61  _dynamic_viscosity(declareProperty<Real>("dynamic_viscosity")),
62 
63  // The momentum components of the inviscid flux Jacobians.
64  _calA(declareProperty<std::vector<RealTensorValue>>("calA")),
65 
66  // "Velocity column" matrices
67  _calC(declareProperty<std::vector<RealTensorValue>>("calC")),
68 
69  // Energy equation inviscid flux matrices, "cal E_{kl}" in the notes.
70  _calE(declareProperty<std::vector<std::vector<RealTensorValue>>>("calE")),
72 
73  // Coupled solution values needed for computing SUPG stabilization terms
74  _u_vel(coupledValue(NS::velocity_x)),
75  _v_vel(_mesh.dimension() >= 2 ? coupledValue(NS::velocity_y) : _zero),
76  _w_vel(_mesh.dimension() == 3 ? coupledValue(NS::velocity_z) : _zero),
77 
78  _temperature(coupledValue(NS::temperature)),
79  _enthalpy(coupledValue(NS::enthalpy)),
80 
81  // Coupled solution values
82  _rho(coupledValue(NS::density)),
83  _rho_u(coupledValue(NS::momentum_x)),
84  _rho_v(_mesh.dimension() >= 2 ? coupledValue(NS::momentum_y) : _zero),
85  _rho_w(_mesh.dimension() == 3 ? coupledValue(NS::momentum_z) : _zero),
86  _rho_E(coupledValue(NS::total_energy)),
87 
88  // Time derivative values
89  _drho_dt(coupledDot(NS::density)),
90  _drhou_dt(coupledDot(NS::momentum_x)),
91  _drhov_dt(_mesh.dimension() >= 2 ? coupledDot(NS::momentum_y) : _zero),
92  _drhow_dt(_mesh.dimension() == 3 ? coupledDot(NS::momentum_z) : _zero),
93  _drhoE_dt(coupledDot(NS::total_energy)),
94 
95  // Gradients
96  _grad_rho(coupledGradient(NS::density)),
97  _grad_rho_u(coupledGradient(NS::momentum_x)),
98  _grad_rho_v(_mesh.dimension() >= 2 ? coupledGradient(NS::momentum_y) : _grad_zero),
99  _grad_rho_w(_mesh.dimension() == 3 ? coupledGradient(NS::momentum_z) : _grad_zero),
100  _grad_rho_E(coupledGradient(NS::total_energy)),
101 
102  // Material properties for stabilization
103  _hsupg(declareProperty<Real>("hsupg")),
104  _tauc(declareProperty<Real>("tauc")),
105  _taum(declareProperty<Real>("taum")),
106  _taue(declareProperty<Real>("taue")),
107  _strong_residuals(declareProperty<std::vector<Real>>("strong_residuals")),
108  _fp(getUserObject<IdealGasFluidProperties>("fluid_properties"))
109 {
110 }
Definition: NS.h:15
const std::string momentum_x
Definition: NS.h:18
MaterialProperty< Real > & _dynamic_viscosity
MaterialProperty< Real > & _hsupg
const unsigned int _mesh_dimension
const VariableValue & _drhow_dt
MaterialProperty< std::vector< RealTensorValue > > & _calC
const std::string velocity_z
Definition: NS.h:25
const VariableGradient & _grad_rho_v
const std::string density
Definition: NS.h:17
const std::string velocity_x
Definition: NS.h:23
const VariableValue & _u_vel
const std::string temperature
Definition: NS.h:27
MaterialProperty< Real > & _tauc
const VariableGradient & _grad_rho_w
const std::string enthalpy
Definition: NS.h:28
MaterialProperty< Real > & _taue
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _calE
std::vector< const VariableGradient * > _vel_grads
const VariableValue & _enthalpy
MaterialProperty< std::vector< RealTensorValue > > & _calA
const VariableValue & _v_vel
MaterialProperty< std::vector< Real > > & _strong_residuals
const VariableGradient & _grad_rho_E
const std::string velocity_y
Definition: NS.h:24
const VariableGradient & _grad_rho_u
const VariableValue & _drhou_dt
const std::string momentum_y
Definition: NS.h:19
const VariableGradient & _grad_rho
const VariableValue & _rho_w
const VariableValue & _w_vel
const VariableGradient & _grad_u
MaterialProperty< RealTensorValue > & _viscous_stress_tensor
const VariableValue & _drhoE_dt
MaterialProperty< Real > & _thermal_conductivity
const std::string total_energy
Definition: NS.h:21
const std::string momentum_z
Definition: NS.h:20
const VariableGradient & _grad_v
Ideal gas fluid properties.
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

Member Function Documentation

◆ computeHSUPG()

void NavierStokesMaterial::computeHSUPG ( unsigned int  qp)
private

Definition at line 178 of file NavierStokesMaterial.C.

Referenced by computeProperties().

179 {
180  // // Grab reference to linear Lagrange finite element object pointer,
181  // // currently this is always a linear Lagrange element, so this might need to
182  // // be generalized if we start working with higher-order elements...
183  // FEBase*& fe(_assembly.getFE(FEType(), _current_elem->dim()));
184  //
185  // // Grab references to FE object's mapping data from the _subproblem's FE object
186  // const std::vector<Real> & dxidx(fe->get_dxidx());
187  // const std::vector<Real> & dxidy(fe->get_dxidy());
188  // const std::vector<Real> & dxidz(fe->get_dxidz());
189  // const std::vector<Real> & detadx(fe->get_detadx());
190  // const std::vector<Real> & detady(fe->get_detady());
191  // const std::vector<Real> & detadz(fe->get_detadz());
192  // const std::vector<Real> & dzetadx(fe->get_dzetadx()); // Empty in 2D
193  // const std::vector<Real> & dzetady(fe->get_dzetady()); // Empty in 2D
194  // const std::vector<Real> & dzetadz(fe->get_dzetadz()); // Empty in 2D
195  //
196  // // Bounds checking on element data
197  // mooseAssert(qp < dxidx.size(), "Insufficient data in dxidx array!");
198  // mooseAssert(qp < dxidy.size(), "Insufficient data in dxidy array!");
199  // mooseAssert(qp < dxidz.size(), "Insufficient data in dxidz array!");
200  //
201  // mooseAssert(qp < detadx.size(), "Insufficient data in detadx array!");
202  // mooseAssert(qp < detady.size(), "Insufficient data in detady array!");
203  // mooseAssert(qp < detadz.size(), "Insufficient data in detadz array!");
204  //
205  // if (_mesh_dimension == 3)
206  // {
207  // mooseAssert(qp < dzetadx.size(), "Insufficient data in dzetadx array!");
208  // mooseAssert(qp < dzetady.size(), "Insufficient data in dzetady array!");
209  // mooseAssert(qp < dzetadz.size(), "Insufficient data in dzetadz array!");
210  // }
211  //
212  // // The velocity vector at this quadrature point.
213  // RealVectorValue U(_u_vel[qp],_v_vel[qp],_w_vel[qp]);
214  //
215  // // Pull out element inverse map values at the current qp into a little dense matrix
216  // Real dxi_dx[3][3] = {{0.,0.,0.}, {0.,0.,0.}, {0.,0.,0.}};
217  //
218  // dxi_dx[0][0] = dxidx[qp]; dxi_dx[0][1] = dxidy[qp];
219  // dxi_dx[1][0] = detadx[qp]; dxi_dx[1][1] = detady[qp];
220  //
221  // // OK to access third entries on 2D elements if LIBMESH_DIM==3, though they
222  // // may be zero...
223  // if (LIBMESH_DIM == 3)
224  // {
225  // /**/ /**/ dxi_dx[0][2] = dxidz[qp];
226  // /**/ /**/ dxi_dx[1][2] = detadz[qp];
227  // }
228  //
229  // // The last row of entries available only for 3D elements.
230  // if (_mesh_dimension == 3)
231  // {
232  // dxi_dx[2][0] = dzetadx[qp]; dxi_dx[2][1] = dzetady[qp]; dxi_dx[2][2] = dzetadz[qp];
233  // }
234  //
235  // // Construct the g_ij = d(xi_k)/d(x_j) * d(xi_k)/d(x_i) matrix
236  // // from Ben and Bova's paper by summing over k...
237  // Real g[3][3] = {{0.,0.,0.}, {0.,0.,0.}, {0.,0.,0.}};
238  // for (unsigned int i = 0; i < 3; ++i)
239  // for (unsigned int j = 0; j < 3; ++j)
240  // for (unsigned int k = 0; k < 3; ++k)
241  // g[i][j] += dxi_dx[k][j] * dxi_dx[k][i];
242  //
243  // // Compute the denominator of the h_supg term: U * (g) * U
244  // Real denom = 0.;
245  // for (unsigned int i = 0; i < 3; ++i)
246  // for (unsigned int j = 0; j < 3; ++j)
247  // denom += U(j) * g[i][j] * U(i);
248  //
249  // // Compute h_supg. Some notes:
250  // // .) The 2 coefficient in this term should be a 1 if we are using tets/triangles.
251  // // .) The denominator will be identically zero only if the velocity
252  // // is identically zero, in which case we can't divide by it.
253  // if (denom != 0.0)
254  // _hsupg[qp] = 2.* sqrt( U.norm_sq() / denom );
255  // else
256  // _hsupg[qp] = 0.;
257 
258  // Simple (and fast) implementation: Just use hmin for the element!
259  _hsupg[qp] = _current_elem->hmin();
260 }
MaterialProperty< Real > & _hsupg

◆ computeProperties()

void NavierStokesMaterial::computeProperties ( )
protectedvirtual

Must be called after the child class computes dynamic_viscocity.

Must be called after the child class computes dynamic_viscosity.

Reimplemented in Air.

Definition at line 116 of file NavierStokesMaterial.C.

Referenced by Air::computeProperties().

117 {
118  for (unsigned int qp = 0; qp < _qrule->n_points(); ++qp)
119  {
120  /******* Viscous Stress Tensor *******/
121  // Technically... this _is_ the transpose (since we are loading these by rows)
122  // But it doesn't matter....
123  RealTensorValue grad_outer_u(_grad_u[qp], _grad_v[qp], _grad_w[qp]);
124 
125  grad_outer_u += grad_outer_u.transpose();
126 
127  Real div_vel = 0.0;
128  for (unsigned int i = 0; i < 3; ++i)
129  div_vel += (*_vel_grads[i])[qp](i);
130 
131  // Add diagonal terms
132  for (unsigned int i = 0; i < 3; ++i)
133  grad_outer_u(i, i) -= 2.0 / 3.0 * div_vel;
134 
135  grad_outer_u *= _dynamic_viscosity[qp];
136 
137  _viscous_stress_tensor[qp] = grad_outer_u;
138 
139  // Tabulated values of thermal conductivity vs. Temperature for air (k increases slightly with
140  // T):
141  // T (K) k (W/m-K)
142  // 273 0.0243
143  // 373 0.0314
144  // 473 0.0386
145  // 573 0.0454
146  // 673 0.0515
147 
148  // Pr = (mu * cp) / k ==> k = (mu * cp) / Pr = (mu * gamma * cv) / Pr.
149  // TODO: We are using a fixed value of the Prandtl number which is
150  // valid for air, it may or may not depend on temperature? Since
151  // this is a property of the fluid, it could possibly be moved to
152  // the FluidProperties module...
153  const Real Pr = 0.71;
154  _thermal_conductivity[qp] = (_dynamic_viscosity[qp] * _fp.cp()) / Pr;
155 
156  // Compute stabilization parameters:
157 
158  // .) Compute SUPG element length scale.
159  computeHSUPG(qp);
160  // Moose::out << "_hsupg[" << qp << "]=" << _hsupg[qp] << std::endl;
161 
162  // .) Compute SUPG parameter values. (Must call this after computeHSUPG())
163  computeTau(qp);
164  // Moose::out << "_tauc[" << qp << "]=" << _tauc[qp] << ", ";
165  // Moose::out << "_taum[" << qp << "]=" << _taum[qp] << ", ";
166  // Moose::out << "_taue[" << qp << "]=" << _taue[qp] << std::endl;
167 
168  // .) Compute strong residual values.
170  // Moose::out << "_strong_residuals[" << qp << "]=";
171  // for (unsigned i=0; i<_strong_residuals[qp].size(); ++i)
172  // Moose::out << _strong_residuals[qp][i] << " ";
173  // Moose::out << std::endl;
174  }
175 }
MaterialProperty< Real > & _dynamic_viscosity
void computeTau(unsigned int qp)
void computeHSUPG(unsigned int qp)
void computeStrongResiduals(unsigned int qp)
std::vector< const VariableGradient * > _vel_grads
const VariableGradient & _grad_u
MaterialProperty< RealTensorValue > & _viscous_stress_tensor
MaterialProperty< Real > & _thermal_conductivity
const VariableGradient & _grad_v
const IdealGasFluidProperties & _fp
const VariableGradient & _grad_w

◆ computeStrongResiduals()

void NavierStokesMaterial::computeStrongResiduals ( unsigned int  qp)
private

Definition at line 353 of file NavierStokesMaterial.C.

Referenced by computeProperties().

354 {
355  // Create storage at this qp for the strong residuals of all the equations.
356  // In 2D, the value for the z-velocity equation will just be zero.
357  _strong_residuals[qp].resize(5);
358 
359  // The timestep is stored in the Problem object, which can be accessed through
360  // the parent pointer of the SubProblem. Don't need this if we are not
361  // approximating time derivatives ourselves.
362  // Real dt = _subproblem.parent()->dt();
363  // Moose::out << "dt=" << dt << std::endl;
364 
365  // Vector object for the velocity
366  RealVectorValue vel(_u_vel[qp], _v_vel[qp], _w_vel[qp]);
367 
368  // A VectorValue object containing all zeros. Makes it easier to
369  // construct type tensor objects
370  RealVectorValue zero(0., 0., 0.);
371 
372  // Velocity vector magnitude squared
373  Real velmag2 = vel.norm_sq();
374 
375  // Debugging: How large are the time derivative parts of the strong residuals?
376  // Moose::out << "drho_dt=" << _drho_dt
377  // << ", drhou_dt=" << _drhou_dt
378  // << ", drhov_dt=" << _drhov_dt
379  // << ", drhow_dt=" << _drhow_dt
380  // << ", drhoE_dt=" << _drhoE_dt
381  // << std::endl;
382 
383  // Momentum divergence
384  Real divU = _grad_rho_u[qp](0) + _grad_rho_v[qp](1) + _grad_rho_w[qp](2);
385 
386  // Enough space to hold three space dimensions of velocity components at each qp,
387  // regardless of what dimension we are actually running in.
388  _calC[qp].resize(3);
389 
390  // Explicitly zero the calC
391  for (unsigned int i = 0; i < 3; ++i)
392  _calC[qp][i].zero();
393 
394  // x-column matrix
395  _calC[qp][0](0, 0) = _u_vel[qp];
396  _calC[qp][0](1, 0) = _v_vel[qp];
397  _calC[qp][0](2, 0) = _w_vel[qp];
398 
399  // y-column matrix
400  _calC[qp][1](0, 1) = _u_vel[qp];
401  _calC[qp][1](1, 1) = _v_vel[qp];
402  _calC[qp][1](2, 1) = _w_vel[qp];
403 
404  // z-column matrix (this assumes LIBMESH_DIM==3!)
405  _calC[qp][2](0, 2) = _u_vel[qp];
406  _calC[qp][2](1, 2) = _v_vel[qp];
407  _calC[qp][2](2, 2) = _w_vel[qp];
408 
409  // The matrix S can be computed from any of the calC via calC_1*calC_1^T
410  RealTensorValue calS = _calC[qp][0] * _calC[qp][0].transpose();
411 
412  // Enough space to hold five (=n_sd + 2) 3*3 calA matrices at this qp, regarless of dimension
413  _calA[qp].resize(5);
414 
415  // 0.) _calA_0 = diag( (gam - 1)/2*|u|^2 ) - S
416  _calA[qp][0].zero(); // zero this calA entry
417  _calA[qp][0](0, 0) = _calA[qp][0](1, 1) = _calA[qp][0](2, 2) =
418  0.5 * (_fp.gamma() - 1.0) * velmag2; // set diag. entries
419  _calA[qp][0] -= calS;
420 
421  for (unsigned int m = 1; m <= 3; ++m)
422  {
423  // Use m_local when indexing into matrices and vectors
424  unsigned int m_local = m - 1;
425 
426  // For m=1,2,3, calA_m = C_m + C_m^T + diag( (1.-gam)*u_m )
427  _calA[qp][m].zero(); // zero this calA entry
428  _calA[qp][m](0, 0) = _calA[qp][m](1, 1) = _calA[qp][m](2, 2) =
429  (1. - _fp.gamma()) * vel(m_local); // set diag. entries
430  _calA[qp][m] += _calC[qp][m_local]; // Note: use m_local for indexing into _calC!
431  _calA[qp][m] += _calC[qp][m_local].transpose(); // Note: use m_local for indexing into _calC!
432  }
433 
434  // 4.) calA_4 = diag(gam - 1)
435  _calA[qp][4].zero(); // zero this calA entry
436  _calA[qp][4](0, 0) = _calA[qp][4](1, 1) = _calA[qp][4](2, 2) = (_fp.gamma() - 1.0);
437 
438  // Enough space to hold the 3*5 "cal E" matrices which comprise the inviscid flux term
439  // of the energy equation. See notes for additional details
440  _calE[qp].resize(3); // Three rows, 5 entries in each row
441 
442  for (unsigned int k = 0; k < 3; ++k)
443  {
444  // Make enough room to store all 5 E matrices for this k
445  _calE[qp][k].resize(5);
446 
447  // Store and reuse the velocity column transpose matrix for the
448  // current value of k.
449  RealTensorValue Ck_T = _calC[qp][k].transpose();
450 
451  // E_{k0} (density gradient term)
452  _calE[qp][k][0].zero();
453  _calE[qp][k][0] = (0.5 * (_fp.gamma() - 1.0) * velmag2 - _enthalpy[qp]) * Ck_T;
454 
455  for (unsigned int m = 1; m <= 3; ++m)
456  {
457  // Use m_local when indexing into matrices and vectors
458  unsigned int m_local = m - 1;
459 
460  // E_{km} (momentum gradient terms)
461  _calE[qp][k][m].zero();
462  _calE[qp][k][m](k, m_local) = _enthalpy[qp]; // H * D_{km}
463  _calE[qp][k][m] += (1. - _fp.gamma()) * vel(m_local) * Ck_T; // (1-gam) * u_m * C_k^T
464  }
465 
466  // E_{k4} (energy gradient term)
467  _calE[qp][k][4].zero();
468  _calE[qp][k][4] = _fp.gamma() * Ck_T;
469  }
470 
471  // Compute the sum over ell of: A_ell grad(U_ell), store in DenseVector or Gradient object?
472  // The gradient object might be more useful, since we are multiplying by VariableGradient
473  // (which is a MooseArray of RealGradients) objects?
474  RealVectorValue mom_resid = _calA[qp][0] * _grad_rho[qp] + _calA[qp][1] * _grad_rho_u[qp] +
475  _calA[qp][2] * _grad_rho_v[qp] + _calA[qp][3] * _grad_rho_w[qp] +
476  _calA[qp][4] * _grad_rho_E[qp];
477 
478  // No matrices/vectors for the energy residual strong form... just write it out like
479  // the mass equation residual. See "Momentum SUPG terms prop. to energy residual"
480  // section of the notes.
481  Real energy_resid =
482  (0.5 * (_fp.gamma() - 1.0) * velmag2 - _enthalpy[qp]) * (vel * _grad_rho[qp]) +
483  _enthalpy[qp] * divU +
484  (1. - _fp.gamma()) * (vel(0) * (vel * _grad_rho_u[qp]) + vel(1) * (vel * _grad_rho_v[qp]) +
485  vel(2) * (vel * _grad_rho_w[qp])) +
486  _fp.gamma() * (vel * _grad_rho_E[qp]);
487 
488  // Now for the actual residual values...
489 
490  // The density strong-residual
491  _strong_residuals[qp][0] = _drho_dt[qp] + divU;
492 
493  // The x-momentum strong-residual, viscous terms neglected.
494  // TODO: If we want to add viscous contributions back in, should this kernel
495  // not inherit from NSViscousFluxBase so it can get tau values? This would
496  // also involve shape function second derivative values.
497  _strong_residuals[qp][1] = _drhou_dt[qp] + mom_resid(0);
498 
499  // The y-momentum strong residual, viscous terms neglected.
500  _strong_residuals[qp][2] = _drhov_dt[qp] + mom_resid(1);
501 
502  // The z-momentum strong residual, viscous terms neglected.
503  if (_mesh_dimension == 3)
504  _strong_residuals[qp][3] = _drhow_dt[qp] + mom_resid(2);
505  else
506  _strong_residuals[qp][3] = 0.;
507 
508  // The energy equation strong residual
509  _strong_residuals[qp][4] = _drhoE_dt[qp] + energy_resid;
510 }
const unsigned int _mesh_dimension
const VariableValue & _drhow_dt
MaterialProperty< std::vector< RealTensorValue > > & _calC
const VariableGradient & _grad_rho_v
const VariableValue & _u_vel
const VariableGradient & _grad_rho_w
MaterialProperty< std::vector< std::vector< RealTensorValue > > > & _calE
const VariableValue & _enthalpy
MaterialProperty< std::vector< RealTensorValue > > & _calA
const VariableValue & _v_vel
MaterialProperty< std::vector< Real > > & _strong_residuals
const VariableGradient & _grad_rho_E
const VariableGradient & _grad_rho_u
const VariableValue & _drhou_dt
const VariableGradient & _grad_rho
const VariableValue & _w_vel
const VariableValue & _drhoE_dt
const VariableValue & _drhov_dt
const IdealGasFluidProperties & _fp
const VariableValue & _drho_dt

◆ computeTau()

void NavierStokesMaterial::computeTau ( unsigned int  qp)
private

Definition at line 263 of file NavierStokesMaterial.C.

Referenced by computeProperties().

264 {
265  Real velmag =
266  std::sqrt(_u_vel[qp] * _u_vel[qp] + _v_vel[qp] * _v_vel[qp] + _w_vel[qp] * _w_vel[qp]);
267 
268  // Moose::out << "velmag=" << velmag << std::endl;
269 
270  // Make sure temperature >= 0 before trying to take sqrt
271  // if (_temperature[qp] < 0.)
272  // {
273  // Moose::err << "Negative temperature "
274  // << _temperature[qp]
275  // << " found at quadrature point "
276  // << qp
277  // << ", element "
278  // << _current_elem->id()
279  // << std::endl;
280  // mooseError("Can't continue, would be nice to throw an exception here?");
281  // }
282 
283  // The speed of sound for an ideal gas, sqrt(gamma * R * T). Not needed unless
284  // we want to use a form of Tau that requires it.
285  // Real soundspeed = _fp.c_from_v_e(_specific_volume[_qp], _internal_energy[_qp]);
286 
287  // If velmag == 0, then _hsupg should be zero as well. Then tau
288  // will have only the time-derivative contribution (or zero, if we
289  // are not including dt terms in our taus!) Note that using the
290  // time derivative contribution in this way assumes we are solving
291  // unsteady, and guarantees *some* stabilization is added even when
292  // u -> 0 in certain regions of the flow.
293  if (velmag == 0.)
294  {
295  // 1.) Tau without dt terms
296  // _tauc[qp] = 0.;
297  // _taum[qp] = 0.;
298  // _taue[qp] = 0.;
299 
300  // 2.) Tau *with* dt terms
301  _tauc[qp] = _taum[qp] = _taue[qp] = 0.5 * _dt;
302  }
303  else
304  {
305  // The element length parameter, squared
306  Real h2 = _hsupg[qp] * _hsupg[qp];
307 
308  // The viscosity-based term
309  Real visc_term = _dynamic_viscosity[qp] / _rho[qp] / h2;
310 
311  // The thermal conductivity-based term, cp = gamma * cv
312  Real k_term = _thermal_conductivity[qp] / _rho[qp] / _fp.cp() / h2;
313 
314  // 1a.) Standard compressible flow tau. Does not account for low Mach number
315  // limit.
316  // _tauc[qp] = _hsupg[qp] / (velmag + soundspeed);
317 
318  // 1b.) Inspired by Hauke, the sum of the compressible and incompressible tauc.
319  // _tauc[qp] =
320  // _hsupg[qp] / (velmag + soundspeed) +
321  // _hsupg[qp] / (velmag);
322 
323  // 1c.) From Wong 2001. This tau is O(M^2) for small M. At small M,
324  // tauc dominates the inverse square sums and basically makes
325  // taum=taue=tauc. However, all my flows occur at low Mach numbers,
326  // so there would basically never be any stabilization...
327  // _tauc[qp] = (_hsupg[qp] * velmag) / (velmag*velmag + soundspeed*soundspeed);
328 
329  // For use with option "1",
330  // (tau_c)^{-2}
331  // Real taucm2 = 1./_tauc[qp]/_tauc[qp];
332  // _taum[qp] = 1. / std::sqrt(taucm2 + visc_term*visc_term);
333  // _taue[qp] = 1. / std::sqrt(taucm2 + k_term*k_term);
334 
335  // 2.) Tau with timestep dependence (guarantees stabilization even
336  // in zero-velocity limit) incorporated via the "r-switch" method,
337  // with r=2.
338  Real sqrt_term = 4. / _dt / _dt + velmag * velmag / h2;
339 
340  // For use with option "2", i.e. the option that uses dt in the definition of tau
341  _tauc[qp] = 1. / std::sqrt(sqrt_term);
342  _taum[qp] = 1. / std::sqrt(sqrt_term + visc_term * visc_term);
343  _taue[qp] = 1. / std::sqrt(sqrt_term + k_term * k_term);
344  }
345 
346  // Debugging
347  // Moose::out << "_tauc[" << qp << "]=" << _tauc[qp] << std::endl;
348  // Moose::out << "_hsupg[" << qp << "]=" << _hsupg[qp] << std::endl;
349  // Moose::out << "velmag[" << qp << "]=" << velmag << std::endl;
350 }
MaterialProperty< Real > & _dynamic_viscosity
MaterialProperty< Real > & _hsupg
const VariableValue & _u_vel
MaterialProperty< Real > & _tauc
MaterialProperty< Real > & _taue
const VariableValue & _v_vel
const VariableValue & _w_vel
MaterialProperty< Real > & _thermal_conductivity
const IdealGasFluidProperties & _fp
const VariableValue & _rho
MaterialProperty< Real > & _taum

Member Data Documentation

◆ _calA

MaterialProperty<std::vector<RealTensorValue> >& NavierStokesMaterial::_calA
protected

Definition at line 62 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _calC

MaterialProperty<std::vector<RealTensorValue> >& NavierStokesMaterial::_calC
protected

Definition at line 66 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _calE

MaterialProperty<std::vector<std::vector<RealTensorValue> > >& NavierStokesMaterial::_calE
protected

Definition at line 71 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _drho_dt

const VariableValue& NavierStokesMaterial::_drho_dt
protected

Definition at line 97 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _drhoE_dt

const VariableValue& NavierStokesMaterial::_drhoE_dt
protected

Definition at line 101 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _drhou_dt

const VariableValue& NavierStokesMaterial::_drhou_dt
protected

Definition at line 98 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _drhov_dt

const VariableValue& NavierStokesMaterial::_drhov_dt
protected

Definition at line 99 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _drhow_dt

const VariableValue& NavierStokesMaterial::_drhow_dt
protected

Definition at line 100 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _dynamic_viscosity

MaterialProperty<Real>& NavierStokesMaterial::_dynamic_viscosity
protected

Definition at line 56 of file NavierStokesMaterial.h.

Referenced by Air::computeProperties(), computeProperties(), and computeTau().

◆ _enthalpy

const VariableValue& NavierStokesMaterial::_enthalpy
protected

Definition at line 87 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _fp

const IdealGasFluidProperties& NavierStokesMaterial::_fp
protected

Definition at line 122 of file NavierStokesMaterial.h.

Referenced by computeProperties(), computeStrongResiduals(), and computeTau().

◆ _grad_rho

const VariableGradient& NavierStokesMaterial::_grad_rho
protected

Definition at line 104 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _grad_rho_E

const VariableGradient& NavierStokesMaterial::_grad_rho_E
protected

Definition at line 108 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _grad_rho_u

const VariableGradient& NavierStokesMaterial::_grad_rho_u
protected

Definition at line 105 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _grad_rho_v

const VariableGradient& NavierStokesMaterial::_grad_rho_v
protected

Definition at line 106 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _grad_rho_w

const VariableGradient& NavierStokesMaterial::_grad_rho_w
protected

Definition at line 107 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _grad_u

const VariableGradient& NavierStokesMaterial::_grad_u
protected

Definition at line 50 of file NavierStokesMaterial.h.

Referenced by computeProperties(), and NavierStokesMaterial().

◆ _grad_v

const VariableGradient& NavierStokesMaterial::_grad_v
protected

Definition at line 51 of file NavierStokesMaterial.h.

Referenced by computeProperties(), and NavierStokesMaterial().

◆ _grad_w

const VariableGradient& NavierStokesMaterial::_grad_w
protected

Definition at line 52 of file NavierStokesMaterial.h.

Referenced by computeProperties(), and NavierStokesMaterial().

◆ _hsupg

MaterialProperty<Real>& NavierStokesMaterial::_hsupg
protected

Definition at line 112 of file NavierStokesMaterial.h.

Referenced by computeHSUPG(), and computeTau().

◆ _mesh_dimension

const unsigned int NavierStokesMaterial::_mesh_dimension
protected

Definition at line 48 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _rho

const VariableValue& NavierStokesMaterial::_rho
protected

Definition at line 90 of file NavierStokesMaterial.h.

Referenced by computeTau().

◆ _rho_E

const VariableValue& NavierStokesMaterial::_rho_E
protected

Definition at line 94 of file NavierStokesMaterial.h.

◆ _rho_u

const VariableValue& NavierStokesMaterial::_rho_u
protected

Definition at line 91 of file NavierStokesMaterial.h.

◆ _rho_v

const VariableValue& NavierStokesMaterial::_rho_v
protected

Definition at line 92 of file NavierStokesMaterial.h.

◆ _rho_w

const VariableValue& NavierStokesMaterial::_rho_w
protected

Definition at line 93 of file NavierStokesMaterial.h.

◆ _strong_residuals

MaterialProperty<std::vector<Real> >& NavierStokesMaterial::_strong_residuals
protected

Definition at line 119 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals().

◆ _tauc

MaterialProperty<Real>& NavierStokesMaterial::_tauc
protected

Definition at line 113 of file NavierStokesMaterial.h.

Referenced by computeTau().

◆ _taue

MaterialProperty<Real>& NavierStokesMaterial::_taue
protected

Definition at line 115 of file NavierStokesMaterial.h.

Referenced by computeTau().

◆ _taum

MaterialProperty<Real>& NavierStokesMaterial::_taum
protected

Definition at line 114 of file NavierStokesMaterial.h.

Referenced by computeTau().

◆ _temperature

const VariableValue& NavierStokesMaterial::_temperature
protected

Definition at line 84 of file NavierStokesMaterial.h.

◆ _thermal_conductivity

MaterialProperty<Real>& NavierStokesMaterial::_thermal_conductivity
protected

Definition at line 55 of file NavierStokesMaterial.h.

Referenced by computeProperties(), and computeTau().

◆ _u_vel

const VariableValue& NavierStokesMaterial::_u_vel
protected

Definition at line 79 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals(), and computeTau().

◆ _v_vel

const VariableValue& NavierStokesMaterial::_v_vel
protected

Definition at line 80 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals(), and computeTau().

◆ _vel_grads

std::vector<const VariableGradient *> NavierStokesMaterial::_vel_grads
protected

Definition at line 75 of file NavierStokesMaterial.h.

Referenced by computeProperties().

◆ _viscous_stress_tensor

MaterialProperty<RealTensorValue>& NavierStokesMaterial::_viscous_stress_tensor
protected

Definition at line 54 of file NavierStokesMaterial.h.

Referenced by computeProperties().

◆ _w_vel

const VariableValue& NavierStokesMaterial::_w_vel
protected

Definition at line 81 of file NavierStokesMaterial.h.

Referenced by computeStrongResiduals(), and computeTau().


The documentation for this class was generated from the following files: