3DEC Modeling • Problem Solving with 3DEC

Reaching Equilibrium

At the end of setting up initial conditions, at the end of the model, and perhaps during multiple intermediate stages, the model must reach “equilibrium.” In the examples so far, this is done with the simple command model solve. But how do you know when you have reached equilibrium? How can you tell if you have performed more cycles than is necessary, or if some specific part of the model still needs computation time?

Calculation cycles in 3DEC can be performed manually with the model step or model cycle command. There are occasions where this is useful, most often when you have prescribed a velocity at a boundary and want to step until a specific total displacement has been reached. Generally, however, you wish to cycle until a certain condition has been reached—for example, a dynamic simulation could cycle until a certain amount of time has passed. The model solve command allows cycling to continue until one or more limits (specified using subsequent keywords) has been met. In a static model, the most common limit of interest is “until the model is in equilibrium.”

Because 3DEC uses an iterative procedure informed by physics to reach static equilibrium, when we refer to being “in equilibrium,” we mean the same thing as being fully converged to a static solution. 3DEC provides a number of different ways to measure and visualize how close you are to being converged or having reached an equilibrium state in your model. But, as always, the more you are aware of the issues involved, the more efficient and accurate your models are likely to be.

Convergence Criteria

3DEC has four different criteria available to determine when equilibrium is reached. Local Force Ratio, Average Force Ratio, Maximum Force Ratio, and Convergence. Each of these options is described below.

Local Maximum Force Ratio

The local force ratio is a unitless measure of local convergence at a gridpoint. The idea is simple: divide the out-of-balance component remaining after summing forces by a measure of the total forces being applied to the gridpoint. For performance purposes, we use the “Manhattan norm” or “Taxicab norm” of a vector, defined as:

\[\langle v_i \rangle = |v_x| + |v_y| + |v_z|\]

or the sum of the absolute values of the components of the vector. The local force ratio can then be defined as:

\[\frac{\langle \sum_i f_i\rangle}{\sum_i \langle f_i \rangle}\]

where \(f_i\) are all the force vectors being added to a gridpoint.

The local maximum force ratio is the maximum local ratio of any gridpoint in the model. This is a very conservative measure of the convergence of a model. However it has a weakness when there is less than two forces being applied to a gridpoint, or when all forces being applied to a gridpoint act in the same direction. In this circumstance, the only “equilibrium” value possible is for both numerator and denominator to converge toward zero; the ratio between them will tend toward 1.

3DEC will detect when only one force vector is being applied to a gridpoint and automatically remove it from consideration. But if multiple aligned forces are acting on a gridpoint, or if the second force is negligible compared to a main driving force, then this check is insufficient. These points will dominate the maximum local force ratio, and the model will seem to never effectively converge.

In practice this does not happen very often. The majority of models never see this effect. But it does happen, and so maximum local force ratio is not used as the default convergence criteria for model solve in order to guarantee that any given model will converge using the default criteria.

The last local force ratio calculated for a gridpoint is stored and is available to FISH as block.gp.ratio and to plot items as the Local Force Ratio.

Use the keyword ratio-local with the model solve command to use the local ratio for the solve criterion.

Average Force Ratio

The average force ratio is defined as the sum of all out-of-balance force components at every gridpoint divided by the sum of all total forces applied at a gridpoint. Or:

\[\frac{ \sum_j^n \langle \sum_i f_i \rangle }{ \sum_j^n \sum_i \langle f_i \rangle }\]

This measure serves as a reliable measure of the overall convergence of the system. For models that are relatively uniform, it performs very well. It is bulletproof enough that the default convergence criteria for model solve is ratio-average 1e-5.

However, for very non-uniform models, the effects of localized convergence problems can be lost in the average, making the user unaware that an area that might be of most interest still needs computation. This can occur either because of disparities in zone size, or because of large differences in stiffness. For large complex models, this measure should be checked carefully before being trusted. See the following discussion.

Use the keyword ratio-average or simply ratio with the model solve command to use the average ratio for the solve criterion.

Maximum Force Ratio

The maximum force ratio is defined as the maximum out-of-balance force divided by the average total force acting on all gridpoints. The average total force is defined as:

\[\frac{\sum_j^n \sum_i \langle f_i^j \rangle}{N}\]

where \(N\) is the total number of gridpoints in the model.

Use the keyword ratio-maximum with the model solve command to use the maximum ratio for the solve criterion.

Convergence

The convergence solve was created as a means of giving the user more control over the convergence criteria desired and overcoming the possible disadvantages of relying on the local force ratio.

First - a new gridpoint field variable called the Local Force Target is created. This value defaults to 1e-4 in every gridpoint and node. It can be changed with the block gridpoint initialize ratio-target command and with the block.gp.ratio.target FISH intrinsic. This value indicates the local force ratio that is to be considered “converged” for that specific gridpoint. This means that the local convergence criteria can vary in different locations of the grid.

