www.mooseframework.org
INSExplicitTimestepSelector.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 
11 
12 #include "libmesh/quadrature.h"
13 
15 
18 {
20 
21  params.addClassDescription("Postprocessor that computes the minimum value of h_min/|u|, where "
22  "|u| is coupled in as an aux variable.");
23 
24  // Coupled variables
25  params.addRequiredCoupledVar("vel_mag", "Velocity magnitude");
26 
27  // Required parameters
28  params.addRequiredParam<Real>("beta",
29  "0 < beta < 1, choose some fraction of the limiting timestep size");
30 
31  // Optional parameters
32  params.addParam<MaterialPropertyName>("mu_name", "mu", "The name of the dynamic viscosity");
33  params.addParam<MaterialPropertyName>("rho_name", "rho", "The name of the density");
34 
35  return params;
36 }
37 
39  : ElementPostprocessor(parameters),
40  _vel_mag(coupledValue("vel_mag")),
41 
42  // Other parameters
43  _beta(getParam<Real>("beta")),
44 
45  // Material properties
46  _mu(getMaterialProperty<Real>("mu_name")),
47  _rho(getMaterialProperty<Real>("rho_name"))
48 {
49 }
50 
52 
53 void
55 {
56  _value = std::numeric_limits<Real>::max();
57 }
58 
59 void
61 {
62  Real h_min = _current_elem->hmin();
63 
64  // The space dimension plays a role in the diffusive dt limit. The more
65  // space dimensions there are, the smaller this limit is.
66  Real dim = static_cast<Real>(_current_elem->dim());
67 
68  for (unsigned qp = 0; qp < _qrule->n_points(); ++qp)
69  {
70  // Don't divide by zero...
71  Real vel_mag = std::max(_vel_mag[qp], std::numeric_limits<Real>::epsilon());
72 
73  // For explicit Euler, we always have to satisfy the Courant condition for stability.
74  Real courant_limit_dt = h_min / vel_mag;
75 
76  // But we also have to obey the diffusive time restriction,
77  // dt <= 1/(2*nu)*(1/h1^2 + 1/h2^2 + 1/h3^2)^(-1) <=
78  // <= h_min^2 / n_dim / (2*nu)
79  Real diffusive_limit_dt = 0.5 * h_min * h_min / (_mu[qp] / _rho[qp]) / dim;
80 
81  // And the "combined" condition, dt <= 2*nu/|u|^2
82  Real combined_limit_dt = 2. * (_mu[qp] / _rho[qp]) / vel_mag / vel_mag;
83 
84  // // Debugging:
85  // Moose::out << "courant_limit_dt = " << courant_limit_dt << "\n"
86  // << "diffusive_limit_dt = " << diffusive_limit_dt << "\n"
87  // << "combined_limit_dt = " << combined_limit_dt
88  // << std::endl;
89 
90  _value = std::min(
91  _value,
92  _beta * std::min(std::min(courant_limit_dt, diffusive_limit_dt), combined_limit_dt));
93  }
94 }
95 
96 Real
98 {
99  return _value;
100 }
101 
102 void
104 {
105  gatherMin(_value);
106 }
107 
108 void
110 {
111  const auto & pps = static_cast<const INSExplicitTimestepSelector &>(uo);
112  _value = std::min(_value, pps._value);
113 }
Real _beta
We can compute maximum stable timesteps based on the linearized theory, but even those timesteps are ...
const MaterialProperty< Real > & _mu
Material properties: the explicit time scheme limit for the viscous problem also depends on the kinem...
static InputParameters validParams()
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
registerMooseObject("NavierStokesApp", INSExplicitTimestepSelector)
unsigned int dim
static InputParameters validParams()
void gatherMin(T &value)
void addRequiredParam(const std::string &name, const std::string &doc_string)
const MaterialProperty< Real > & _rho
Postprocessor that computes the minimum value of h_min/|u|, where |u| is coupled in as an aux variabl...
const VariableValue & _vel_mag
Velocity magnitude. Hint: Use VectorMagnitudeAux in Moose for this.
virtual Real getValue() const override
void addRequiredCoupledVar(const std::string &name, const std::string &doc_string)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const QBase *const & _qrule
const Elem *const & _current_elem
INSExplicitTimestepSelector(const InputParameters &parameters)
Real _value
The value of dt (NOTE: _dt member variable is already defined)
void addClassDescription(const std::string &doc_string)
virtual void threadJoin(const UserObject &uo) override