DirectCentralDifference

Implementation of Explicit/Forward Euler without invoking any of the nonlinear solver

Overview

DirectCentralDifference applies a time integrator for central difference in which the acceleration used for the solution update is calculated directly from the residual forces.

The formulation assumes a constant acceleration between midpoints. An average between the old and current time step is used to increment midpoint velocity to account for changing time steps, which is the same method used in Abaqus.

For example if, then,

When using Dirichlet BCs, one must use the (DirectDirichletBC,DirectFunctionDirichletBC) variations to enforce Dirichlet BC's properly.

Additionally, the time integrator must be used with MassMatrix, with a properly tagged mass matrix.

Example Input File Syntax

An example input file is shown below:



[GlobalParams]
    displacements = 'disp_x disp_y'
[]

[Kernels]
    [DynamicSolidMechanics]
        displacements = 'disp_x disp_y'
    []
    [massmatrix]
        type = MassMatrix
        density = 1
        matrix_tags = 'system'
        variable = disp_x
    []
    [massmatrix_y]
        type = MassMatrix
        density = 1
        matrix_tags = 'system'
        variable = disp_y
    []
[]

[Executioner]
    type = Transient

    [TimeIntegrator]
        type = CentralDifferenceDirect
        mass_matrix_tag = 'system'
    []
[]

Input Parameters

  • mass_matrix_tagmassThe tag for the mass matrix

    Default:mass

    C++ Type:TagName

    Unit:(no unit assumed)

    Controllable:No

    Description:The tag for the mass matrix

  • 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

    Unit:(no unit assumed)

    Controllable:No

    Description:If set to true, will only compute the mass matrix in the first time step, and keep using it throughout the simulation.

  • variablesA subset of the variables that this time integrator should be applied to

    C++ Type:std::vector<VariableName>

    Unit:(no unit assumed)

    Controllable:No

    Description:A subset of the variables that this time integrator should be applied to

Optional Parameters

  • control_tagsAdds user-defined labels for accessing object parameters via control logic.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    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

    Unit:(no unit assumed)

    Controllable:No

    Description:Set the enabled status of the MooseObject.

Advanced Parameters