INSADTauMaterialTempl< T > Class Template Reference

#include <INSADTauMaterial.h>

Inheritance diagram for INSADTauMaterialTempl< T >:

Public Member Functions
	INSADTauMaterialTempl (const InputParameters &parameters)

Static Public Member Functions
static InputParameters	validParams ()

Protected Member Functions
virtual void	computeProperties () override
virtual void	computeQpProperties () override
void	computeHMax ()
	Compute the maximum dimension of the current element.
void	computeViscousStrongResidual ()
	Compute the viscous strong residual.
void	viscousTermRZ ()
	Compute the strong form corresponding to RZ pieces of the viscous term.
bool	doVelocityDerivatives () const
	Whether to seed with respect to velocity derivatives.

Protected Attributes
const Real	_alpha
ADMaterialProperty< Real > &	_tau
ADMaterialProperty< RealVectorValue > &	_viscous_strong_residual
	Strong residual corresponding to the momentum viscous term.
ADMaterialProperty< RealVectorValue > &	_momentum_strong_residual
	The strong residual of the momentum equation.
ADReal	_hmax
const VectorMooseVariable *const	_velocity_var
	The velocity variable.
const FEBase *const &	_scalar_lagrange_fe
	A scalar Lagrange FE data member to compute the velocity second derivatives since they're currently not supported for vector FE types.
std::vector< ADRealTensorValue >	_d2u
	Containers to hold the matrix of second spatial derivatives of velocity.
std::vector< ADRealTensorValue >	_d2v
std::vector< ADRealTensorValue >	_d2w
const unsigned int	_vel_number
	The velocity variable number.
const unsigned int	_vel_sys_number
	The velocity system number.
ADReal	_speed
	The speed of the medium.

Detailed Description

template<typename T>
class INSADTauMaterialTempl< T >

Definition at line 32 of file INSADTauMaterial.h.

Constructor & Destructor Documentation

INSADTauMaterialTempl()

template<typename T >

INSADTauMaterialTempl< T >::INSADTauMaterialTempl

(

const InputParameters &

parameters

)

Definition at line 152 of file INSADTauMaterial.h.

  : T(parameters),
    _alpha(this->template getParam<Real>("alpha")),
    _tau(this->template declareADProperty<Real>("tau")),
    _viscous_strong_residual(
        this->template declareADProperty<RealVectorValue>("viscous_strong_residual")),
    _momentum_strong_residual(
        this->template declareADProperty<RealVectorValue>("momentum_strong_residual")),
    _velocity_var(getVectorVar("velocity", 0)),
    _scalar_lagrange_fe(_assembly.getFE(FEType(_velocity_var->feType().order, LAGRANGE)
                                            .set_p_refinement(_velocity_var->feType().p_refinement),
                                        _mesh.dimension())),
    _vel_number(_velocity_var->number()),
    _vel_sys_number(_velocity_var->sys().number())
{
  _scalar_lagrange_fe->get_d2phi();
}

Member Function Documentation

computeHMax()

template<typename T >

void INSADTauMaterialTempl< T >::computeHMax ( )

protected

Compute the maximum dimension of the current element.

Definition at line 180 of file INSADTauMaterial.h.

{
  if (_disp_x_num == libMesh::invalid_uint || !ADReal::do_derivatives)
  {
    _hmax = _current_elem->hmax();
    return;
  }
 
  _hmax = 0;
  std::array<unsigned int, 3> disps = {_disp_x_num, _disp_y_num, _disp_z_num};
  std::array<unsigned int, 3> disp_sys_nums = {_disp_x_sys_num, _disp_y_sys_num, _disp_z_sys_num};
 
  for (unsigned int n_outer = 0; n_outer < _current_elem->n_vertices(); n_outer++)
    for (unsigned int n_inner = n_outer + 1; n_inner < _current_elem->n_vertices(); n_inner++)
    {
      VectorValue<ADReal> diff = (_current_elem->point(n_outer) - _current_elem->point(n_inner));
      for (const auto i : index_range(disps))
      {
        const auto disp_num = disps[i];
        if (disp_num == libMesh::invalid_uint)
          continue;
        const auto sys_num = disp_sys_nums[i];
 
        // Here we insert derivatives of the difference in nodal positions with respect to the
        // displacement degrees of freedom. From above, diff = outer_node_position -
        // inner_node_position
        diff(i).derivatives().insert(
            _current_elem->node_ref(n_outer).dof_number(sys_num, disp_num, 0)) = 1.;
        diff(i).derivatives().insert(
            _current_elem->node_ref(n_inner).dof_number(sys_num, disp_num, 0)) = -1.;
      }
 
      _hmax = std::max(_hmax, diff.norm_sq());
    }
 
  using std::sqrt;
  _hmax = sqrt(_hmax);
}

