14.5. Automatic Time Stepping

The method of automatic time stepping (or automatic loading) is one in which the time step size and/or the applied loads are automatically determined in response to the current state of the analysis under consideration. This method (accessed with AUTOTS,ON) may be applied to structural, thermal, electric, and magnetic analyses that are performed in the time domain (using the TIME command), and includes static (or steady state) (ANTYPE,STATIC) and dynamic (or transient) (ANTYPE,TRANS) situations.

An important point to be made here is that automatic loading always works through the adjustment of the time step size; and that the loads that are applied are automatically adjusted if ramped boundary conditions are activated (using KBC,0). In other words the time step size is always subjected to possible adjustment when automatic loading is engaged. Applied loads and boundary conditions, however, will vary according to how they are applied and whether the boundary conditions are stepped or ramped. That is why this method may also be thought of as automatic loading.

There are two important features of the automatic time stepping algorithm. The first feature concerns the ability to estimate the next time step size, based on current and past analysis conditions, and make proper load adjustments. In other words, given conditions at the current time, tn, and the previous time increment, Δtn, the primary aim is to determine the next time increment, Δtn+1. Since the determination of Δtn+1 is largely predictive, this part of the automatic time stepping algorithm is referred to as the time step prediction component.

The second feature of automatic time stepping is referred to as the time step bisection component. Its purpose is to decide whether or not to reduce the present time step size, Δtn, and redo the substep with a smaller step size. For example, working from the last converged solution at time point tn-1, the present solution begins with a predicted time step, Δtn. Equilibrium iterations are performed; and if proper convergence is either not achieved or not anticipated, this time step is reduced to Δtn/2 (i.e., it is bisected), and the analysis begins again from time tn-1. Multiple bisections can occur per substep for various reasons (discussed later).

14.5.1. Time Step Prediction

At a given converged solution at time, tn, and with the previous time increment, Δtn, the goal is to predict the appropriate time step size to use as the next substep. This step size is derived from the results of several unrelated computations and is most easily expressed as the minimization statement:

(14–57)

where:

Δteq = time increment which is limited by the number of equilibrium iterations needed for convergence at the last converged time point. The more iterations required for convergence, the smaller the predicted time step. This is a general measure of all active nonlinearities. Increasing the maximum number of equilibrium iterations (using the NEQIT command) will tend to promote larger time step sizes.
Δt1 = time increment which is limited by the response eigenvalue computation for 1st order systems (e.g., thermal transients) (input on the TINTP command).
Δt2 = time increment which is limited by the response frequency computation for 2nd order systems (e.g., structural dynamics). The aim is to maintain 20 points per cycle (described below). Note when the middle step criterion is used, this criterion can be turned off.
Δtg = time increment that represents the time point at which a gap or a nonlinear (multi-status) element will change abruptly from one condition to another (status change). KEYOPT(7) allows further control for the CONTAC elements.
Δtc = time increment based on the allowable creep strain increment (described below).
Δtp = time increment based on the allowable plastic strain increment. The limit is set at 5% per time step (described below).
Δtm = time increment which is limited by the middle step residual tolerance (described below) for 2nd order systems (e.g., structural dynamics) (input on the MIDTOL command). When it is enabled, the Δt2 criterion can be turned off.

Several trial step sizes are calculated, and the minimum one is selected for the next time step. This predicted value is further restricted to a range of values expressed by

(14–58)

and

(14–59)

where:

F = increase/decrease factor. F = 2, if static analysis; F = 3, if dynamic (see the ANTYPE and TIMINT commands)
Δtmax = maximum time step size (DTMAX from the DELTIM command or the equivalent quantity calculated from the NSUBST command)
Δtmin = minimum time step size (DTMIN from the DELTIM command or the equivalent quantity calculated from the NSUBST command)

In other words, the current time step is increased or decreased by at most a factor of 2 (or 3 if dynamic), and it may not be less than Δtmin or greater than Δtmax.

14.5.2. Time Step Bisection

When bisection occurs, the current substep solution (Δtn) is removed, and the time step size is reduced by 50%. If applied loads are ramped (KBC,0), then the current load increment is also reduced by the same amount. One or more bisections can take place for several reasons, namely:

  1. The number of equilibrium iterations used for this substep exceeds the number allowed (NEQIT command).

  2. It appears likely that all equilibrium iterations will be used.

  3. A negative pivot message was encountered in the solution, suggesting instability.

  4. The largest calculated displacement DOF exceeds the limit (DLIM on the NCNV command).

  5. An illegal element distortion is detected (e.g., negative radius in an axisymmetric analysis).

  6. For transient structural dynamics, when the middle step residual is greater than the given tolerance. This check is done only when the middle step residual check is enabled by the MIDTOL command.

More than one bisection may be performed per substep. However, bisection of the time-step size is limited by the minimum size (defined by DTMIN input on the DELTIM command or the equivalent NSUBST input).

14.5.3. The Response Eigenvalue for 1st Order Transients

The response eigenvalue is used in the computation of Δt1 and is defined as:

(14–60)

where:

λr = response eigenvalue (item RESEIG for POST26 SOLU command and *GET command)
{Δu} = substep solution vector (tn-1 to tn)
[KT] = the Dirichlet matrix. In a heat transfer or an electrical conduction analysis this matrix is referred to as the conductivity matrix; in magnetics this is called the magnetic “stiffness”. The superscript T denotes the use of a tangent matrix in nonlinear situations
[C] = the damping matrix. In heat transfer this is called the specific heat matrix.

