2. Definition of a problem adapted to transitory dynamics#

In this chapter, we will offer advice for developing the modeling of a mechanical problem for which we want to simulate a nonlinear transient numerical simulation with Code_Aster. The next chapter will deal with discretization in space and time, while in this chapter only the continuous model will be discussed.

The nonlinearities in transient dynamics with DYNA_NON_LINE are: * major movements and major transformations, * nonlinear constitutive relationships (plasticity, damage, viscoplasticity…) ), * interface nonlinearities: touch-friction, collisions…

2.1. Model#

The first step is the definition of a model consistent with the hypothesis of transitory evolution.

Thus, as any mechanical system should not have a mode with zero kinetic energy, it is important to ensure that the density is defined at every point of the continuous model. Likewise, if we want to introduce discrete elements to simulate solid bodies, we must associate a mass with them.

Likewise, for finite elements such as shells or plates, the user may have to ensure that all degrees of freedom, in particular rotations, have an associated mass term. To verify this condition, it is recommended to specify the option INER_ROTA =” OUI “when calling the [AFFE_CARA_ELEM] command for all COQUE models.

Some models that do not have a mass matrix (this calculation has not been programmed) cannot be used directly with DYNA_NON_LINE, such as non-local models. This difficulty can be circumvented by superimposing on them a classical model to represent kinetic energy (mass matrix): the two models based on the same nodes.

The use of certain common devices in quasistatics, such as very steep zones (fictional material with a very large Young’s modulus) to take into account reinforcements that one does not want to represent finely geometrically, for example, can cause disturbances in dynamics. In fact, this very stiff material can generate high-frequency oscillations and wave propagations whose speed is non-physical. Moreover, with an explicit time pattern, these very steep areas risk causing the value of the critical time step to drop (stability condition CFL, cf. [R5.05.05]).

In addition to the inertial contribution, the system will have a dissipative, and therefore damping, contribution. In non-linear mode, dissipation may be due to the behavioral relationship (plasticity…), to friction in the connections…