computeProperties()

template<typename T >

void INSADTauMaterialTempl< T >::computeProperties ( )

overrideprotectedvirtual

Definition at line 221 of file INSADTauMaterial.h.

{
  computeHMax();
  computeViscousStrongResidual();
 
  T::computeProperties();
}

computeQpProperties()

template<typename T >

void INSADTauMaterialTempl< T >::computeQpProperties ( )

overrideprotectedvirtual

Reimplemented in INSADStabilized3Eqn.

Definition at line 345 of file INSADTauMaterial.h.

{
  T::computeQpProperties();
 
  using std::sqrt;
 
  const auto nu = _mu[_qp] / _rho[_qp];
  const auto transient_part = _has_transient ? 4. / (_dt * _dt) : 0.;
  _speed = NS::computeSpeed<ADReal>(_relative_velocity[_qp]);
  _tau[_qp] = _alpha / sqrt(transient_part + (2. * _speed / _hmax) * (2. * _speed / _hmax) +
                            9. * (4. * nu / (_hmax * _hmax)) * (4. * nu / (_hmax * _hmax)));
 
  _momentum_strong_residual[_qp] =
      _advective_strong_residual[_qp] + _viscous_strong_residual[_qp] + _grad_p[_qp];
 
  if (_has_transient)
    _momentum_strong_residual[_qp] += _td_strong_residual[_qp];
 
  if (_has_gravity)
    _momentum_strong_residual[_qp] += _gravity_strong_residual[_qp];
 
  if (_has_boussinesq)
    _momentum_strong_residual[_qp] += _boussinesq_strong_residual[_qp];
 
  if (_has_advected_mesh)
    _momentum_strong_residual[_qp] += _advected_mesh_strong_residual[_qp];
 
  if (_has_coupled_force)
    _momentum_strong_residual[_qp] += _coupled_force_strong_residual[_qp];
 
  // // Future addition
  // if (_object_tracker->hasMMS())
  //   _momentum_strong_residual[_qp] += _mms_function_strong_residual[_qp];
}

Referenced by INSADStabilized3Eqn::computeQpProperties().

computeViscousStrongResidual()

template<typename T >

void INSADTauMaterialTempl< T >::computeViscousStrongResidual ( )

protected

Compute the viscous strong residual.

Definition at line 231 of file INSADTauMaterial.h.

