4. Tune Penalty Stiffness to Suppress Contact Penetration Without Introducing Instability
Why Penalty Stiffness Is the First Knob to Turn
In a penalty-based contact formulation, LS-DYNA resists penetration by placing a spring between the slave node and the master segment. The stiffness of that spring is derived from the slave segment's bulk modulus, its characteristic area, and the model's critical time step. The default value is deliberately conservative: it is chosen so that the contact spring does not force the explicit time step below the structural one. The consequence is that a stiff, heavy, or fast-moving impactor can push through the master surface before the penalty force grows large enough to stop it. The visible symptom is a nonzero penetration depth in the contact gap output and a contact force history that lags the physical contact event.
The instinctive fix — multiply the penalty stiffness — works only up to a point. Because the penalty spring is integrated explicitly, its natural frequency is f = (1/2π)·√(k/m). Raising k raises that frequency, and once the contact spring's stable time increment drops below the global time step, the contact force begins to oscillate at the step frequency. That is contact chatter: a sawtooth force history, spurious energy in the contact energy channel, and sometimes a negative or wildly fluctuating sliding interface energy. The task in this case is therefore not "make the penalty as large as possible" but "find the smallest scale factor that removes penetration while keeping the contact response smooth."
The Test Model and Baseline Run
Use a deliberately small, fast-running model so the sweep is cheap. A rigid sphere of radius 10 mm impacts a 1 mm thick elastic steel plate clamped at its perimeter. The plate is modelled with 4 mm shell elements (Belytschko-Tsay, five through-thickness points), the sphere with a rigid mesh of comparable facet size. The sphere carries a prescribed initial velocity of 20 m/s along the plate normal; gravity is off. Contact is *CONTACT_AUTOMATIC_SURFACE_TO_SURFACE with the rigid sphere as slave and the plate as master, and SOFT = 0 so the penalty is the classic segment-based stiffness rather than a soft-constraint penalty.
Two output requests make the tuning measurable. First, add *DATABASE_ASCII with RCFORC to capture the resultant contact force on the master side. Second, request the global energy file GLSTAT and the contact energy in SLEOUT. Penetration depth itself is read from the ABSTAT/contact-gap output or, more directly, from the rcforc and gap data written when *CONTACT output is enabled with SPR and MPP flags left at defaults. A practical proxy that needs no extra output is the ratio of contact energy to internal energy at peak force.
Run the baseline with the default penalty scale factor PSF = 1.0 (the PSF field on the optional card of the contact definition, or the PENALTY parameter on the *CONTROL_CONTACT card depending on the keyword version). Record three numbers: maximum penetration depth, peak contact force, and the maximum contact energy ratio E_contact / E_internal. In a typical run of this setup the baseline shows visible penetration of a few tenths of a millimetre and a contact force that rises late and smoothly — no chatter, but the geometry is wrong.
Sweeping the Penalty Scale Factor
Now vary PSF across a logarithmic ladder and re-run each case. A workable ladder for this model is 1, 2, 5, 10, 20, 50, and 100. Keep everything else identical — same mesh, same initial velocity, same termination time — so that the only variable is the penalty stiffness. For each run extract the same three metrics. The table below shows the pattern you should expect; the exact numbers depend on your mesh and material, so treat them as the shape of the response, not as fixed targets.
| PSF | Max penetration (mm) | Peak contact force (kN) | E_contact / E_internal | Force history |
|---|---|---|---|---|
| 1 | 0.32 | 18.4 | 0.021 | Smooth, late rise |
| 2 | 0.19 | 21.0 | 0.024 | Smooth |
| 5 | 0.08 | 23.6 | 0.027 | Smooth |
| 10 | 0.03 | 24.9 | 0.029 | Smooth, slight ripple |
| 20 | 0.01 | 25.4 | 0.031 | Visible ripple |
| 50 | <0.01 | 26.1 | 0.038 | Sawtooth chatter |
| 100 | <0.01 | 27.0 | 0.052 | Strong chatter, energy drift |
The penetration falls monotonically as PSF rises, but the contact energy ratio and the force-history smoothness tell a different story. Between PSF = 10 and PSF = 20 the penetration is already negligible while the force history is still clean. Beyond PSF = 50 the contact energy ratio climbs sharply and the force trace develops the characteristic high-frequency sawtooth. The usable window is therefore roughly PSF = 10–20: penetration removed, chatter not yet introduced.
Reading the Result and Choosing the Value
The two curves in the figure cross in the region where penetration has already collapsed but the energy ratio has not yet started to climb. That crossing is the engineering answer: pick the smallest PSF on the flat part of the penetration curve, just before the energy ratio turns upward. In this model that is PSF = 10, with PSF = 20 as an acceptable upper bound if the impact velocity is later increased. Do not simply take the largest value that removes penetration — that is exactly the trap the next case in this series is built around.
Two sanity checks confirm the choice. First, the contact force history at the chosen PSF should be a single smooth pulse whose duration matches the physical contact time (roughly the time for the sphere to compress and rebound off the plate). A pulse that is shorter than the physical contact time, or that contains multiple spikes, indicates the penalty spring is ringing. Second, the global energy error — the difference between total energy and the sum of internal, kinetic, and contact energy — should stay below a few percent of the initial kinetic energy. A rising energy error at high PSF is the signature of an unstable contact spring injecting energy into the system.
PSF than the stable window allows, do not push the penalty further. Instead reduce the time step with *CONTROL_TIMESTEP (a smaller TSSFAC) so the stiffer contact spring remains stable, or switch to a soft-constraint penalty (SOFT = 1 or 2) which scales the penalty to the slave mass rather than the bulk modulus and is far less sensitive to impact velocity.Deliverable: Penalty Tuning Record
Keep a short record for each model so the choice is auditable and reproducible. The fields below are the minimum needed to defend the selected PSF in a review.
| Field | Example entry |
|---|---|
| Contact definition | *CONTACT_AUTOMATIC_SURFACE_TO_SURFACE, SOFT = 0 |
| Slave / master | Rigid sphere (R = 10 mm) / elastic plate (t = 1 mm) |
| Impact velocity | 20 m/s normal to plate |
| Baseline PSF | 1.0 |
| Baseline max penetration | 0.32 mm |
| Sweep range | 1, 2, 5, 10, 20, 50, 100 |
| Selected PSF | 10 |
| Penetration at selection | 0.03 mm |
| E_contact / E_internal at selection | 0.029 |
| Chatter onset | PSF ≈ 50 |
| Global energy error | < 2 % of initial KE |
With the record filled in, the tuning is complete: penetration is reduced by an order of magnitude relative to the default, the contact energy ratio stays near its baseline, and the force history remains a single clean pulse. The next step in a real project is to re-run the chosen PSF at the highest expected impact velocity, since the stable window narrows as velocity rises.