- solve_typeconsistentThe way to solve the system. A 'consistent' solve uses the full mass matrix and actually needs to use a linear solver to solve the problem. 'lumped' uses a lumped mass matrix with a simple inversion - incredibly fast but may be less accurate. 'lump_preconditioned' uses the lumped mass matrix as a preconditioner for the 'consistent' solve
Default:consistent
C++ Type:MooseEnum
Description:The way to solve the system. A 'consistent' solve uses the full mass matrix and actually needs to use a linear solver to solve the problem. 'lumped' uses a lumped mass matrix with a simple inversion - incredibly fast but may be less accurate. 'lump_preconditioned' uses the lumped mass matrix as a preconditioner for the 'consistent' solve
- use_constant_massFalseIf set to true, will only compute the mass matrix in the first time step, and keep using it throughout the simulation.
Default:False
C++ Type:bool
Description:If set to true, will only compute the mass matrix in the first time step, and keep using it throughout the simulation.
ActuallyExplicitEuler
This time integrator is an implementation of forward/explicit Euler based on ExplicitTimeIntegrator, which avoids a nonlinear solve.
Input Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Options:
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Options:
Description:Set the enabled status of the MooseObject.
Advanced Parameters
Input Files
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-1d-quadratic.i)
- (test/tests/time_integrators/actually_explicit_euler/diverged.i)
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-1d-linear.i)
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-2d-linear-adapt.i)
- (test/tests/time_integrators/actually_explicit_euler/actually_explicit_euler_lumped.i)
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-ode.i)
- (test/tests/time_integrators/actually_explicit_euler/actually_explicit_euler.i)
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-1d-quadratic-neumann.i)
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-2d-linear.i)
- (test/tests/time_integrators/actually_explicit_euler/actually_explicit_euler_lump_preconditioned.i)
- (test/tests/time_integrators/actually_explicit_euler_verification/ee-2d-quadratic.i)
- (modules/navier_stokes/test/tests/finite_volume/cns/scalar_advection/mass-frac-advection.i)
Child Objects
(test/tests/time_integrators/actually_explicit_euler_verification/ee-1d-quadratic.i)
[Mesh]
type = GeneratedMesh
dim = 1
xmin = -1
xmax = 1
nx = 20
elem_type = EDGE3
[]
[Functions]
[./ic]
type = ParsedFunction
value = 0
[../]
[./forcing_fn]
type = ParsedFunction
value = x*x-2*t
[../]
[./exact_fn]
type = ParsedFunction
value = t*x*x
[../]
[]
[Variables]
[./u]
order = SECOND
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
[../]
[./diff]
type = Diffusion
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
preset = false
boundary = '0 1'
function = exact_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
l_tol = 1e-12
start_time = 0.0
num_steps = 20
dt = 0.00005
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
max_rows = 10
[../]
[]
(test/tests/time_integrators/actually_explicit_euler/diverged.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 10
ny = 10
[]
[Variables]
[./u]
[../]
[]
[Kernels]
[./diff]
type = CoefDiffusion
variable = u
coef = 0.1
[../]
[./nan]
type = NanKernel
variable = u
timestep_to_nan = 4
[../]
[./time]
type = TimeDerivative
variable = u
[../]
[]
[BCs]
[./left]
type = DirichletBC
variable = u
boundary = 'left'
value = 0
[../]
[./right]
type = DirichletBC
variable = u
boundary = 'right'
value = 1
[../]
[]
[Executioner]
type = Transient
num_steps = 10
dt = 0.001
l_tol = 1e-12
dtmin = 1e-8
[./TimeIntegrator]
type = ActuallyExplicitEuler
solve_type = lump_preconditioned
[../]
[]
[Outputs]
exodus = false
[]
(test/tests/time_integrators/actually_explicit_euler_verification/ee-1d-linear.i)
[Mesh]
type = GeneratedMesh
dim = 1
xmin = -1
xmax = 1
nx = 200
elem_type = EDGE2
[]
[Functions]
[./ic]
type = ParsedFunction
value = 0
[../]
[./forcing_fn]
type = ParsedFunction
value = x
[../]
[./exact_fn]
type = ParsedFunction
value = t*x
[../]
[]
[Variables]
[./u]
order = FIRST
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
lumping = true
[../]
[./diff]
type = Diffusion
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
preset = false
boundary = '0 1'
function = exact_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
start_time = 0.0
num_steps = 20
dt = 0.00005
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
max_rows = 10
[../]
[]
(test/tests/time_integrators/actually_explicit_euler_verification/ee-2d-linear-adapt.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 1
ymin = -1
ymax = 1
nx = 10
ny = 10
elem_type = QUAD4
[]
[Functions]
[./ic]
type = ParsedFunction
value = 0
[../]
[./forcing_fn]
type = ParsedFunction
value = (x+y)
[../]
[./exact_fn]
type = ParsedFunction
value = t*(x+y)
[../]
[]
[Variables]
[./u]
order = FIRST
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
[../]
[./diff]
type = Diffusion
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
preset = false
boundary = '0 1 2 3'
function = exact_fn
[../]
[]
[Adaptivity]
steps = 1
marker = box
max_h_level = 2
[./Markers]
[./box]
bottom_left = '-0.4 -0.4 0'
inside = refine
top_right = '0.4 0.4 0'
outside = do_nothing
type = BoxMarker
[../]
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
start_time = 0.0
num_steps = 4
dt = 0.005
l_tol = 1e-12
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
max_rows = 10
[../]
[]
(test/tests/time_integrators/actually_explicit_euler/actually_explicit_euler_lumped.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 10
ny = 10
[]
[Variables]
[./u]
[../]
[]
[Kernels]
[./diff]
type = CoefDiffusion
variable = u
coef = 0.1
[../]
[./time]
type = TimeDerivative
variable = u
[../]
[]
[BCs]
[./left]
type = DirichletBC
variable = u
preset = false
boundary = 'left'
value = 0
[../]
[./right]
type = DirichletBC
variable = u
preset = false
boundary = 'right'
value = 1
[../]
[]
[Executioner]
type = Transient
num_steps = 10
dt = 0.001
[./TimeIntegrator]
type = ActuallyExplicitEuler
solve_type = lumped
[../]
[]
[Outputs]
exodus = true
[]
(test/tests/time_integrators/actually_explicit_euler_verification/ee-ode.i)
# Tests that ActuallyExplicitEuler works with scalar variables.
#
# The ODE and IC used are the following:
# du/dt = 2, u(0) = 0
# Thus the solution is u(t) = 2*t.
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
[]
[Variables]
[./u]
family = SCALAR
order = FIRST
initial_condition = 0
[../]
[]
[ScalarKernels]
[./time]
type = ODETimeDerivative
variable = u
[../]
[./source]
type = ParsedODEKernel
variable = u
function = -2
[../]
[]
[Executioner]
type = Transient
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
dt = 1
num_steps = 5
[]
[Outputs]
csv = true
[]
(test/tests/time_integrators/actually_explicit_euler/actually_explicit_euler.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 10
ny = 10
[]
[Variables]
[./u]
[../]
[]
[Kernels]
[./diff]
type = CoefDiffusion
variable = u
coef = 0.1
[../]
[./time]
type = TimeDerivative
variable = u
[../]
[]
[BCs]
[./left]
type = DirichletBC
variable = u
preset = false
boundary = 'left'
value = 0
[../]
[./right]
type = DirichletBC
variable = u
preset = false
boundary = 'right'
value = 1
[../]
[]
[Executioner]
type = Transient
num_steps = 10
dt = 0.001
l_tol = 1e-12
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[]
(test/tests/time_integrators/actually_explicit_euler_verification/ee-1d-quadratic-neumann.i)
[Mesh]
type = GeneratedMesh
dim = 1
xmin = -1
xmax = 1
nx = 10
elem_type = EDGE3
[]
[Functions]
[./ic]
type = ParsedFunction
value = 0
[../]
[./forcing_fn]
type = ParsedFunction
value = x*x-2*t+t*x*x
[../]
[./exact_fn]
type = ParsedFunction
value = t*x*x
[../]
[./left_bc_fn]
type = ParsedFunction
value = -t*2*x
[../]
[./right_bc_fn]
type = ParsedFunction
value = t*2*x
[../]
[]
[Variables]
[./u]
order = SECOND
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./td]
type = TimeDerivative
variable = u
[../]
[./diff]
type = Diffusion
variable = u
[../]
[./abs]
type = Reaction
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
[../]
[]
[BCs]
[./left]
type = FunctionNeumannBC
variable = u
boundary = '0'
function = left_bc_fn
[../]
[./right]
type = FunctionNeumannBC
variable = u
boundary = '1'
function = right_bc_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
l_tol = 1e-12
start_time = 0.0
num_steps = 10
dt = 0.001
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
max_rows = 10
[../]
[]
(test/tests/time_integrators/actually_explicit_euler_verification/ee-2d-linear.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 1
ymin = -1
ymax = 1
nx = 10
ny = 10
elem_type = QUAD4
[]
[Functions]
[./ic]
type = ParsedFunction
value = 0
[../]
[./forcing_fn]
type = ParsedFunction
value = (x+y)
[../]
[./exact_fn]
type = ParsedFunction
value = t*(x+y)
[../]
[]
[Variables]
[./u]
order = FIRST
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
[../]
[./diff]
type = Diffusion
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
preset = false
boundary = '0 1 2 3'
function = exact_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
start_time = 0.0
num_steps = 20
dt = 0.00005
l_tol = 1e-12
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
max_rows = 10
[../]
[]
(test/tests/time_integrators/actually_explicit_euler/actually_explicit_euler_lump_preconditioned.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 10
ny = 10
[]
[Variables]
[./u]
[../]
[]
[Kernels]
[./diff]
type = CoefDiffusion
variable = u
coef = 0.1
[../]
[./time]
type = TimeDerivative
variable = u
[../]
[]
[BCs]
[./left]
type = DirichletBC
variable = u
preset = false
boundary = 'left'
value = 0
[../]
[./right]
type = DirichletBC
variable = u
preset = false
boundary = 'right'
value = 1
[../]
[]
[Executioner]
type = Transient
num_steps = 10
dt = 0.0001
l_tol = 1e-12
[./TimeIntegrator]
type = ActuallyExplicitEuler
solve_type = lump_preconditioned
[../]
[]
[Outputs]
exodus = true
[]
(test/tests/time_integrators/actually_explicit_euler_verification/ee-2d-quadratic.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 1
ymin = -1
ymax = 1
nx = 10
ny = 10
elem_type = QUAD9
[]
[Functions]
[./ic]
type = ParsedFunction
value = 0
[../]
[./forcing_fn]
type = ParsedFunction
value = ((x*x)+(y*y))-(4*t)
[../]
[./exact_fn]
type = ParsedFunction
value = t*((x*x)+(y*y))
[../]
[]
[Variables]
[./u]
order = SECOND
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
[../]
[./diff]
type = Diffusion
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
preset = false
boundary = '0 1 2 3'
function = exact_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
l_tol = 1e-13
start_time = 0.0
num_steps = 20
dt = 0.00005
[./TimeIntegrator]
type = ActuallyExplicitEuler
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
max_rows = 10
[../]
[]
(modules/navier_stokes/test/tests/finite_volume/cns/scalar_advection/mass-frac-advection.i)
rho_initial=1.29
p_initial=1.01e5
T=273.15
gamma=1.4
e_initial=${fparse p_initial / (gamma - 1) / rho_initial}
et_initial=${e_initial}
rho_et_initial=${fparse rho_initial * et_initial}
v_in=1
[GlobalParams]
fp = fp
[]
[Mesh]
[cartesian]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 1
nx = 2
ymin = 0
ymax = 10
ny = 20
[]
[]
[Modules]
[FluidProperties]
[fp]
type = IdealGasFluidProperties
[]
[]
[]
[Variables]
[rho]
type = MooseVariableFVReal
initial_condition = ${rho_initial}
[]
[rho_u]
type = MooseVariableFVReal
initial_condition = 1e-15
[]
[rho_v]
type = MooseVariableFVReal
initial_condition = 1e-15
[]
[rho_et]
type = MooseVariableFVReal
initial_condition = ${rho_et_initial}
scaling = 1e-5
[]
[mass_frac]
type = MooseVariableFVReal
initial_condition = 1e-15
[]
[]
[AuxVariables]
[U_x]
type = MooseVariableFVReal
[]
[U_y]
type = MooseVariableFVReal
[]
[pressure]
type = MooseVariableFVReal
[]
[temperature]
type = MooseVariableFVReal
[]
[courant]
type = MooseVariableFVReal
[]
[]
[AuxKernels]
[U_x]
type = ADMaterialRealAux
variable = U_x
property = vel_x
execute_on = 'timestep_end'
[]
[U_y]
type = ADMaterialRealAux
variable = U_y
property = vel_y
execute_on = 'timestep_end'
[]
[pressure]
type = ADMaterialRealAux
variable = pressure
property = pressure
execute_on = 'timestep_end'
[]
[temperature]
type = ADMaterialRealAux
variable = temperature
property = T_fluid
execute_on = 'timestep_end'
[]
[courant]
type = Courant
variable = courant
u = U_x
v = U_y
[]
[]
[FVKernels]
[mass_time]
type = FVPorosityTimeDerivative
variable = rho
[]
[mass_advection]
type = PCNSFVKT
variable = rho
eqn = "mass"
[]
[momentum_time_x]
type = FVTimeKernel
variable = rho_u
[]
[momentum_advection_and_pressure_x]
type = PCNSFVKT
variable = rho_u
eqn = "momentum"
momentum_component = 'x'
[]
[momentum_time_y]
type = FVTimeKernel
variable = rho_v
[]
[momentum_advection_and_pressure_y]
type = PCNSFVKT
variable = rho_v
eqn = "momentum"
momentum_component = 'y'
[]
[energy_time]
type = FVPorosityTimeDerivative
variable = rho_et
[]
[energy_advection]
type = PCNSFVKT
variable = rho_et
eqn = "energy"
[]
[mass_frac_time]
type = PCNSFVDensityTimeDerivative
variable = mass_frac
rho = rho
[]
[mass_frac_advection]
type = PCNSFVKT
variable = mass_frac
eqn = "scalar"
[]
[]
[Functions]
[ud_in]
type = ParsedVectorFunction
value_x = '0'
value_y = '${v_in}'
[]
[]
[FVBCs]
[rho_bottom]
type = PCNSFVStrongBC
boundary = 'bottom'
variable = rho
superficial_velocity = 'ud_in'
T_fluid = ${T}
eqn = 'mass'
[]
[rho_u_bottom]
type = PCNSFVStrongBC
boundary = 'bottom'
variable = rho_u
superficial_velocity = 'ud_in'
T_fluid = ${T}
eqn = 'momentum'
momentum_component = 'x'
[]
[rho_v_bottom]
type = PCNSFVStrongBC
boundary = 'bottom'
variable = rho_v
superficial_velocity = 'ud_in'
T_fluid = ${T}
eqn = 'momentum'
momentum_component = 'y'
[]
[rho_et_bottom]
type = PCNSFVStrongBC
boundary = 'bottom'
variable = rho_et
superficial_velocity = 'ud_in'
T_fluid = ${T}
eqn = 'energy'
[]
[mass_frac_bottom]
type = PCNSFVStrongBC
boundary = 'bottom'
variable = mass_frac
superficial_velocity = 'ud_in'
T_fluid = ${T}
scalar = 1
eqn = 'scalar'
[]
[rho_top]
type = PCNSFVStrongBC
boundary = 'top'
variable = rho
pressure = ${p_initial}
eqn = 'mass'
[]
[rho_u_top]
type = PCNSFVStrongBC
boundary = 'top'
variable = rho_u
pressure = ${p_initial}
eqn = 'momentum'
momentum_component = 'x'
[]
[rho_v_top]
type = PCNSFVStrongBC
boundary = 'top'
variable = rho_v
pressure = ${p_initial}
eqn = 'momentum'
momentum_component = 'y'
[]
[rho_et_top]
type = PCNSFVStrongBC
boundary = 'top'
variable = rho_et
pressure = ${p_initial}
eqn = 'energy'
[]
[mass_frac_top]
type = PCNSFVStrongBC
boundary = 'top'
variable = mass_frac
pressure = ${p_initial}
eqn = 'scalar'
[]
[momentum_x_walls]
type = PCNSFVImplicitMomentumPressureBC
variable = rho_u
boundary = 'left right'
momentum_component = 'x'
[]
[momentum_y_walls]
type = PCNSFVImplicitMomentumPressureBC
variable = rho_v
boundary = 'left right'
momentum_component = 'y'
[]
[]
[Materials]
[var_mat]
type = PorousConservedVarMaterial
rho = rho
rho_et = rho_et
superficial_rhou = rho_u
superficial_rhov = rho_v
fp = fp
porosity = porosity
[]
[porosity]
type = GenericConstantMaterial
prop_names = 'porosity'
prop_values = '1'
[]
[]
[Executioner]
type = Transient
[TimeIntegrator]
type = ActuallyExplicitEuler
[]
steady_state_detection = true
steady_state_tolerance = 1e-12
abort_on_solve_fail = true
dt = 5e-4
num_steps = 25
[]
[Outputs]
[out]
type = Exodus
execute_on = 'initial timestep_end'
[]
[dof]
type = DOFMap
execute_on = 'initial'
[]
[]
[Debug]
show_var_residual_norms = true
[]
(framework/include/timeintegrators/CentralDifference.h)
// This file is part of the MOOSE framework
// https://www.mooseframework.org
//
// All rights reserved, see COPYRIGHT for full restrictions
// https://github.com/idaholab/moose/blob/master/COPYRIGHT
//
// Licensed under LGPL 2.1, please see LICENSE for details
// https://www.gnu.org/licenses/lgpl-2.1.html
#pragma once
#include "ActuallyExplicitEuler.h"
// Forward declarations
class CentralDifference;
template <>
InputParameters validParams<CentralDifference>();
/**
* Implements a truly explicit (no nonlinear solve) Central Difference time
* integration scheme.
*/
class CentralDifference : public ActuallyExplicitEuler
{
public:
CentralDifference(const InputParameters & parameters);
virtual void initialSetup() override;
virtual int order() override { return 2; }
virtual void computeTimeDerivatives() override;
void computeADTimeDerivatives(DualReal & ad_u_dot,
const dof_id_type & dof,
DualReal & ad_u_dotdot) const override;
protected:
/// solution vector for \f$ {du^dotdot}\over{du} \f$
Real & _du_dotdot_du;
/// The older solution
const NumericVector<Number> & _solution_older;
const NumericVector<Number> & _solution_old_old_old;
/**
* Helper function that actually does the math for computing the time derivative
*/
template <typename T, typename T2, typename T3, typename T4, typename T5>
void computeTimeDerivativeHelper(T & u_dot,
T2 & u_dotdot,
const T3 & u_old,
const T4 & u_old_old,
const T5 & u_old_old_old) const;
};
template <typename T, typename T2, typename T3, typename T4, typename T5>
void
CentralDifference::computeTimeDerivativeHelper(T & u_dot,
T2 & u_dotdot,
const T3 & u_old,
const T4 & u_old_old,
const T5 & u_old_old_old) const
{
// computing first derivative
// using the Central Difference method
// u_dot_old = (first_term - second_term) / 2 / dt
// first_term = u
// second_term = u_older
u_dot = u_old;
u_dot -= u_old_old_old; // 'older than older' solution
u_dot *= 1.0 / (2.0 * _dt);
// computing second derivative
// using the Central Difference method
// u_dotdot_old = (first_term - second_term + third_term) / dt / dt
// first_term = u
// second_term = 2 * u_old
// third_term = u_older
u_dotdot = u_old;
u_dotdot -= u_old_old;
u_dotdot -= u_old_old;
u_dotdot += u_old_old_old;
u_dotdot *= 1.0 / (_dt * _dt);
}