{
  auto resize_and_zero = [this](auto & d2vel)
  {
    d2vel.resize(_qrule->n_points());
    for (auto & d2qp : d2vel)
      d2qp = 0;
  };
  resize_and_zero(_d2u);
  resize_and_zero(_d2v);
  resize_and_zero(_d2w);
 
  auto get_d2 = [this](const auto i) -> std::vector<ADRealTensorValue> &
  {
    switch (i)
    {
      case 0:
        return _d2u;
      case 1:
        return _d2v;
      case 2:
        return _d2w;
      default:
        mooseError("invalid value of 'i'");
    }
  };
 
  // libMesh does not yet have the capability for computing second order spatial derivatives of
  // vector bases. Lacking that capability, we can compute the second order spatial derivatives "by
  // hand" using the scalar field version of the vector basis, e.g. LAGRANGE instead of
  // LAGRANGE_VEC. Adding this implementation allows us to be fully consistent with results from a
  // scalar velocity field component implementation of Navier-Stokes
  const auto & vel_dof_indices = _velocity_var->dofIndices();
  for (const auto i : index_range(vel_dof_indices))
  {
    // This may not work if the element and spatial dimensions are different
    mooseAssert(_current_elem->dim() == _mesh.dimension(),
                "Below logic only applicable if element and mesh dimension are the same");
    const auto dimensional_component = i % _mesh.dimension();
    auto & d2vel = get_d2(dimensional_component);
    const auto dof_index = vel_dof_indices[i];
    ADReal dof_value = (*_velocity_var->sys().currentSolution())(dof_index);
    if (doVelocityDerivatives())
      dof_value.derivatives().insert(dof_index) = 1;
    const auto scalar_i_component = i / _mesh.dimension();
    for (const auto qp : make_range(_qrule->n_points()))
      d2vel[qp] += dof_value * _scalar_lagrange_fe->get_d2phi()[scalar_i_component][qp];
  }
 
  // Now that we have the second order spatial derivatives of velocity, we can compute the strong
  // form of the viscous residual for use in our stabilization calculations
  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
  {
    _viscous_strong_residual[_qp](0) = -_mu[_qp] * _d2u[_qp].tr();
    _viscous_strong_residual[_qp](1) =
        _mesh.dimension() >= 2 ? -_mu[_qp] * _d2v[_qp].tr() : ADReal(0);
    _viscous_strong_residual[_qp](2) =
        _mesh.dimension() == 3 ? -_mu[_qp] * _d2w[_qp].tr() : ADReal(0);
    if (_viscous_form == NS::ViscousForm::Traction)
    {
      _viscous_strong_residual[_qp] -= _mu[_qp] * _d2u[_qp].row(0);
      if (_mesh.dimension() >= 2)
        _viscous_strong_residual[_qp] -= _mu[_qp] * _d2v[_qp].row(1);
      if (_mesh.dimension() == 3)
        _viscous_strong_residual[_qp] -= _mu[_qp] * _d2w[_qp].row(2);
    }
    if (_coord_sys == Moose::COORD_RZ)
      viscousTermRZ();
  }
}

doVelocityDerivatives()

template<typename T >

bool INSADTauMaterialTempl< T >::doVelocityDerivatives ( ) const

protected

Whether to seed with respect to velocity derivatives.

Definition at line 172 of file INSADTauMaterial.h.

{
  return ADReal::do_derivatives &&
         (_vel_sys_number == _fe_problem.currentNonlinearSystem().number());
}

validParams()

template<typename T >

InputParameters INSADTauMaterialTempl< T >::validParams ( )

static

Definition at line 142 of file INSADTauMaterial.h.

{
  InputParameters params = T::validParams();
  params.addClassDescription(
      "This is the material class used to compute the stabilization parameter tau.");
  params.addParam<Real>("alpha", 1., "Multiplicative factor on the stabilization parameter tau.");
  return params;
}

Referenced by INSADStabilized3Eqn::validParams().

viscousTermRZ()

template<typename T >

void INSADTauMaterialTempl< T >::viscousTermRZ ( )

protected

Compute the strong form corresponding to RZ pieces of the viscous term.

Definition at line 304 of file INSADTauMaterial.h.

