ExplicitMidpoint
Time integration using the explicit midpoint method.
The explicit midpoint method is second-order accurate in time. It is a two-step method and a special case of the 2nd-order Runge-Kutta method.
Description
With , the vector of nonlinear variables, and , a nonlinear operator, we write the PDE of interest as:
Using for the current time step, and for the previous step, the explicit midpoint integration scheme can be written:
This method can be expressed as a Runge-Kutta method with the following Butcher Tableau:
All kernels except time-(derivative)-kernels should have the parameter implicit=false
to use this time integrator.
ExplicitRK2-derived TimeIntegrators (ExplicitMidpoint, Heun, Ralston) and other multistage TimeIntegrators are known not to work with Materials/AuxKernels that accumulate 'state' and should be used with caution.
Input Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:No
Description:Set the enabled status of the MooseObject.
Input Files
- (test/tests/time_integrators/rk-2/2d-quadratic.i)
- (test/tests/time_integrators/scalar/scalar.i)
- (modules/rdg/test/tests/advection_1d/1d_aefv_square_wave.i)
- (modules/rdg/test/tests/advection_1d/block_restrictable.i)
- (test/tests/time_integrators/rk-2/1d-linear.i)
- (test/tests/dgkernels/1d_advection_dg/1d_advection_dg.i)
(test/tests/time_integrators/rk-2/2d-quadratic.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 1
ymin = -1
ymax = 1
nx = 20
ny = 20
elem_type = QUAD9
[]
[Functions]
[./ic]
type = ParsedFunction
expression = 0
[../]
[./forcing_fn]
type = ParsedFunction
expression = 2*t*((x*x)+(y*y))-(4*t*t)
[../]
[./exact_fn]
type = ParsedFunction
expression = t*t*((x*x)+(y*y))
[../]
[]
[Variables]
[./u]
order = SECOND
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
implicit = true
[../]
[./diff]
type = Diffusion
variable = u
implicit = false
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
implicit = false
[../]
[]
[BCs]
active = 'all'
[./all]
type = FunctionDirichletBC
variable = u
boundary = '0 1 2 3'
function = exact_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
[./TimeIntegrator]
type = ExplicitMidpoint
[../]
solve_type = 'LINEAR'
start_time = 0.0
num_steps = 10
dt = 0.0001
l_tol = 1e-8
[]
[Outputs]
exodus = true
perf_graph = true
[]
(test/tests/time_integrators/scalar/scalar.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 1
ny = 1
elem_type = QUAD4
[]
[Variables]
[./n]
family = SCALAR
order = FIRST
initial_condition = 1
[../]
[]
[ScalarKernels]
[./dn]
type = ODETimeDerivative
variable = n
[../]
[./ode1]
type = ParsedODEKernel
expression = '-n'
variable = n
# implicit = false
[../]
[]
[Executioner]
type = Transient
[./TimeIntegrator]
# type = ImplicitEuler
# type = BDF2
type = CrankNicolson
# type = ImplicitMidpoint
# type = LStableDirk2
# type = LStableDirk3
# type = LStableDirk4
# type = AStableDirk4
#
# Explicit methods
# type = ExplicitEuler
# type = ExplicitMidpoint
# type = Heun
# type = Ralston
[../]
start_time = 0
end_time = 1
dt = 0.001
dtmin = 0.001 # Don't allow timestep cutting
solve_type = 'PJFNK'
nl_max_its = 2
nl_abs_tol = 1.e-12 # This is an ODE, so nl_abs_tol makes sense.
[]
[Functions]
[./exact_solution]
type = ParsedFunction
expression = exp(t)
[../]
[]
[Postprocessors]
[./error_n]
# Post processor that computes the difference between the computed
# and exact solutions. For the exact solution used here, the
# error at the final time should converge at O(dt^p), where p is
# the order of the method.
type = ScalarL2Error
variable = n
function = exact_solution
# final is not currently supported for Postprocessor execute_on...
# execute_on = 'final'
[../]
[]
[Outputs]
csv = true
[]
(modules/rdg/test/tests/advection_1d/1d_aefv_square_wave.i)
############################################################
[GlobalParams]
order = CONSTANT
family = MONOMIAL
u = u
slope_limiting = lslope
implicit = false
[]
############################################################
[Mesh]
type = GeneratedMesh
dim = 1
xmin = 0
xmax = 1
nx = 100
[]
############################################################
[Functions]
[./ic_u]
type = PiecewiseConstant
axis = x
direction = right
xy_data = '0.1 0.5
0.6 1.0
1.0 0.5'
[../]
[]
############################################################
[UserObjects]
[./lslope]
type = AEFVSlopeLimitingOneD
execute_on = 'linear'
scheme = 'none' #none | minmod | mc | superbee
[../]
[./internal_side_flux]
type = AEFVUpwindInternalSideFlux
execute_on = 'linear'
[../]
[./free_outflow_bc]
type = AEFVFreeOutflowBoundaryFlux
execute_on = 'linear'
[../]
[]
############################################################
[Variables]
[./u]
[../]
[]
############################################################
[ICs]
[./u_ic]
type = FunctionIC
variable = 'u'
function = ic_u
[../]
[]
############################################################
[Kernels]
[./time_u]
implicit = true
type = TimeDerivative
variable = u
[../]
[]
############################################################
[DGKernels]
[./concentration]
type = AEFVKernel
variable = u
component = 'concentration'
flux = internal_side_flux
[../]
[]
############################################################
[BCs]
[./concentration]
type = AEFVBC
boundary = 'left right'
variable = u
component = 'concentration'
flux = free_outflow_bc
[../]
[]
############################################################
[Materials]
[./aefv]
type = AEFVMaterial
block = 0
[../]
[]
############################################################
[Executioner]
type = Transient
[./TimeIntegrator]
type = ExplicitMidpoint
[../]
solve_type = 'LINEAR'
l_tol = 1e-4
nl_rel_tol = 1e-20
nl_abs_tol = 1e-8
nl_max_its = 60
start_time = 0.0
num_steps = 4 # 4 | 400 for complete run
dt = 5e-4
dtmin = 1e-6
[]
[Outputs]
[./Exodus]
type = Exodus
file_base = 1d_aefv_square_wave_none_out
time_step_interval = 2
[../]
perf_graph = true
[]
(modules/rdg/test/tests/advection_1d/block_restrictable.i)
############################################################
[GlobalParams]
order = CONSTANT
family = MONOMIAL
u = u
slope_limiting = lslope
implicit = false
[]
############################################################
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 1
xmin = 0
xmax = 1
nx = 100
[]
[./subdomain1]
type = SubdomainBoundingBoxGenerator
bottom_left = '0.5 0 0'
block_id = 1
top_right = '1.0 1.0 0'
input = gen
[../]
[./interface]
type = SideSetsBetweenSubdomainsGenerator
primary_block = '0'
paired_block = '1'
new_boundary = 'primary0_interface'
input = subdomain1
[../]
[./interface_again]
type = SideSetsBetweenSubdomainsGenerator
primary_block = '1'
paired_block = '0'
new_boundary = 'primary1_interface'
input = interface
[../]
[]
############################################################
[Functions]
[./ic_u]
type = PiecewiseConstant
axis = x
direction = right
xy_data = '0.1 0.5
0.4 1.0
0.5 0.5'
[../]
[]
############################################################
[UserObjects]
[./lslope]
type = AEFVSlopeLimitingOneD
execute_on = 'linear'
scheme = 'superbee' #none | minmod | mc | superbee
block = 0
[../]
[./internal_side_flux]
type = AEFVUpwindInternalSideFlux
execute_on = 'linear'
[../]
[./free_outflow_bc]
type = AEFVFreeOutflowBoundaryFlux
execute_on = 'linear'
[../]
[]
############################################################
[Variables]
[./u]
block = 0
[../]
[./v]
block = 1
family = LAGRANGE
order = FIRST
[../]
[]
############################################################
[ICs]
[./u_ic]
type = FunctionIC
variable = 'u'
function = ic_u
[../]
[]
############################################################
[Kernels]
[./time_u]
implicit = true
type = TimeDerivative
variable = u
block = 0
[../]
[./diff_v]
implicit = true
type = Diffusion
variable = v
block = 1
[../]
[./time_v]
implicit = true
type = TimeDerivative
variable = v
block = 1
[../]
[]
############################################################
[DGKernels]
[./concentration]
type = AEFVKernel
variable = u
component = 'concentration'
flux = internal_side_flux
block = 0
[../]
[]
############################################################
[BCs]
[./concentration]
type = AEFVBC
boundary = 'left primary0_interface'
variable = u
component = 'concentration'
flux = free_outflow_bc
[../]
[./v_left]
type = DirichletBC
boundary = 'primary1_interface'
variable = v
value = 1
[../]
[./v_right]
type = DirichletBC
boundary = 'right'
variable = v
value = 0
[../]
[]
############################################################
[Materials]
[./aefv]
type = AEFVMaterial
block = 0
[../]
[./dummy_1]
type = GenericConstantMaterial
block = 1
prop_names = ''
prop_values = ''
[../]
[]
############################################################
[Executioner]
type = Transient
[./TimeIntegrator]
type = ExplicitMidpoint
[../]
solve_type = 'LINEAR'
l_tol = 1e-4
nl_rel_tol = 1e-20
nl_abs_tol = 1e-8
nl_max_its = 60
start_time = 0.0
num_steps = 4 # 4 | 400 for complete run
dt = 5e-4
dtmin = 1e-6
[]
[Outputs]
[./out]
type = Exodus
time_step_interval = 2
[../]
perf_graph = true
[]
(test/tests/time_integrators/rk-2/1d-linear.i)
[Mesh]
type = GeneratedMesh
dim = 1
xmin = -1
xmax = 1
nx = 20
elem_type = EDGE2
[]
[Functions]
[./ic]
type = ParsedFunction
expression = 0
[../]
[./forcing_fn]
type = ParsedFunction
expression = x
[../]
[./exact_fn]
type = ParsedFunction
expression = t*x
[../]
[]
[Variables]
[./u]
order = FIRST
family = LAGRANGE
[./InitialCondition]
type = FunctionIC
function = ic
[../]
[../]
[]
[Kernels]
[./ie]
type = TimeDerivative
variable = u
implicit = true
[../]
[./diff]
type = Diffusion
variable = u
implicit = false
[../]
[./ffn]
type = BodyForce
variable = u
function = forcing_fn
implicit = false
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
boundary = '0 1'
function = exact_fn
[../]
[]
[Postprocessors]
[./l2_err]
type = ElementL2Error
variable = u
function = exact_fn
[../]
[]
[Executioner]
type = Transient
[./TimeIntegrator]
type = ExplicitMidpoint
[../]
solve_type = 'LINEAR'
start_time = 0.0
num_steps = 10
dt = 0.001
l_tol = 1e-15
[]
[Outputs]
exodus = true
perf_graph = true
[]
(test/tests/dgkernels/1d_advection_dg/1d_advection_dg.i)
[Mesh]
type = GeneratedMesh
dim = 1
nx = 100
xmin = 0
xmax = 1
[]
[Functions]
[ic_u]
type = PiecewiseConstant
axis = x
direction = right
xy_data = '0.1 0.0
0.6 1.0
1.0 0.0'
[]
[]
[Variables]
[u]
order = FIRST
family = MONOMIAL
[]
[]
[Kernels]
[time_u]
type = TimeDerivative
variable = u
[]
[adv_u]
implicit = false
type = ConservativeAdvection
variable = u
velocity = '1 0 0'
[]
[]
[DGKernels]
[dg_advection_u]
implicit = false
type = DGConvection
variable = u
velocity = '1 0 0'
[]
[]
[ICs]
[u_ic]
type = FunctionIC
variable = u
function = ic_u
[]
[]
[Executioner]
type = Transient
[TimeIntegrator]
type = ExplicitMidpoint
[]
solve_type = 'LINEAR'
num_steps = 4
dt = 2e-4
[]
[Outputs]
exodus = true
[]