The product of the response eigenvalue and the previous time step (Δtn) has been employed by Hughes([145]) for the evaluation of 1st order explicit/implicit systems. In Hughes([145]) the quantity Δtnλ is referred to as the “oscillation limit”, where λ is the maximum eigenvalue. For unconditionally stable systems, the primary restriction on time-step size is that the inequality Δtnλ >> 1 should be avoided. Hence it is very conservative to propose that Δtnλ = 1.

Since the time integration used employs the trapezoidal rule (Equation 15–40), all analyses of 1st order systems are unconditionally stable. The response eigenvalue supplied by means of Equation 14–60 represents the dominate eigenvalue and not the maximum; and the time-step restriction above is restated as:

(14–61)

This equation expresses the primary aim of automatic time stepping for 1st order transient analyses. The quantity Δtnλr appears as the oscillation limit output during automatic loading. The default is f = 1/2, and can be changed (using OSLM and TOL on the TINTP command). The quantity Δt1 is approximated as:

(14–62)

14.5.3.1. Additional Considerations for Thermal Transient Analysis

If the thermal solution reaches steady-state, the time step is allowed to increase by a factor of 2. Steady-state is determined if the maximum temperature change is less than 3 degrees for 3 time steps in a row. Use the OPNCONTROL command to modify these defaults, if desired.

For the Quasi option (THOPT command), the time step is also controlled by the property change tolerance such that the chosen time step keeps the property change below the tolerance. See the THOPT command for more details.

14.5.4. The Response Frequency for Structural Dynamics

The response frequency is used in the computation of Δt2 and is defined as (Bergan([105])):

(14–63)

where:

fr = response frequency (item RESFRQ for POST26 SOLU command and *GET command)
{Δu} = substep solution vector (tn-1 to tn)
[KT] = tangent stiffness matrix
[M] = mass matrix

This equation is a nonlinear form of Rayleigh's quotient. The related response period is:

(14–64)

Using Tr the time increment limited by the response frequency is:

(14–65)

When the middle step criterion is used, this criterion can be turned off.

14.5.5. Creep Time Increment

The time step size may be increased or decreased by comparing the value of the creep ratio Cmax (Rate-Dependent Plasticity (Including Creep and Viscoplasticity)) to the creep criterion Ccr. Ccr is equal to .10 unless it is redefined (using the CRPLIM command). The time step estimate is computed as:

(14–66)

Δtc is used in Equation 14–57 only if it differs from Δtn by more than 10%.

14.5.6. Plasticity Time Increment

The time step size is increased or decreased by comparing the value of the effective plastic strain increment (Equation 4–26) to 0.05 (5%). The time step estimate is computed as:

(14–67)

Δtp is used in Equation 14–57 only if it differs from Δtn by more than 10%.

14.5.7. Midstep Residual for Structural Dynamic Analysis

The midstep residual is used in the computation of Δtm. The midstep residual for the determination of the time step is based on the following consideration. The solution of the structural dynamic analysis is carried out at the discrete time points, and the solution at the intermediate time remains unknown. However, if the time step is small enough, the solution at the intermediate time should be close enough to an interpolation between the beginning and end of the time step. If so, the unbalanced residual from the interpolation should be small. On the other hand, if the time step is large, the interpolation will be very different from the true solution, which will lead to an unbalanced residual that is too large. The time step is chosen to satisfy the criterion set by the user (e.g. MIDTOL command).

Refer to the discussion in Newton-Raphson Procedure. The residual force at any time between the time step n and n+1 can be written as:

(14–68)

where:

ν = intermediate state between the time step n and n+1 (0 < ν < 1)
{R} = residual force vector
= vector of the applied load at n + ν
= vector of the restoring load corresponding to the element internal load at n + ν, which depends on the intermediate state of displacement at n + ν, and also the velocity and acceleration at n + ν. This intermediate state is approximately calculated based on the Newmark assumption.

A measure of the magnitude of {R} is established in a manner similar to the convergence check at the end of the time step (see Convergence). After the solution has converged at the end of the time step (n+1), the midstep residual force is compared to the reference value:

(14–69)

where:

||{R}|| = magnitude (vector norm) of residual force vector
Rref = reference force (see Convergence)

The convergence criterion for the midstep residual is defined by the value of τb (input as TOLERB on MIDTOL command):

If τb > 0, the value is used as a tolerance. If τb = 0 is specified or τb is not specified, then a default positive value is used as a tolerance. The midstep residual is assumed to have converged if its value is within the desired tolerance (ε τb ). Depending on how well the convergence criterion is satisfied the time step size for the next increment is increased or kept unchanged.

If the midstep residual hasn’t converged (ε > τb), the time step is repeated with a smaller increment:

(14–70)

where:

Δtn = old time step size
τb = midstep residual tolerance(TOLERB on MIDTOL command)

If τb < 0, the value is used as a reference force (reference moment is computed from reference force value) for midstep convergence check. A procedure similar to the one described above is followed with modified definition of time step size:

(14–71)


Release 18.2 - © ANSYS, Inc. All rights reserved.