Transient analysis using contact elements

[miran]

Dear STKO Team,

I read in a previous post that you are preparing a tutorial on soil-pile interaction. Will you including transient solutions?

Meanwhile, I have a few issues with my soil-pile analysis that I hope you can comment on. An appoligy in advance for what became a long post


1) Issues with small (!) time step in implicit integration scheme.
I have a pile-soil model (based on your example, but with absorbing boundaries), where the interaction is modeled with the ContactBeam3D elements. Load is applied at the pile top.

Static analysis works well.

When I run my model using a implicit time scheme with a time step of 0.005 s (Modified newton, TRBDF2), it seems to produce fairly good results. When I decease the time step to 0.0005 s, the nodes corresponding to solid elements that are in contact get very large and very noisy deformations, but the deformations in all other nodes in model, including the beams (piles), remain seemingly reasonable.

I have tried both time steps for a variety of soil and beam stiffness's (elastic), and the same thing happens, i.e. good results for dt 0.005 s and large and noisy results for dt 0.0005 s.

As far as I can see, the literature is pretty consistent with the fact that contact problems should benefit from small time steps due to the high non-linearity of the contact formulation (implicit schemes included). This is also very logical to me, and thus my problem is seems somewhat strange.

Am I missing something obvious here, and is there a logical explanation to why this occurs?

NB! In order to decrease the non-linearity that occurs from contact, I decreased the stiffness parameter G by a factor of 10 in ConctacMaterial3D. The idea was to prolonged the elastic stick and thus smooth out some of the non-linearity. It did not help. I also tried increase G just to see what happened, and the results got even worse, which I guess was to expect.

2) Issues with explicit integration scheme.
Since I have a relatively small time step and small a small time frame (1 s), an explicit scheme should be fine. I used the central difference method with time step 0.0005 s. Element sizes and time steps are well below critical values. Everything else is as for the implicit scheme.

However, this does not work at all. The displacements go towards infinity, and is impossible to visualize anything in the post-processor. Explicit schemes are common in high-velocity contact-impact problems, so it seems strange that it does not work here.

Do you agree the explicit schemes should work if all the critical values and parameters are as they should be? What could the problem be?

2) Issues with ContactBeam3Dp
It is my understanding that the contact problems constrained with Lagrange multipliers tend to create noise for high-velocity problems. I do not consider my problem high-velocity in that sense, but I thought it might be worth while to try ContactBeam3Dp. I used your example model posted in on of the previous posts, and simply changed the contact element from ContactBeam3D to ContactBeam3Dp. I am not an expert on choosing penalty values, but I could not get this element to work for a variety of penalty values. The analysis runs, but the results are terrible in the sense that beam displacements are off the charts. There is no information on this element in OpenSees Wiki as far as I can see.

Do you have experience with this element, and is there something in particular to be aware of?

[STKO Team]

Dear Miran,

I've read something in this regard, but did not test it deeply.
On the website there is this note about instability caused by the Lagrange Multipliers:

In transient analysis, the presence of the contact constraints can effect the stability of commonly-used time integration methods in the HHT or Newmark family (e.g., Laursen, 2002). For this reason, use of alternative time-integration methods which numerically damp spurious high frequency behavior may be required. The TRBDF2 integrator is an effective method for this purpose. The Newmark integrator can also be effective with proper selection of the gamma and beta coefficients. The trapezoidal rule, i.e., Newmark with gamma = 0.5 and beta = 0.25, is particularly prone to instability related to the contact constraints and is not recommended.

Regarding the penalty version: I have never tested it, but a common rule for penalty paramter is:

• Evaluate (approximately) the largest entry on the global (assembled) stiffness matrix diagonal... not the exact value, just an approximation is fine. Call it Kmax
• Get the order of magnitue: O = ceil(log10(Kmax))
• Penalty = 10^(O+8). Where 8 is half the significant number of digits for a double-precision floating point number

Side note: When using penatly approaches, either in contant or in the constraint handlers, NEVER use a test on the residual (Unbalance)!

[miran]

Thank you for the reply,

I have in fact used the TRBDF2 scheme as OpenSees recommends, and it does seems to work for the larger time step. However, small time steps produce large noise in the displacements. This is very strange.

Regarding the convergence criteria.

