Bulk Viscosity for Shells and Solids in LS-DYNA

In any wave propagation code, such as LS-DYNA, which belongs to a family of ‘Hydro’ codes, bulk viscosity is essential to treat shocks. Smooth initial data can lead into shock discontinuities and if left untreated can result in severe instabilities. LS-DYNA has the capability (performed by default) to automatically detect the shocks and treat them by adding a pressure term that is based on the element dimension, density, strain-rate and some user-defined coefficients (Q1 and Q2). The addition of the pressure term guarantees unperturbed solution away from the shock and satisfies Hugoniot jump conditions that requires the conservation of mass, momentum and energy across the shock front. The automatic detection of the shock was originally performed by looking at the divergence of the velocity field which is now replaced with the trace of the strain-rate tensor in multidimensional problems. The divergence of a vector field such as a velocity shows us the amount of flux or density that is flowing in (sink or compression or convergence) or out of (source or expansion or divergence) of a certain area. The divergence of a vector field provides us with a scalar function whose values are positive for out-flow (source, expansion,divergence) and negative (sink, compression, convergence). Fortunately, in multi-dimensions, the trace of a strain-rate tensor provides us this information very easily which is computed every cycle by LS-DYNA. If the trace of the strain-rate tensor is negative (compression, convergence, sink), then LS-DYNA automatically adds a pressure them that is based on the element dimension (square-root of Area for 2D and cube-root of Volume for 3D) while a positive trace is ignored.

Bulk Viscosity Coefficients
LS-DYNA provides two types of bulk viscosity coefficients namely Q1 and Q2. Q1 is called as the quadratic term that helps to smear the shocks and also helps in preventing the element from collapsing under high velocities where the particle velocity is exceeds the sounds speed for the material. Q2 which is called as the linear term, helps to rapidly damp out the oscillations, often called as the ‘ringing’. By default, these coefficients are fixed at 1.5 (Q1) and 0.06 (Q2) and are both active for solid elements. Starting 970, LS-DYNA allows the application of the bulk viscosity pressure for shell elements which can be invoked by setting TYPE = -2 in *CONTROL_BULK_VISCOSITY. Use of bulk viscosity for shells and solids are highly recommended to improve simulation stability. The internal energy that is dissipated by viscosity is computed and included in the overall energy balance.

Simple Elastic Impact Example
To demonstrate the effects of the adding the bulk viscosity, two slender bars are set up to impact against each other as shown in figure below. The bars are meshed using shells and a single element away from the point of impact is monitored for the uni-axial stress time history. The figures below show the stress time history for no bulk viscosity, quadratic bulk viscosity and quadratic+linear bulk viscosity.

(Click image to enlarge)

(Click image to enlarge)

(Click image to enlarge)

(Click image to enlarge)

References:
1. LS-DYNA Theory Manual, 2007, Livermore Software Technology Corp.
2. LS-DYNA Keyword Manual, 2007, Livermore Software Technology Corp.

Acknowledgement:
Jim Day for reviewing this post and providing valuable comments and suggestions

Limitations of Using LCSR for Strain-Rate Inclusion in MAT_024

The widely popular material model *MAT_PIECEWISE_LINEAR_PLASTICITY (MAT_024) supports various methods to include the strain-rate effects. One of the methods is the scaling of the quasi-static stress-strain curves using a load-curve, LCSR, that defines a scale-factor as a function of strain-rate. This method works by first finding the yield-stress, SIG_QUASI_STATIC, as a function of the effective plastic strain, that may be defined using a load curve LCSS or using discrete points using ES-EPS, and then multiplies this value by the scale-factor by looking up the curve LSCR to give a ’scaled’ yield stress, SIG_SCALED. Plasticity of the material is then determined by comparing the trial stresses with the scaled stress, SIG_SCALED. If the trial stress is greater than SIG_SCALED, then the plastic strain increment is computed and the trial stress is scaled back to the yield surface using the radial-return algorithm. The primary limitation of this method is in its assumption that the stress-strain curves at higher strain-rates is a simply y-scaling of the quasi-static curves which is not the case in many materials. It is well known that at higher strain-rates, the materials not only show an increase in the yield strength but they also exhibit brittleness that is not captured by the LCSR method. This is graphically demonstrated in the figure below. The black curve is the quasi-static curve that can be input using either EPS-EPSS or LCSS parameters. The dark colored curves, red and blue, show the scaled curves as a function of strain-rates using the LCSR parameters. The lighter counterparts of red and blue curves show the actual behavior that is often seen in the physical tests where some amount of brittleness is observed.

(Click image to enlarge)

Modeling Loading and Unloading Behavior in Seatbelt Materials

Seatbelt constitutive model, invoked by using *MAT_SEATBELT, in LS-DYNA provides features to model the loading and unloading characteristics from a uni-axial test. Parameter LLCID provides ability to model the loading curve which allows the definition of force as a function of engineering strain. Parameter ULCID, provides ability to model the unloading curve again allowing the definition of force as a function engineering strain. Inherently, LS-DYNA models the seatbelt as an elastic-plastic where the first non-zero point on the loading curve, LLCID, is taken as the yield force/stress. Unloading always following the loading curve when the force level is below the first non-zero value. Unloading int the post-yield region is done by first shifting the unloading curve, ULCID, along the x-axis (striain) until it intersects the current yield point on the loading curve before unloading along the unloading curve. This is shown in the following figure.

(Click image to enlarge)

In certain cases when the unloading is extremely nonlinear, then the shifting of the unloading curve could yield conditions where the unloading curve could intersect the loading curve at more that one location. In such cases the unloading occurs as shown below.

(Click image to enlarge)

Identifying Necking in Metals and Plastics

When characterizing materials such as Metals and Plastics in LS-DYNA, most constitutive models provide a yield criteria that accounts for a 3D state of stress which reduces to a uniaxial yield stress in 1D. This allows us to directly input the true stress-strain curve from a one-dimensional state of stress testing such as in uniaxial testing which then can be used for any 3D problems. One important consideration in the usage of such one-dimensional testing data is the phenomenon called necking which occurs in both Metals and Plastics. It is well known that post-necking stress state is no longer uni-axial but is a 3D state of stress. This limits us to use uniaxial information only upto necking and use some sort of iterative process to characterize the post-necking behavior. The first task therefore is to identify the necking point on a uni-axial true stress-strain curve. This is done by plotting the true stress-strain curve against its derivative. The intersection point of these curves is the necking point. When this curve overlay is done using the raw test data, which is sampled at a high frequency, the derivative may be extremely noisy and this may provide multiple intersection points. A sample of such noisy overlay is shown here:

(Click image to enlarge)

To overcome this, it is usually recommended to digitize the raw curve using anywhere from 10-100 equi-distant points which provides a very smooth derivative of the true-stress-strain curve ensuring a unique point of intersection (necking point). A digitized curve using the raw test data is shown here.

(Click image to enlarge)