Simulations of a stretching bar with the STZ model

This webpage describes some preliminary results of simulations carried out at Santa Barbara, using the STZ description of amorphous plasticity to model a 2D bar being stretched.


The simulations are based on a two dimensional elastoplastic model. The horizontal and vertical components of the velocity are denoted by u and v, and the stress tensor is

The stress tensor

where p is the pressure, and s and τ are the deviatoric components. The deviatoric stress tensor, where the pressure is subtracted, is denoted by σ0. In addition, we also consider the effective temperature χ, which is an internal measure of molecular disorder that affects the total amount of plastic deformation. In rescaled units, the system of equations can be written as

The rescaled STZ system of equations

For details about rescaling, see this PDF document. In this system, the yield stress of the material is scaled to unity. The function q is defined as

The STZ equation for q

so that plastic yielding only sets in when the deviatoric stress exceeds the yield stress. Based on values for real materials, we set μ=30, c0=1, K=60, and χ=0.15. This leaves us with three parameters to vary:

Computational implementation

The simulations are carried out using a finite difference method on a square grid, and the details are similar to those in previous simulations. The boundary of the bar is tracked using a level set library provided by Frederic Gibou. Finite difference calculations are carried out within the object, and when the finite difference calculations need to reference values outside of the object, these values are extrapolated out from the defined field values.

A very small amount of viscosity is employed in the simulations to help remove unwanted high frequency elastic modes. The inclusion of these terms is also partly based on physical reasoning, since these high frequency modes are damped out in practice.

The level set routines used here, like many others, have some problems resolving sharp corners in the boundary. The level set tends to smooth out any angles before they can appear, effectively creating a small surface tension. Thus, in the current implementation, it may be impossible for us to properly simulate a crack.

A typical simulation run

The simulations are based on a rectangular bar in the region |x|<0.5, |y|<0.25 with a small dimple located at (x,y)=(0,0.25). The bar is attached to vertical walls at x=0.5 that move to x=0.9 at constant velocity. The image below shows a typical simulation in progress, using parameter values of ν=1000, χ0=0.074, with a bar velocity of 0.004.

Shear bands in a typical simulation

Shear bands, where the effective temperature is higher, nucleate diagonally from the dimple. The bar deforms in a ductile fashion and eventually splits into two. See below for movies of this situation.

Varying the size of ν

Three simulation snapshots for different values of nu

The three images to the right are snapshots taken at t=40 using ν=10 (top), ν=1000 (middle), and ν=104.5 (bottom) with all other parameters the same as above. Despite changing ν by more than three orders of magnitude, the simulations exhibit qualitatively the same behavior, with diagonal shear bands forming from the dimple. The onset of the shear band formation is only slightly different between the three simulations. At late times, the ν=10 simulation exhibits different behavior, with the triangle of material between the shear bands moving upwards as an almost-rigid region. This may be indicative that the bar may want to crack along the shear bands, but that the numerical method is incapable of resolving this. Without cracking, the triangle of material is forced to slide upwards.

Varying other parameters

Other simulations were carried out to investigate the effect of the initial temperature. The movies below show two simulation runs with ν=10, for initial temperatures of χ0=0.074 (top) and χ0=0.068 (bottom). The run with lower initial temperature exhibits more crack-like behavior, with the shear bands taking slightly longer to form, but there is no significant qualitative change.

The effect of the strain rate was also investigated, and the movies below show the typical process (top), a simulation five times as slow (middle), and a simulation ten times as slow (bottom). For each simulation, the time is rescaled so that the walls appear to move at the same speed. The slower simulations exhibit more ductile behavior, with the bar taking longer to split into two.

An off-center initial dimple was also considered. Despite the break in symmetry, the shear bands still both form at approximately the same angles.

Higher resolution

A higher resolution simulation

The image above was taken at t=40 in a simulation with χ0=0.074 and ν=10, where the grid spacing was halved. Reducing the grid spacing by a factor of two requires that the time interval be reduced by a factor of four. Thus the simulation takes sixteen times as long to run than the previous ones. The simulation provides cleaner results, and it is clearer than ever that the triangle of material between the two shear bands behaves like a solid body. However, the method is still unable to resolve any possible cracking. Since the computation time scales like the fourth power of the grid spacing, it suggests that in order to resolve cracks, we may have to make use of a more advanced technique such as an adaptive grid, which will only resolve more detail in the vicinity of the crack.

An alternative geometry

The simulation below is of a tall bar of width 2.5 being stretched. For these simulations, the bar is allowed to slide up and down the wall, as opposed to being fixed as in the previous simulations; this resolves some problems at the join between the wall and the edge of the bar.

A tall simulation

Shear bands nucleate from the initial dimple, but now appear to reflect off the vertical walls forming a criss-cross pattern. The bands also appear to become more diffuse, and it is unclear what causes this.

Future work