Usually, the detailed lack of knowledge of all the dissipative mechanisms in the system is compensated by a simplified representation that makes it possible to define global damping. In Code_Aster, we have two global viscous damping models (which can be associated in modeling with dampers on discrete elements, see command [AFFE_CARA_ELEM] (keyword factor DISCRET or DISCRET_2D, etc…): the Rayleigh model and modal damping.

In practice, in both cases, it is essential to have carried out a prior modal calculation. Indeed, modal damping is defined on modes that are then arguments of DYNA_NON_LINE. For Rayleigh damping, the only way to simply understand its physical meaning is to adjust it to modal damping values (which may come from experimental results). We’ll come back to Rayleigh damping in more detail in the next chapter Discretization of the ongoing problem.

In general, it should be borne in mind that all global damping models have been designed for linear cases and that they are then used to take into account certain dissipative phenomena such as friction or imperfect connections. Therefore, if we add to the model the consideration of these non-linearities, while keeping the global damping values adjusted to a linear approach, this can lead to too much dissipation in the non-linear system, especially at low frequencies.

Ideally, taking into account all non-linearities, combined with a fine discretization of the system, should make it possible to avoid having to use a fixed global damping in addition. In reality, we are forced to neglect certain non-linear aspects for reasons of cost CPU and/or lack of knowledge about the mechanisms implemented, and in this case, global depreciation has its role to play. The [U2.06.03] documentation provides more details on how to use it.

2.2. Boundary conditions and excitations#

The first decision to be taken to solve the transient is the choice of the configuration frame of reference (the field of movements): the general case is a Galilean frame of reference (absolute frame of reference). In frequent cases where the system is involved in a rigid translational movement, it is possible to solve the transient in a relative frame of reference. The calculation will be done on a built-in basis, see Mono-support. If the system is stressed by various supports excited differently, we are in case Multi-supports and we will have to solve the transient in absolute reference frame (unless the behavior is completely in HPP linear elastic).

In quasistatics, it is common to impose boundary conditions whose evolutions over time are functions that are simple to define, such as ramps or continuous functions that are refined piecewise.

In dynamics, it is common for excitations to be described by the product of a function of space and a ramp or function of time (or time signal). We will thus use the keyword factor: EXCIT.

It is recommended to set this signal in time so as to start and finish at zero. It may also be useful for the list of moments chosen for the transient with DYNA_NON_LINE to last longer than the imposed signal in order to let the system under study return to rest due to the effect of the damping that has been chosen.

In dynamics, if we want to avoid parasitic oscillations in the solution, it is absolutely necessary to impose quantities that are sufficiently regular in time: therefore at least continuously differentiable, piecewise. To achieve this simply, we can choose to define polynomial time functions instead of affine functions. Any irregularities in the imposed loads can be partially compensated by the use of an adapted damping model (or in particular thanks to an integration diagram in dissipative time of type HHT, see Time integration diagram). However, excessive nonlinear oscillations can irreparably compromise the rest of the calculations because the solution (including the internal variables) depends on its history.

In diagram u2.06.13-lissage we show a possible regularization of an affine curve (in dashed lines) by a curve adapted to dynamic calculation. We can note the need, on this example, to have to start dynamic evolution in a negative time, if we want to keep comparable values for positive moments (outside the regulation zones).

_images/Shape1.gif

In all cases where the spectrum of imposed signals can be calculated, the Shannon rule for defining sampling is not strict enough for nonlinear calculations. In practice, it is recommended to have between 4 and 10 steps of time in the shortest period of all the imposed signals.

Any discretization that is insufficient in time risks resulting in high-frequency over-oscillations of the solution. There is then a palliative to this by choosing certain schemes in time, but at the risk of introducing too much amortization into the system.

2.2.1. Mono-support#

Mono-support loads (driven by a rigid translational movement), which are commonly used during seismic calculations, are characterized by particularities that should be noted here.

First of all, the resolution, according to the single-support hypothesis, is done natively in the relative frame of reference linked to the support (write off…). In post-processing, it is therefore necessary to recombine with the training movement to obtain the absolute kinematic fields. For acceleration, it is generally immediate because training acceleration is the input data known in the form of the accelerograms used.

On the other hand, to reconstruct the absolute displacement, it is necessary to have the training displacement. In the classical case where we only have accelerograms, we have to integrate twice as much time to calculate the corresponding displacement (this integration can be done with the operator [CALC_FONCTION], keyword INTEGRE). The user must then be very attentive to the problems of the digital integration of sampled signals such as integration constants that can induce drift or oscillations at very low frequency. The [CALC_FONCTION] operator has options for correcting these drifts (keyword CORR_DEPL and « CORR_ACCE »), but it is very strongly recommended to carefully check the quality of the final integrated result because, in some cases, these corrections are insufficient, including the correction by high-pass filtering of [CALC_FONCTION] (with « CORR_ACCE », « METHODE = “FILTRAGE “ »). An additional specific correction is then essential. Going through FFT to do the frequency integration can then be a good alternative allowing to obtain integrated signals that are more correct or at least easier to correct (just a shift of a constant, for example, but without drift).

2.2.2. Multi-supports#

Among the types of loads that may require an adaptation to DYNA_NON_LINE, we can mention the case of multi-supports (in accelerations imposed on the supports) for seismic studies.

This method, originally developed for transient calculations under linear hypotheses, is based on the definition of static elastic modes at the supports. In non-linear mode, as long as these modes remain relevant, i.e. the nonlinearities are sufficiently weak to not require them to be updated, the usual multi-support approach remains acceptable.

On the other hand, in the general case, the use of static elastic modes calculated initially must be called into question. In this case, to avoid their use, multi-support solicitation must be taken into account by imposing on the supports the movements corresponding to the integration of accelerograms, which are input data for the usual linear multi-support method.

The method is then rigorous, regardless of the type of non-linearity, but the integration of the accelerograms must be carried out with care, in order not to result in incorrect signals.

Apart from sampling, it is essential, among other things, to check the absence of signal drift (a correction can be made with the option « CORR_ACCE » from [CALC_FONCTION]). In Code_Aster, you can easily do this integration using the « FFT » that is available in [CALC_FONCTION]. It is also an opportunity to analyze the frequency content of the signal. We then have the relationship, in the frequency domain: \(U=-A/{\omega }^{2}\).

You can also use the « INTEGRE » option from [CALC_FONCTION]. For more information on the implementation of single or multi-support requests, reading the documentation [U2.06.09] and [U2.06.10] is ideal.

2.2.3. Contact and collisions#

The charges resulting from contact by penalization can also disturb the solution by high frequency oscillations, which are linked to the value of the penalty coefficient. These artifacts can be reduced with the use of a high-frequency dissipative time integration scheme.

We can also have a list of calculation moments that manage difficulties by division, cf. the [DEFI_LIST_INST] command with the keyword factor « ECHEC = _F (EVENEMENT = “DELTA_GRANDEUR”,…) », of the « event-driven » type, by testing whether the increment of a component of a field (speed for example at the point of contact) exceeds a fixed threshold.

If you want to use an explicit time integration diagram, an excessively large penalty coefficient will cause the critical time step to fall. Compared to statics, it may prove mandatory to lower the penalty coefficient, with the disadvantage of increasing the interpenetration of solids during contact.

You can consult the test cases: * sdnv100, impact of a beam on a rigid wall [V5.03.100], * sdnv105, balancing a block on a table [V5.03.105], * sdnv301, radio collision of two elastic balls [V5.03.301].

note

For explicit time diagrams () Time integration diagram), if we Imposes boundary conditions on the move that change over time, you have to keep taking into account the fact that these conditions are in fact imposed in acceleration in explicit (because it is the primal unknown). That means we have to enter DYNA_NON_LINE the derivative second of the moving signal that you want to impose. This evolution of displacement Imposed must therefore be differentiable at least twice in time.