Like you stated, the right hand side (residual) of the equilibrium (static or dynamic) equation will contain penalty values and is clearly not suited for convergence. However, when you apply Lagrangian constraints, the Lagrangian multipliers appear both in the right hand side of the equation (residual) and on the left hand side as the Lagrangian field. Since the Lagrangian multipliers equal the traction (normal and tangential nodal forces with respect to contact plane), this should be fine for linear systems, but I am not sure how straight forward the implementation is for high non-linear problems such as contact formulations.

I encourage anyone who has any experience with transient analysis in combination with contact formulations to please comment this thread. It would be much appreciated.

[STKO Team]

Dear milan,

Regarding the use of the lagrange multipliers, you are right, I never use lagrange multipliers in OpenSees for one main reason: When OpenSees does a test on the displacement or residual vectors, it takes the norm of the entire vector, and I don't like it, because it is not consistent... I mean, with lagrange multipliers you put forces in the displacement vector (lagrange multipliers) and displacement constraints (in case of contant the gap) in the residual... when you then take the norm of the entire vector, you are mixing things with different dimensions...

If also the TRBDF2 shows the same issue with smaller time steps, then probably it is not due to stability issues of the integrator scheme...

There is something else worth trying... Are you doing a quasi static analysis (transient with very large dT) for gravity loads, before the fast dynamic? I have seen many examples where people do slow (quasti static) transient analysis (let's say dt = 500s) for gravity, and then they suddenly run a fast dynamic (say dt = 0.0001s)... Is this the case?

In dynamic analysis it is not a good idea to drasticallyl change the dT, since the time integration schemes typically assume a regular grid...
If this is the case, I can send you a custom TCL script (I did it for some works in the past) where you do "relaxation" analysis (no loads) between the gravity (large dT) and the fast dynamic (small dT). In this "relaxation" analysis the TCL script will basically run N steps with a smoothly decreasing dT, just to damp out any undesired acceleration.

[miran]

Regarding static analysis. No, I did a linear static analysis with one increment for the constant gravity loads (linear because the non-linear contact elements are added at a later stage and all other elements are linear).

Also, I used the initial state analysis, so that the displacements before the transient analysis are zero.

I have tried to run the dynamic analysis both with and without a gravity analysis. The same issues occur.

Basically, this is the implementation in short:

1) Add solids (soil), beam (pile), appropriate fixed constraints and gravity loads.
2) Run initial state analysis (static, one pseudo time step)
3) Run another static analysis (tried this according to your video SSI video)
4) Remove prior constraints and add permanent constraints, LK boundaries, dynamic loads and reaction forces.
5) Add contact elements.
6) Add damping.
7) Run dynamic analysis

I have also noticed noise at the end of the analysis (regardless of the time step). Some time after the dynamic loading stops and everything is damped out, some substantial acceleration and velocity values can be seen at some of the solid elements connected to the slave nodes. See screenshot below.

I am sorry to bother you with these issues since this is not an STKO-problem. I am not expecting this to solved here, but any comment or tip is very much appropriated.

[STKO Team]

Dear Miran,

Can you post your SCD file (zipped) here? or send it by e-mail. So that we can have a look and do some debug in the opensees code to unserstand what's going on

[miran]

Dear STKO Team,

Did you get chance to have a look?

[STKO Team]

Dear Miran,

I confirm that it's not easy to get rid of those spurious oscillations at the contact nodes due to lagrange multipliers. Reading some literature, it is a well known problem. I am now testing your model with the BeamContact3Dp (penalty). At a first try, I noticed that the spurious oscillations are gone (when using a relaxation intermediate analysis). However another problem occurs (conpenetration sometimes fails). I will give you some updates tomorrow

[miran]

Sure, thank you!

Clearly, Lagrange works fine for static problems. From what the literature tells us, I was expecting noise for high-frequency excitation. However, it seems noise occurs for dynamic problems in general when using implicit integrations.

Contact formulations, both penalty and Lagrange, will cause problems in the implicit integration scheme due to high non-linearity in combination with Newtons algorithms since we suddenly add/remove terms in the K-matrix. In my opinion, it is therefore logical use explicit solutions. Even more so when the duration of load is relatively short and we already use a small time step. I will play around some more with this and let you know if I get decent results.

[STKO Team]

This is for sure, Contact conditions leads to strong discontinuities in the solution, very challenging for implicit schemes. We are currently implementing another node-to-node contact element in OpenSees, that uses the so-called IMPL-EX scheme at the constitutive level. IMPL-EX is a mixed implicit explicit integration scheme at the constitutive level, that basically turns a nonlinear problem into a step-wise linear one. As soon as it is ready I will tell you