We define a new quantity, unimaginatively called the convergence, that is the ratio of the current local force ratio to the target local force ratio. A convergence of 1.0 is therefore considered “converged”. The maximum convergence is then a measure of overall model convergence tied to the convergence keyword in model solve. The convergence value is available to plot items, and to FISH as the block.gp.convergence intrinsic. By default, a convergence of 1.0 is exactly the same as a ratio-local limit of 1e-4.

In the event that the local ratio calculation is failing in corner cases, the ratio-target can be set to a high value in that gridpoint or node, effectively removing them from the overall convergence criteria. In addition, it is not unusual for a small region of the model to have an overly large effect on the time it takes to reach convergence—either because they are very small, or because they are very stiff, or both. For example, it is not unusual for a small number of very stiff structural supports to dominate the convergence time of a very large model. The ratio-target for those nodes could be set to 1e-3 or even higher, relaxing the criteria for those nodes once it is determined that this will have negligible effect on the model results of primary interest.

Use the keyword convergence with the model solve command to use the convergence criterion.

Choosing Convergence Criteria

The default convergence criteria if none is specified for a model solve command is a ratio of 1e-5. In our experience, an average force ratio of 1e-5 corresponds to a local force ratio limit of around 1e-4 for a model of average complexity. This is why 1e-4 was chosen as a default ratio-target.

A maximum local ratio of 1e-4 is generally sufficient, even conservative, for most engineering applications. A maximum local ratio of 1e-3, or even as high as 1e-2 may be sufficient for many purposes. This is especially true for intermediate stages (excavation sequences, etc.) of a model, where taking the time to reach “full” convergence after every step may have negligible effects on the final state of a model.

Evaluating Equilibrium

For the purposes of illustration, we will demonstrate with a cartoonishly simple model of a surface excavation. The basic geometry is shown in Figure 1.

Figure 1: Cartoon model of a trench excavation.

The full data file for the first demonstration is below.

model new
model large-strain off

; Create trench
block create brick 0 30
block cut joint-set origin 0,0,24 join
block hide range position-z 0 24
block cut joint-set dip 90 dip-direction 90 origin 3,0,0 join
block hide range position-x 3 30
block cut joint-set dip 90 origin 0,6,0 join
block group 'Trench' range position-y 0 6
block hide off

; zone the model with smaller zones near the trench
block zone size isolated-point 0 0 30 edgelength 0.5
block zone generate-new minimum-edge 0.5 maximum-edge 5 gradation-surface 0.1

; Constitutive model and properties
block zone cmodel assign mohr-coulomb
block zone property bulk 3e6 shear 2e6 density 1000
block zone property cohesion 1e5 friction 25 tension 1e3

; Boundary conditions
block gridpoint apply velocity-x 0 range position-x 0
block gridpoint apply velocity-x 0 range position-x 30
block gridpoint apply velocity-y 0 range position-y 0
block gridpoint apply velocity-y 0 range position-y 30
block gridpoint apply velocity-z 0 range position-z 0

; Initial conditions
model gravity 0 0 -10
block insitu topography ratio-x 0.5 ratio-y 0.5
model solve

; Histories
model history mechanical ratio-average
block history displacement-x position (3.1,0,30)
block history displacement-z position (3.1,0,30)

; Excavate
block delete range group 'Trench'

; Solve
model solve
model save 'simple'

;; Solve with convergence
;model solve convergence 1
;model save 'simple-convergence'

There are several things that should always be done to check to see if your model has actually reached sufficient equilibrium. First, take histories of displacements and/or velocities in regions of primary interest. If these displacements are still changing relatively rapidly when cycling completes, the model is not done adjusting. Second, make plots (contours and isosurfaces) of Local Force Ratio or Convergence. An isosurface plot of Convergence is particularly useful to tell you where the model is taking the longest to converge, and where you should look closer in a contour plot (perhaps using cut-planes to see into the interior of the model).

Figure 2 shows the histories of averaged force ratio and displacements in the top of the trench. Note that the average force ratio is plotted on a log scale. For a well-behaved model, this sort of linear logarithmic (or exponential) convergence is expected. Failure or other non-linear behavior will show disturbances in convergence, and progressive uncontained failure will make the curve level out. The displacement history, again, shows a typical slightly underdamped response. It is clear that the top of the trench wall has stopped moving significantly and the model has effectively reached equilibrium.

Figure 2: Histories of average unbalanced ratio and displacement for the simple trench.

Figure 3 shows the result of examining contours of convergence. While the bulk of the model has a convergence well below 1, a few gridpoints near the surface are still much higher. This is fairly typical of this sort of problem; the surface gridpoints have the same inertial mass resisting motion and smaller driving forces due to stress. The gridpoints at the bottom tend to (and need to) converge first, and small adjustments propagate against gravity, so the surface gridpoints converge last. If a stricter convergence criteria of convergence 1 was used, the model would take half again as many steps to reach equilibrium. How important this is depends on the significance of relatively small adjustments in the surface gridpoints. If exact surface settlement values were critical, this might be a necessary change.

Figure 3: Contours of gridpoint convergence.