2.3. Initial conditions#

In the same logic as for boundary conditions, it is recommended to avoid singular initial conditions, such as what imposed loads that evolve over time such as Dirac or Heaviside can generate.

We will favor regular conditions therefore with derivatives (first and if possible second) in zero \({t}_{0}\) time.

If you want to start from a prestressed (non-virgin) state, it is entirely possible and recommended to combine a dynamic calculation to a quasistatic calculation. In this case, the initial state (keyword keyword factor ETAT_INIT) must then be perfectly balanced (displacements, constraints, internal variables, forces applied at this moment) and, subject to regularity in the conditions imposed at the time of the static to dynamic transition, the dynamic excitation having to start with a zero signal, the transitory solution must not present non-physical oscillations. This method of linking quasistatics to dynamics allows easy consideration of all preloads, such as gravity. For example, in a release simulation, with an initial displacement imposed that is not zero, then you have to fill in the constraints and internal variables corresponding to this moment.

If one had wanted to take into account this initial prestressed state with only dynamic resolution, it would have been necessary to introduce a high initial damping and wait for the oscillations to have dissipated in order to have the prestressed static state. Then, we can continue the calculation with physical amortization. One is obliged to proceed in this way when one uses a calculation code that does not treat nonlinear quasistatics (cf. Europlexus, see [bib2] _).

In certain situations, it is necessary to combine a simulation in statics, with blocked boundary conditions, then dynamically under the action of an incident wave field with boundary conditions that are absorbent and released. We will consult [U2.04.08].

During the initialization of the integration diagram for the nonlinear transient calculation, an attempt is therefore made to invert the mass matrix. If it is singular, then a message alerts the user and the initial acceleration is arbitrarily set to zero. We can then fear the presence of parasitic oscillations at the beginning of the transient; we can reduce them by adopting a smaller time step at the very beginning of the transient. In general, the non-inversibility of the mass matrix should lead the user to check his model, unless this is voluntary.

During prosecutions, two aspects are sensitive.

On the one hand, it is better to avoid having too sudden variations in time steps: the coherence of the solution during the pursuit would be affected. A strong change in time steps can be seen as a filter. A large time step constitutes a low-pass filter of a solution initially calculated with a fine time step.

On the other hand, if you want to change the type of schema over time, certain rules must be respected. If we want to go from an implicit time integration method (quasistatic or dynamic) to an explicit diagram, the pursuit will be mathematically valid because the initial state will be balanced (to the nearest residue). On the other hand, the direct inverse flip-flop introduces an error because the initial state, coming from an explicit calculation, will not respect the balance in the implicit sense. In fact, this operator will seek to recheck the balance: the equilibrium is solved by inverting the mass matrix, which gives the initial acceleration. If this field is non-zero, this reflects the non-verification, in the static sense, of the balance, at the initial moment. In practice, this imprecision can cause oscillations in the solution (overshooting phenomenon). In a more subtle way, switching from the \(\theta\) -schema to a second-order schema like the average acceleration « NEWMARK » or « HHT » is going to introduce a small error on the acceleration terms, which can disrupt the numerical solution. If, in spite of everything, we want to make these flip-flops that induce errors, we can minimize the inaccuracies by choosing to make these flip-flops during phases where the excitation and therefore the solution does not change much. In summary:

Tableau 2.1 Enchaînements de calculs transitoires#

Origin

Destination

Possible

STAT_NON_LINE

DYNA_NON_LINE

yes

DYNA_NON_LINE with implicit time schema

DYNA_NON_LINE

yes

DYNA_NON_LINE with explicit time schema

DYNA_NON_LINE with explicit time schema

yes

DYNA_NON_LINE with explicit time schema

DYNA_NON_LINE with implicit time schema

no