{
  // To understand the code immediately below, visit
  // https://en.wikipedia.org/wiki/Del_in_cylindrical_and_spherical_coordinates.
  // The u_r / r^2 term comes from the vector Laplacian. The -du_i/dr * 1/r term comes from
  // the scalar Laplacian. The scalar Laplacian in axisymmetric cylindrical coordinates is
  // equivalent to the Cartesian Laplacian plus a 1/r * du_i/dr term. And of course we are
  // applying a minus sign here because the strong form is -\nabala^2 * \vec{u}
  //
  // Another note: libMesh implements grad(v) as dvi/dxj or more obviously:
  //
  // grad(v) = (vx_x vx_y)
  //           (vy_x vy_y)
  //
  // so, the gradient of the velocity with respect to the radial coordinate will correspond to a
  // *column* slice
 
  const auto & r = _ad_q_point[_qp](_rz_radial_coord);
 
  {
    // Do the "Laplace" form. This will be present in *both* Laplace and Traction forms
    ADRealVectorValue rz_term;
    for (const auto i : make_range((unsigned int)2))
      rz_term(i) = -_mu[_qp] * _grad_velocity[_qp](i, _rz_radial_coord) / r;
    rz_term(_rz_radial_coord) += _mu[_qp] * _velocity[_qp](_rz_radial_coord) / (r * r);
    _viscous_strong_residual[_qp] += rz_term;
  }
  if (_viscous_form == NS::ViscousForm::Traction)
  {
    ADRealVectorValue rz_term;
    for (const auto i : make_range((unsigned int)2))
      // This is the transpose of the above
      rz_term(i) = -_mu[_qp] * _grad_velocity[_qp](_rz_radial_coord, i) / r;
    // This is the same as above (since the transpose of the diagonal is the diagonal)
    rz_term(_rz_radial_coord) += _mu[_qp] * _velocity[_qp](_rz_radial_coord) / (r * r);
    _viscous_strong_residual[_qp] += rz_term;
  }
}

Member Data Documentation

_alpha

template<typename T >

const Real INSADTauMaterialTempl< T >::_alpha

protected

Definition at line 63 of file INSADTauMaterial.h.

_d2u

template<typename T >

std::vector<ADRealTensorValue> INSADTauMaterialTempl< T >::_d2u

protected

Containers to hold the matrix of second spatial derivatives of velocity.

Definition at line 83 of file INSADTauMaterial.h.

_d2v

template<typename T >

std::vector<ADRealTensorValue> INSADTauMaterialTempl< T >::_d2v

protected

Definition at line 84 of file INSADTauMaterial.h.

_d2w

template<typename T >

std::vector<ADRealTensorValue> INSADTauMaterialTempl< T >::_d2w

protected

Definition at line 85 of file INSADTauMaterial.h.

_hmax

template<typename T >

ADReal INSADTauMaterialTempl< T >::_hmax

protected

Definition at line 73 of file INSADTauMaterial.h.

_momentum_strong_residual

template<typename T >

ADMaterialProperty<RealVectorValue>& INSADTauMaterialTempl< T >::_momentum_strong_residual

protected

The strong residual of the momentum equation.

Definition at line 71 of file INSADTauMaterial.h.

_scalar_lagrange_fe

template<typename T >

const FEBase* const& INSADTauMaterialTempl< T >::_scalar_lagrange_fe

protected

A scalar Lagrange FE data member to compute the velocity second derivatives since they're currently not supported for vector FE types.

Definition at line 80 of file INSADTauMaterial.h.

Referenced by INSADTauMaterialTempl< T >::INSADTauMaterialTempl().

_speed

template<typename T >

ADReal INSADTauMaterialTempl< T >::_speed

protected

The speed of the medium.

This is the norm of the relative velocity, e.g. the velocity minus the mesh velocity, at the current _qp

Definition at line 95 of file INSADTauMaterial.h.

_tau

template<typename T >

ADMaterialProperty<Real>& INSADTauMaterialTempl< T >::_tau

protected

Definition at line 64 of file INSADTauMaterial.h.

_vel_number

template<typename T >

const unsigned int INSADTauMaterialTempl< T >::_vel_number

protected

The velocity variable number.

Definition at line 88 of file INSADTauMaterial.h.

_vel_sys_number

template<typename T >

const unsigned int INSADTauMaterialTempl< T >::_vel_sys_number

protected

The velocity system number.

Definition at line 91 of file INSADTauMaterial.h.

_velocity_var

template<typename T >

const VectorMooseVariable* const INSADTauMaterialTempl< T >::_velocity_var

protected

The velocity variable.

Definition at line 76 of file INSADTauMaterial.h.

_viscous_strong_residual

template<typename T >

ADMaterialProperty<RealVectorValue>& INSADTauMaterialTempl< T >::_viscous_strong_residual

protected

Strong residual corresponding to the momentum viscous term.

This is only used by stabilization kernels

Definition at line 68 of file INSADTauMaterial.h.

---

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

• navier_stokes/include/materials/INSADTauMaterial.h