- Email: [email protected]
A reduction method for nonlinear structural dynamic analysis
COMPUTER METHODS NORTH-HOLLAND
IN APPLIED
MECHANICS
AND ENGINEERING
49 (1985) 253-279
A REDUCTION METHOD FOR NONLINEAR STRUCTURAL DYNAMIC ANALYSIS* Sergio R. IDELSOHN
and AIberto CARDONA
Mechanics Laboratory, Institute of Technological Development for the Chemical Industry, Universidad National de! Litoral, and National Council for Scientific and Technological Research of Argentina, 3MxI Santa Fe, Argentina Manuscript received 8 August 1983 Revised manuscript received 30 April 1984
A computational algorithm for predicting the nonlinear dynamic response of a structure is presented. The nonlinear system of ordinary differential equations restthing from the finite etement discretization is highly reduced by means of a Rayleigh-Ritz analysis. The basis vectors are chosen to be the current tangent eigenmodes together with some modal derivatives that indicate the way in which the spectrum is changing. Only a few basis updatings are required during the whole time integration. The truncation error introduced at every change of basis is pointed out as the cause for a divergence-type behaviour, and some means for eliminating it are discussed. Results for examples involving large displacements are shown and compared to the results obtained by integrating the complete system of equations.
1. Introduction The numerical solution of structural static and dynamic problems is currently being performed by means of the popularized finite element method. Nevertheless, the computational effort required to analyze large-scale structures quite often exceeds a reasonable amount of computer time, even on the large computers available today. This is especially the case when handling nonlinear problems. Considerable research is being devoted to reducing the computer time based on the recognized fact that the geometry of the structure very often determines the number of degrees of freedom of a structural finite element model. We can englobe these techniques as those being directed towards the improvement of the discretization before proceeding to the solution itself. Within this context, the following techniques can be mentioned: subst~cturing [l-4], block-iteration schemes [5-71, definition of master and slave nodes [S] and those techniques embodied as reduction methods [9-221. A reduction method, in particular, can be thought of as two-step discretization. The first step, is the classic finite element discretization whereas the second is based on the computation of some basis vectors in order to perform a Rayteigh-Ritz analysis. Clearly, the success of the method depends, to a great extent, on the proper selection of the basis vectors. *The work reported here received financial support from the National Council for Scientific and Technological Research of Argentina (CONICET) and from Universidad National de1 Litoral de Santa Fe, Argentina (U.N.L.). 00457825/85/$3.30
@ 1985, Elsevier Science Publishers B.V. worth-Holland)
2.54
S.R. Idelsohn, A. Cardona,
Nonlinear structural dynamic
analysis
The modal superposition technique has already been applied in nonlinear structural dynamics anaiysis [17-201 by following the actual tangent spectrum of the structure at every instant. However, the effectiveness of the employed techniques in practical analysis is questionable due to the excessive effort involved in solving the tangent eigenprobiem at every time step. Another approach has been used: a constant stiffness with a nonlinear internal-forces vector [21-221. The basis vectors are computed only once as the modes of the initial eigenspectrum of the structure. This procedure can be applied accurateIy in the analysis of only slightly nonlinear systems or systems with only local nonlinearities. The reduction method developed by Noor and co-workers, to be applied in nonlinear static structural problems [l&16], is recognized as well-suited for treating highly nonlinear behaviours. It is based on choosing as basis vectors those commonly used in the static path perturbation technique, that is to say a nonlinear solution and its various-order derivatives. Based on Noor’s idea, we propose a new choice of the basis vectors to be used in nonlinear structural dynamics, which we have called the ‘modal derivatives method’ [23]. It combines some of the main features of the above-mentioned papers. In it the modes of the actual spectrum are employed together with some modal derivatives that indicate the way in which the spectrum is changing. By doing so, a more complete knowledge of the nonlinear solution path is gained and the number of required basis updatings is highly reduced. When a basis change is performed the incompatibility between the old basis and the new one obliges the projection of the computed velocities and accelerations onto the new spanned subspace, thus resulting in a truncation error. We emphasize the importance of this truncation error as it introduces a continuously growing lack of equilibrium that finally deteriorates the solution. This spurious effect is illustrated by Example 5.1. The basis proposed in this paper, does not require such frequent updating due to its power in following the solution path, thereby considerably reducing the associated truncation error. A strategy for controlling the basis vectors updating process is proposed based on maintaining the residual norm of the complete system within a prescribed tolerance all along the time history. The basis updating has been implemented so as to diminish the incompat~biljty between both bases, thus reducing the truncation error. The reduced system of equations is integrated by means of Newmark’s implicit algorithm. Nevertheless, we think that more elaborate schemes, which incorporate some features such as automatic step-size control, can be economically used in view of the reduced size of the final system of equations. Several examples dealing with geometrically nonlinear problems are shown. No difficulty exists in extending the concept to material nonlinearity. 2. Equations of motion Let Ui denote the three components of a displacement field that satisfy ‘a priori’ the compatibility equations, related to the stress and strain fields, by the constitutive equations. The structural dynamics equilibrium equations are then: trij,j+fi = piij
in V,
O;,?Zj= ti
on r,,
(2.1)
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
25.5
where uii is the stress tensor, fi the body forces, p the density, ti the imposed surface tractions and nj the outward pointing unit normal component to the boundary surface r,. Dots are used to denote differentiation with respect to time. Damping effects were not taken into account, but the generalization for their introduction is obvious. ‘First-stage’ finite element discretization may be written in general as follows: NN
(2.2) where Lii represents the discrete approximation to the displacement field Ui, h;(xj) are the interpolation functions adopted to represent the displacement field, a,(t) characterizes the amplitude assumed in time by the interpolation functions and NN is the number of nodes of the discretization. By making use of Hamilton’s variational principle or the weighted-residual technique equation (2.1) may be expressed in a discrete manner by the N-dimensional system of differential equations: G(a)+ Mti = F,
(2.3)
where G(a) is the internal-forces vector, which is a nonlinear function of the unknown parameter a, M is the mass matrix and F is the applied-forces vector. With given initial conditions at t = to, (2.3) can be integrated to produce the time-history response of the structure. An incremental form of the discretized equilibrium equations with respect to a deformed position a0 can be obtained by defining the tangent stiffness matrix as:
Then the discretized incremental K(uo)Au + Mii = F-
and higher-order
equilibrium
G(uo),
where
equation reads: Au = a - a,
terms in Au were neglected.
3. Second-stage discretization A Rayleigh-Ritz technique is used to replace (2.3) by a reduced system of equations. This can be done by approximating Au by the linear combination of R linearly independent vectors: Au=Aa=+y,
(3.1)
256
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
where $ is the (N x R) matrix formed with R basis vectors and y is the (R unknowns. The system of differential equations to be solved now reads:
G(ao+ J/y)+ATj
=
x 1)
vector of new
F,
(3.2)
where (liy) = W@o+
&o+
$y),
ii? = $‘Mt,b,
F = t,VF.
(3.3)
The success in the performed approximation depends, to a great extent, on the proper choice of the basis vectors. These vectors must satisfy the following basic requirements [15]: (a) linear independence; (b) sufficient characterization of the nonlinear dynamic response of the structure in the neighbourhood of the considered point; (c) good approximation properties, in the sense that they can approximate the solution in a significant time interval. The first idea that appears, in order to generate a convenient basis, is to consider the incremental equilibrium equations (2.5). By assuming small increments about the position a0 we can accept that the system behavior will not differ significantly from that of a linear system. We can then conveniently express the displacement increments in terms of the homogeneous incremental linearized problem modes [17]. Then, by solving (K(ao) - w;M)#’
= 0
and assuming that the R lower frequencies response, we can express
(3.4) and their corresponding
modes
govern
the
R
Aa=C4’z,=cCly,
(3.5)
r=l
where z, are the modal participation factors which are directly identified with the generalized displacements y. The displacement increments should be small enough to maintain the validity of the linearization hypothesis. Numerical experiences show that this fact forces us to update the basis too frequently (specially when treating geometrical nonlinearities), and so the method becomes expensive. In order to decrease the number of updatings we will generate some new vectors that complete the previously used basis. By considering equation (3.4) this one gives a function dependence of the eigenmodes 4’ on the instantaneous deformed position a. The displacement vector a is changing continuously and, taking this fact into consideration, equation (3.5) may be rewritten as Aa =
5 qb’(a)~.
,=I
This equation expresses the displacements
(3.6) a as nonlinear function of the modal amplitudes
257
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
z,. If we now assume that Au can be developed position ao, we get
Aai = Aai(z,)=
[210z,+
[s],T+
into a Taylor series around the reference
(3.7)
* * * .
r
It is readily seen in (3.6) that the first term in the Taylor expansion is simply the eigenmode, and we realize that in (3.5) only these first derivatives were taken into account. Nevertheless, in order to obtain a better basis we will also include the higher-order terms. For the purpose of computing the second and higher derivatives, a first idea could be to differentiate (3.4). This procedure leads to a linear system of equations in which the unknowns are the desired derivatives; unfortunately, the coefficients matrix is singular, and special questions arise as to whether the independent term belongs to the image or not. So, instead of doing that, we propose an approximate computation as follows. We will assume that the accelerations can be expanded around the reference acceleration term ii0 in terms of the same parameters z,. We also assume that in this expansion the secondand higher-order terms can be neglected. Thus: ..
ii = iOi +
[aai az,Ioz’*
(3.8)
This proposal is in accordance with the well-known fact that the accelerations approximated to a lesser degree than the displacements. By differentiating the equilibrium equation (2.3) with respect to z,, we obtain:
can be
(3.9) This equation, computed
at instant to, gives the following expression:
= -CO@;.
(3.10)
Now, by differentiating (3.9) and noting that the second-order terms vanish, (3.8), we obtain aGi a2a. --+ &Zj
82,
J’Gi
a.2,aaj
ask
&j auk --=
a.2, a.2,
0.
derivatives
of the inertia
(3.11)
An expression which, evaluated at to, reads: (3.12)
258
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
In geometrically nonlinear problems the internal-forces vector is a cubic function of the displacements and, therefore, its derivatives of higher order than three vanish identically. Then, by consecutively differentiating equation (3.11) we can obtain recurrence relations as:
(3.13) Then, the displacement increment will be written as the linear combination number of modal displacement derivatives with respect to the parameters 2,:
of a given
(3.14) where y is the R-dimensional vector of generalized displacements in the new basis (R G N). It may be seen that this basis is made up of the tangent eigenmodes at to and a certain number of derivatives of these modes with respect to the modal amplitude parameters. The assumption that the first modes govern the response still holds. Therefore, these will exhibit a largely nonlinear behavior and so will require the computation of derivatives of higher order. This criterion helps the user to choose the number of basis vectors to be computed and how to calculate them. Finally, it should be mentioned that by slightly modifying the basis vector’s generation, (3.13), the scheme employed in [14-161 for the basis vector’s generation when solving static non-linear problems can be obtained. Then, both techniques can be easily programmed in a single computer code.
4. Evaluation of the internal-forces
vector derivatives
In order to calculate the modal derivatives, second and third derivatives of the internalforces vector have to be computed, (3.13). These are normally not evaluated in a standard finite element code and so, special routines for the element’s computations should be provided. Next we will summarize the main computations to be carried out when using a displacement-based finite element in which geometrical nonlinearities are considered. Although some advantages were reported when employing mixed models in conjunction with reduction methods (in static nonlinear problems), we have used here a displacement’s one for the sake of simplicity. The vector of strains within the element is expressed, as usual, as Ei =
Bjj
(Uk)Uj
(4.1)
5
where Ei is the strains vector, Uj is the vector of elemental
nodal displacements,
and Bij is a
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
259
matrix which relates both of them being a linear function of the displacements. The variation of strains is denoted as: S&i = Bi ( Uk)stJ
(4.2)
where the matrix B: relates variation of strains with variation of displacements. function of the displacements, too. The elemental internal-forces vector Gy is then computed as BXiuk dV,
G;=
It is a linear
(4.3)
where V, is the volume of the element and ok is the stress vector. By differentiating (4.3) the elemental tangent stiffness matrix Kfj is obtained: BtiGk dV.
(4.4)
If we take into account that s
(4.5)
= D,,B; , I
where Dkl is the Hookean matrix, and defining the three-dimensional
constant array B:, as (4.6)
then matrix K; is computed K&z
as
I
(BziDk,BE + BiijCTk)dV. VC
(4.7)
The elemental array of second derivatives of the internal forces Htk is computed by differentiating the last equation. By noting that: aB:,/aU,,, = 0
afterwards
(4.8)
we obtain (4.9) If we differentiate internal forces F;kl:
again, we can get the four-dimensional
array of third derivatives of the
S.R. Idelsohn, A. Cardona, No&near
260
structural dynamic analysis
It should be pointed out that the arrays H& and F&, are never explicitly computed. Directly, the code evaluates the contribution of each element to the right-hand side of (3.13). So, the involved computational cost is of the same order as in an evaluation of the internal-forces vector times the number of terms of the right-hand side. For instance, for the evaluation of f~3ai/&:]0 the following system of equations has to be solved: (4.11) Each term between brackets in the right-hand side is directly evaluated at the element level and then assembled.
5. Numerical time integration of the reduced equations of motion Either implicit or explicit time-integration methods may be conveniently employed for the integration of the reduced equations of motion. Since the reduction process increases the least period of the discretized model, the time step is not limited seriously by stability considerations when using an explicit scheme. But also, as we are integrating a reduced system, we are able to choose more elaborate schemes without having problems with excessive storage or CPU time consumption. In this paper, we have worked only with the Newmark’s time-integration scheme. This was implemented as follows: Step 1. Computations
associated with the basis change: ji, = J/‘Mii .
jn = @‘MA,
Step 2. Initial time-step computations: “1 Yti+i= -(l/j3 At)& + (I$+I= 1
y,+1=
1/2P)j;, ,
jn + (1- y)Atji + y Atjt+, , 0,
1 a,+~
=
a,,
pn+,
=
t+VF,+,- A?&;+, .
Step 3. Iterative cycle: G- i,,+, = JI’G;+&z:+J,
&+I = F,+j - G’,+, - (l/p At2)n;iy;+~,
If:
@t+,lj< TOL, +
end of the iterations,
If:
llRa+,ll< TOL, j
the stiffness matrix is not recomputed ,
(5.1)
S.R. Idelsohn, A. Cardona, Nonlinear structuraldynamic analysis
si = (,+&pqi, y
i+l_ n+l
-
y
i n+l
Ay’,
+
Ayi =
261
(S’)-lRi+l,
azl, = a;+,+ q!Jby’.
Step 4. End of the iterative cycle:
jn+1= j;+, + (r/P Af)yn,t,
jL+, = $+I + (l/P At2)yn+, .
If the velocities and accelerations (in+, = $j,+,,
are to be computed:
H,+* = $jL+1 *
As it may be seen in the described algorithm, the incompatibility between the old and the new basis introduces a truncation error in velocities and accelerations every time a basis change is made. This error produces a continuously growing lack of equilibrium in the complete system that affects the final results. The following example illustrates this effect. EXAMPLE 5.1. A clamped beam with a sinusoidal load on its end, was modelled by using five equally spaced beam elements (Fig. 1). The reduced system was built by using as shape functions the eight first instantaneous eigenfunctions of the system, and these were changed at every time step. The vertical displacement obtained for the end versus time is shown in Fig. 2 together with the correct answer given by integrating the complete system. Also shown (Fig. 3) is the weighted-residual norm evaluated as:
&= where
7i
IIRII
(5.2)
(pq + IIMtill) ’
R=F-G(a)-ML
(5.3)
r
I
2
A= 70.
L = 100.
I * 585.
p4.0070
E = 16000.
F 8 2500. w= 9.5
Fig. 1. Cantilever
beam. Problem description.
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
262
s x
Fig. 2. Vertical appears.
displacement domplete
at node 6. When the basis is changed at every reduced model (8 modes). system; --
time
step
I
I
a divergent
behaviour
8d
Sd
%m Fig. 3. Error measure. norm (8 modes).
1
I
0.40
The divergent
1
I
0.80
I
I
I
1.20
I
1.80
I
2.00
I
I
2.40
I
2.80
T behaviour
is clearly
pointed
out by the error
measure.
-
-
Residual
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
263
This residual norm exhibits a clearly divergence-type behavior that we attribute mainly to the truncation error committed at each time step. The same phenomenon occurs when employing some modal derivatives as shape functions and changing them at every time step. It will be shown later in the examples, that this problem can be solved accurately by using three modes and three modal derivatives, without changing the basis during the whole time period considered (the first 4sec. were solved). Nevertheless, once the system’s characteristics have changed sufficiently, a new basis needs to be computed. Therefore, some means for eliminating this error are needed. A first idea is to perform the velocities and accelerations computations over the complete system. While doing so, the displacement increment is computed as the linear combination of some pre-established basis vectors. This procedure leads to the formulation of a reduced system of nonlinear equations to be solved at every time step. Some details of such an integrator are given in Appendix A. Unfortunately, the analysis made in the appendix shows that this scheme is unstable in velocities and accelerations. As a result of this analysis, it is also concluded that any particular expression for the velocities and accelerations after performing the change of basis should be contained exactly in the new basis to avoid any exponential growth, and the truncation seems to be the most adequate form of computing them. Bearing in mind the idea of reducing the truncation error associated with the change of basis, we propose the general strategy illustrated in Fig. 4. In this scheme, the initial computed basis is retained during the whole integration. When an error measure exceeds a prescribed value, the last m steps are rejected and then, some new vectors are computed, added to the initial
t GENERATE
FIRST BASIS
+ ADVANCE
I STEP I
COMPUTE
ERROR MEASURE E I
1 REJECT
LAST
’
m STEPS
I
1
Fig. 4. Selected strategy for the basis updating process.
264
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
computed basis, and orthonormalized. The user has to choose the vectors in the initial basis, the vectors to be added, the number of steps to reject, and the error tolerance. The whole strategy is based on obtaining a computed response for which the error measure remains beneath the prescribed tolerance for the complete time history. The adopted error measure is the weighted-residual norm defined in equation (5.3).
6. Some comments about the selection of the basis vectors The strategy for the computation of a suitable basis, where the nonlinear response is to be expanded, can be seen as formed by two different phases: (1) The computation of a set of basis vectors capable of solving the linearized dynamic problem, i.e. the tangent eigenmodes. (2) The computation of a collection of vectors giving the modes of the new tangent eigenproblems that appear when the system deviates from the current position, i.e. the modal derivatives. Concerning the determination of which vectors have to be included in the basis the criterion for their selection should consider, side by side, these previously mentioned aspects. Thus, the two following criteria have to be verified: (1) Convergence from the point of view of a dynamic problem. (2) Convergence from the point of view of a nonlinear problem. The convergence of the dynamic problem can be assured by following similar criteria to those commonly employed for linear dynamic problems. For this reason, a recommended practice is to verify that the projection of the load onto the orthogonal subspace to the span of the basis does not contain any important component [24]. After having defined which modes are to be included, we should analyze which modal derivatives shall be computed and added to the basis. The most desirable formulation would be one in which derivatives of every selected mode are included. However, this criterion could lead to an extremely large number of degrees of freedom in the reduced system. It is not difficult to recognize that only those modes with large enough participation factors will develop large values in the associated derivative terms of the Taylor’s expansion. Thus, we can select derivatives of the modes of higher participation, i.e. the modes of smallest associated eigenfrequency. The enunciated general criteria serve only as a guide for the selection of the basis vectors. Once a particular basis has been chosen, the system is integrated by following the standard Rayleigh-Ritz procedure. However, if the basis has been selected badly, we will not be able to maintain the error norm below an acceptable tolerance, although the basis is recomputed too frequently. Thus, we can see that the process itself indicates when additional vectors must be included in the basis in order to achieve a good approximation to the true response.
7. Numerical examples Three geometrically nonlinear examples were tested: a clamped beam and two arches. All of them were modelled by using a 8-DOF (6 connected, 2 internal) straight-beam element with
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
265
cubic interpolation for both the axial and transversal displacement fields f25]. The reduced solutions were compared with the solution obtained by integrating the complete system of nonlinear equations with the Newmark’s time integrator. Both integrations were performed using the same time step. It is not our intention to discuss whether the results are physically meaningful, but rather to compare both numerical solutions, taking as exact the complete system solution. 7.1. Clamped beam with a sinusoidal varying load In Example 5.1, a clamped beam was solved by using a basis formed by the eight first tangent modes. This example showed that although these modes were updated at each time step, the system’s nonlinear response could not be predicted accurately due to the divergent behaviour associated with the change of basis. Now, the same beam was solved by using the six following basis vectors:
that is to say, the three first modes, the first derivative of the first mode with respect to its own amplitude, the first derivative of the first mode with respect to the second mode’s amplitude and the first derivative of the second mode with respect to its own amplitude. This basis was computed only once, at the initial time instant, and no changes of basis were required during the whole integration. The nonlinear system was integrated up to t = 4 using a time step of 0.025 (160 steps). Figs. 5 to 8 show the computed vertical and horizontal displacement in time of nodes 6 and 3 compared with the exact solution obtained by integrating the complete system of equations. We can see that, simply by including three modal derivatives in the basis, an almost complete agreement with the full system solution was obtained without recomputing the basis, Previously, when eight modes were taken as generalized degrees of freedom, the results diverged even though these modes were recomputed at each time step (Fig. 2). Fig. 9 displays the evolution in time of the error measure, its maximum value being under 5.7 per cent. When eight modes were employed, the error measure reached peaks of 120 per cent, as indicated in Fig. 3. We should remark that not only the vertical but also the horizontal displacement components which have a relatively high frequency content have been approximated correctly. 7.2. Clamped arch with a suddenly appi~ed load One half of the clamped arch shown in Fig. 10 was modelled by 14 equally spaced beam elements leading to 73 degrees of freedom. The example in Section 7.1 displayed a case for which the inclusion of the modal derivatives in the first computed basis allowed to predict the whole nonlinear response. The present example is intended to show that, for most structures, the basis also has to be updated periodically to conform the system’s local behaviour. We also show that the time instant in which the basis needs to be updated is clearly indicated by the error measure. In order to see the need of performing the basis updating, the arch was first solved by
266
8,
‘0.00
t
I
I
I
,
1.20
0.60
I
I
I .80
t
i
2.40
3.00
i
I
3.60
I
4.20
T
Fig. 5. Cantilever beam. Vertical displacement at node 6. The three first modes together with three modal derivatives were employed in the basis. The same set was used throughout the response analysis. ~ reduced model (3 modes, 3 derivatives, no basis change). Complete system; -8 *;
8
6
0.00
’
1
0.60
Fig. 6. Horizontal displacement derivatives, no basis change).
I
I
1.20
at node 6.
I
I
1.80
I
2.40
I
T Complete
I
1
3.00
system; --
I
3.0
I
,
I
4.20
reduced
model (3 modes, 3
267
gl] 0.00
0.60
1.80
1.20
2.40
3.00
s-60
4.20
T Fig. 7. Vertical displacement derivatives, no basis change).
I
Complete
at node 3.
I 0.00
1
I
1.20
1
1
1.80
,
I
2.40
system; --
1
I
3.00
reduced
1
I
3.60
I
model (3 modes, 3
1 4.20
T Fig. 8. Horizontal displacement derivatives, no basis change).
at node 3.
Complete
system; --
reduced
model (3 modes, 3
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
268
“I
00.00
I
I
I
0.60
1
1.20
I
I
I
I.80
I
I
2.40
I
3.00
1
1
I
I
3.60
4.20
T
Fig. 9. The error measure remained derivatives, no basis change).
t
below 5.7% during
the whole
E. 6.895
2
I = 1.0925
I
analysis.
x IO”
h= 5.08
x lO-3
Fig. 10. Clamped
arch. Problem
Residual
x IO- 2
8 =0.245 fib 2.607x
R * 2.54 PO’ 4.7034
--
x IO’
description.
rod IO3
norm (3 modes,
3
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
269
leaving the first computed basis intentionally fixed. That is, no tolerance was set to the error measure and the basis remained unchanged although the error measure reached high values. The first mode and its first derivative were used as basis vectors and a total of 50 time steps of 0.0002 were computed. Fig. 11 presents the response obtained at the apex vertical displacement compared with the complete solution. As it was expected, reduced solution does not agree with the complete one, and the error measure indicates that the results were not being computed correctly (Fig. 12). The same structure was again solved by following the next procedure (Fig. 4): each time the error measure exceeded 0.05, the last two steps computed were rejected and the first tangent mode and its first derivative, were calculated and added to the initial basis made of the first mode and its first derivative. Only two changes of basis were made at times 0.0016 and 0.0034 and 50 time steps of 0.0002 were solved. The vertical displacement of nodes 15 and 8 and the horizontal displacement of node 8 are plotted together with the complete system responses in Figs. 13 to 15. Fig. 16 shows the evolution in time of the residual norm. We can appreciate that now the 4-DOF reduced solution and the complete system solution have a close agreement. 7.3. Clamped arch with an unsymmetrical suddenly applied load Here the same arch as in Section 7.2 was loaded with an unsymmetrical load. The entire was now modelled with 14 elements, leading to 73 degrees of freedom (Fig. 17).
Fig. 11. Clamped arch. Vertical displacement at node 15. The basis formed by the first mode and its derivative was left fixed intentionally. Thus, the results were highly inaccurate. Complete system; -reduced model (1 mode, 1 derivative, no basis change).
Fig. 12. The error measure indicated clearly that the basis should have been recalculated. (1 mode, 1 derivative, no basis change).
-
-
Residual norm
t 3
8 d ::
9 L*
Cn?
;;P 0
4 t 9 c
90.00
0.16
0.32
0*46
0.64
T
0.60
0.96
1.12
*10-2
Fig. 13. Clamped arch. Vertical displacement at node 15. The first tangent mode and its derivative were added to reduced model (1 Complete system; -the first step computed basis at times 0.0016 and 0.0034. mode, 1 derivative + 1 mode, 1 derivative). 270
?
P
Fig.
14. Vertical
derivative
displacement
+ 1 mode,
at node
8. -
Complete
system;
--
reduced
model
(1 mode,
1
1 derivative).
Fig. 15. Horizontal displacement derivative + 1 mode, 1 derivative).
at node
8. -
Complete
271
system;
--
reduced
model
(1 mode,
1
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
272
8 d
8d
Sd ‘;
0
?\ f\ I \I \ I \
A II II
-
a-/
/1
I
I
\
I_ir;pi
f
I
&~I, / i I-,I”1 I 1 ,i’” I (1
\/,I
E-1 4,
,
I
Oo.00 0.16
I
I
I
0.32
I
I
0.46 1
Fig. 16, When the error measure -Residual norm (1 mode,
I
0.64
I
llo-2
I
0.60
I
I
0.96
exceeded 0.05, a basis updating was performed 1 derivative + 1 mode, 1 derivative).
I
(at t = 0.0016 and 0.0034).
I E=6.895
x IOK)
h=S.O8xIO+
I = 1.0925 x 10-s
8 ~0.245
R=2.54
fJ -2.607x
rad IO’
PO’ lx107
Fig. 17. Unsymmetrically
loaded
clamped
arch. Problem
I
1.12
description.
273
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
Various reduced solutions were computed and compared with the exact complete one. The example was intended to show how the response is more closely approximated, when the number of generalized degrees of freedom is increased. All integrations were performed up to t = 0.1, with a time step of 0.0002 (50 steps). First of all, a reduced model made up by three modes and the first three derivatives of the two first modes was solved without performing, intentionally, any basis updating during the whole integration. Figs. 18 and 19 show the computed responses for the vertical component of displacements at nodes 5 and 8. We clearly see that the reduced solution is badly computed, and that the error norm (Fig. 20) indicates the need of changing the basis. Then, an error tolerance of 0.045 was set. The first mode and its first derivative were added to the previously mentioned basis, and the last three time steps were rejected, every time a basis change was performed. A total of two changes of basis were made at times 0.0016 and 0.0024. An evident improvement in the results can be noted (Figs. 18 and 19). At last, the same error tolerance of 0.045 was imposed while the first two modes and their first three derivatives were added and the last three steps were rejected at each basis change. Only one change of basis was required at time 0.0016. A high degree of accuracy was obtained evidencing a clear asymptotic behaviour, when the number of degrees of freedom in the reduced model is increased. Finally, we should comment that the horizontal component of displacements, which in this case has an extremely high frequency content, was estimated within an acceptable margin of
tit .
.
0.32
0.42
T
0.64
llo-2
0.20
0.26
1.12
Fig. 18. Unsymmetrically loaded clamped arch. Vertical displacement at node 8. Complete system; -reduced model (3 modes, 3 derivatives, no basis change); ----reduced model (3 modes, 3 derivatives + 1 mode, 1 derivative); --reduced model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
8 d
1
$1 I 0.00
1
0.16I
I
0.32I
1
0.48 I
I
1
0.64 1
1
I
I
0.80
1
0.99
1
1
1.12
a10-2
Fig. 19. Vertical displacement at node 5. ~ Complete system; derivatives, no basis change); ----reduced model (3 modes, 3 derivatives model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
-+ 1 mode,
reduced model (3 modes, 3 1 derivative); --reduced
(Y c;
0 d
? I I 00.000.16
I
0.22
I
1 0.481 I 1
0.94
I
llo-2
1 0.80r I
0.99
I
I
1
1.12
Fig. 20. Error measure. -Residual norm (3 modes, 3 derivatives, no basis change); ----residual norm (3 residual norm (3 modes, 3 derivatives + 2 modes, 3 derivatives). modes, 3 derivatives + 1 mode, 1 derivative); ---
274
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
275
error (Figs. 21 and 22). It must be remarked that these displacements have an amplitude several orders of magnitude below the vertical displacements amplitude, thus their influence in predicting the overall structure response is minimal. The same behavior may be found in similar linear problems solved by the standard mode superposition method.
8. Concluding
remarks
A computational algorithm for predicting the nonlinear dynamic response of a structure by means of a reduction scheme has been described. In it, the nonlinear system of ordinary differential equations obtained from the finite element discretization, is reduced by employing a Rayleigh-Ritz analysis. Some new aspects have been discussed, i.e. the employment of the modal derivatives as basis vectors, the truncation error introduced at each change of basis, a specially oriented implementation of the Newmark’s time-integration scheme, and a reductionerror measure. The use of the tangent modes as basis vectors and the need for performing a periodic basis updating in nonlinear dynamics have been pointed out by several authors. However, the cost of such method seriously questions its applicability, specially when treating severe nonlinearities, due to the large number of required basis updatings. Besides, in this paper, we have
R d
a d
h
c: 0
d
‘g
d i% ;; 6 d 8 + n
40.00
0.16
0.32
0.40
0.64
1
llo-2
0.00
0.m
1.12
Fig. 21. Horizontal displacement at node 8. Complete system; -reduced model (3 modes, 3 derivatives, no basis change); ----reduced model (3 modes, 3 derivatives + 1 mode, 1 derivative); --reduced model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
276
S.R. Idelsohn, A. Cardona, Nonlineur structural dynamic analysis
Fig. 22. Horizontal displacement at node 5. reduced model (3 modes, 3 Complete system; -derivatives, no basis change); ----reduced model (3 modes, 3 derivatives f 1 mode, 1 derivative); --reduced model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
shown the important role played by the truncation error introduced at each change of basis. We have pointed out that it produces a divergence in the solution when too many updatings are performed. We think that this error can be more drastic in those schemes where not only velocities and accelerations but also displacements are re-expanded at each change of basis. We are now proposing to employ some modal derivatives in the basis together with the tangent modes. These derivatives give to the basis a valuable knowledge on how the system nonlinearity is being developed. By employing them, the number of required basis updatings is highly reduced and the nonlinear dynamic solution path can be followed using only a few basis vectors as generalized degrees of freedom. The proposed method of time integration and basis updating was developed, bearing in mind the purpose of preventing the separation of the reduced solution from the true solution path at any time. Their separation is efficiently detected by the error measure and whenever it occurs the last computed steps are rejected so that the solution path is reached again. Then, the basis is changed and the time integration is continued. The additional computational effort the method requires is negligible compared with the computer time gained during the reduced systems time integration. This is particularly true due to the limited number of basis updatings needed in the proposed formulation. Results concerning the application of this technique to some nonlinear undamped dynamic problems were presented, which showed the power of the method.
S.R. Idelsohn, A. Cardona, Nonlinear structuraf dynamic analysis
277
Appendix A An alternative integration scheme to that presented can be thought of as follows:
in (5.1) in which no truncation
is made
Step 1. Initial time-step computations: 1
art+1 = 48 t 1 = - (l/p At)& + (1 - l/2/3)&, a*aPI+1
‘1 a,+1 = h, +
1 yn+1=0,
(1- y)At& + y Atii:,,
,
(A4
F,,, = Cft’F,,,, - $‘iWi:+~.
Step 2. Iterative cycle: The same as in the scheme (5.1). Step 3. End of the iterative cycle:
This integrator would lead to a continuously growing error in velocities and accelerations as it can be seen from the following analysis. The velocities and accelerations at time t,, can be written as hi, = (1 - #)Li,-I+ &, =
(-l/p
At{1 - 7/2/3)ii,-, + (r/p At)Aa, ,
At)&
64.2)
+ (1- 1/2@&_1+ (l/p At2)Aa,.
By applying the iatter repeatedly
and writing it in matrix form:
+ 5 A’-‘L Aai, It] =Af;:l i=l
(A-3)
where
and A” denotes the nth power of A. Suppose that at instant to a basis change was performed. read:
So, the velocities and accelerations
(A-4)
278
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
The error terms appear because the velocities and accelerations can not be exactly re-expanded in the new basis. By replacing in equation (A.3) we obtain:
44 e (4
].
(A.9
The last term introduces an error that grows exponentially with II. It cannot be balanced out by the first terms because they lie in a subspace orthogonal to the space spanned by e(&) and e(&). So, it is necessary to make the corrections on velocities and accelerations indicated in Step 1 in (5.1).
References [II R.W. Clough and E.L. Wilson, Dynamic
analysis of large structural systems with local nonlinearities, Comput. Meths. Appl. Mech. Engrg. 17/18 (1979) 107-129. in linear and nonlinear analysis, Internat. J. Numer. Meths. PI R.H. Dodds and L.A. Lopez, Substructuring Engrg. 15 (1980) 583-597. techniques in material nonlinear analysis, Internat. J. [31 D.R.J. Owen and O.I.A. Goncalvez, Substructuring Comput. & Structures 15 (1982) 205-213. [41 L. Meirovitch and A.L. Hale, On the substructure synthesis method, AIAA J. 19 (1981) 940-947. in: K.J. Bathe, J.T. Oden and W. Wunderlich, eds., 151 E.L. Wachspress, Two-level finite element computation, Formulations and Computational Algorithms in Finite Element Analysis (MIT, Cambridge, MA, 1977) 877-913. scheme for the solution of systems of equations resulting bl L.C. Wellford Jr. and B. Vahdani, A block-iteration from linear and nonlinear finite element models, Comput. Meths. Appl. Mech. Engrg. 26 (1981) 33-52. An iteration procedure for reducing the expenses of static, elastoplastic and eigenvalue [71 W. Dirschmid, problems in finite element analysis, Comput. Meths. Appt. Mech. Engrg. 35 (1982) 15-33. in structural analysis, Internat. J. Numer. Meths. PI Y.T. Leung, An accurate method of dynamic condensation Engrg. 12 (1978) 1705-1715. 191 R. Guyan, Reduction of stiffness and mass matrices, AIAA J. 3 (1965) 380. [lOI B.M. Irons, Structural eigenvalue problem: elimination of unwanted variables, AIAA J. 3 (1965) 961-962. of Ritz vectors, [Ill E.L. Wilson, M.W. Yuan and Y.M. Dickens, Dynamic analysis by direct superposition Earthquake Engrg. Struct. Dyn. 10 (1982) 813-821. nonlinear finite element behavior using buckling mode superposition, [121 D.A. Nagy and M. Konig, Geometrically Comput. Meths. Appl. Mech. Engrg. 19 (1979) 447-484. [13] E.G. Carnoy, Asymptotic study of the elastic post-buckling behavior of structures by the finite element method, Comput. Meths. Appl. Mech. Engrg. 29 (1981) 147-173. [14] A.K. Noor and J.M. Peters, Reduced basis technique for nonlinear analysis of structures, AIAA J. 18 (1980) 455-462. [15] A.K. Noor, Recent advances in reduction methods for nonlinear problems, Internat. J. Comput. & Structures 13 (1981) 314. [16] A.K. Noor, On making large nonlinear problems small, Comput. Meths. Appl. Mech. Engrg. 34 (1982) 955-985. [17] R.E. Nickell, Nonlinear dynamics by mode superposition, Comput. Meths. Appl. Mech. Engrg. 7 (1976) 107-129. [18] N.F. Morris, The use of modal superposition in nonlinear dynamics, Internat. J. Comput. & Structures 7 (1977) 65-72. [19] S.N. Remseth, Nonlinear static and dynamic analysis of framed structures, Internat. J. Comput. & Structures 10 (1979) 879-897.
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
279
[20] L. Landau, N.F. Favilla Ebecken and E.C. Prates de Lima, Comportamiento dinamico nHo-linear de estruturas pelo metodo de superposicao modal, Proc. I Cong. Lat. Amer. Met. Comp. Engng., Porto Alegre, Brasil (1979) 129-151. [21] V.N. Shah, G.I. Bohm and A.N. Nahavandi, Modal superposition method for computationally economical nonlinear structural analysis, J. Press. Vessel Tech. 101 (1979) 134-141. [22] K.J. Bathe and S. Gracewski, On nonlinear dynamic analysis using substructuring and mode superposition, Internat. J. Comput. & Structures 13 (1981) 699-707. [23] A. Cardona, S. Idelsohn and G. Sander, Two stage finite element discretization in structural dynamics, Proc. Third Internat. Symp. Numerical Methods Engineering, March 1983, Paris (1983) 437-445. [24] A. Huck, Application de l’analyse modale par elements finis au calcul de la reponse transitoire et forcee des structures elastiques, Rapport VF-14, Laboratoire de techniques aeronautiques et spatiales, Universite de Liege, Octobre 1972. [25] H.R. Mimer, Accurate finite element analysis of large displacements in skeletal frames, Internat. J. Comput. & Structures 14 (1981) 205-210.
IN APPLIED
MECHANICS
AND ENGINEERING
49 (1985) 253-279
A REDUCTION METHOD FOR NONLINEAR STRUCTURAL DYNAMIC ANALYSIS* Sergio R. IDELSOHN
and AIberto CARDONA
Mechanics Laboratory, Institute of Technological Development for the Chemical Industry, Universidad National de! Litoral, and National Council for Scientific and Technological Research of Argentina, 3MxI Santa Fe, Argentina Manuscript received 8 August 1983 Revised manuscript received 30 April 1984
A computational algorithm for predicting the nonlinear dynamic response of a structure is presented. The nonlinear system of ordinary differential equations restthing from the finite etement discretization is highly reduced by means of a Rayleigh-Ritz analysis. The basis vectors are chosen to be the current tangent eigenmodes together with some modal derivatives that indicate the way in which the spectrum is changing. Only a few basis updatings are required during the whole time integration. The truncation error introduced at every change of basis is pointed out as the cause for a divergence-type behaviour, and some means for eliminating it are discussed. Results for examples involving large displacements are shown and compared to the results obtained by integrating the complete system of equations.
1. Introduction The numerical solution of structural static and dynamic problems is currently being performed by means of the popularized finite element method. Nevertheless, the computational effort required to analyze large-scale structures quite often exceeds a reasonable amount of computer time, even on the large computers available today. This is especially the case when handling nonlinear problems. Considerable research is being devoted to reducing the computer time based on the recognized fact that the geometry of the structure very often determines the number of degrees of freedom of a structural finite element model. We can englobe these techniques as those being directed towards the improvement of the discretization before proceeding to the solution itself. Within this context, the following techniques can be mentioned: subst~cturing [l-4], block-iteration schemes [5-71, definition of master and slave nodes [S] and those techniques embodied as reduction methods [9-221. A reduction method, in particular, can be thought of as two-step discretization. The first step, is the classic finite element discretization whereas the second is based on the computation of some basis vectors in order to perform a Rayteigh-Ritz analysis. Clearly, the success of the method depends, to a great extent, on the proper selection of the basis vectors. *The work reported here received financial support from the National Council for Scientific and Technological Research of Argentina (CONICET) and from Universidad National de1 Litoral de Santa Fe, Argentina (U.N.L.). 00457825/85/$3.30
@ 1985, Elsevier Science Publishers B.V. worth-Holland)
2.54
S.R. Idelsohn, A. Cardona,
Nonlinear structural dynamic
analysis
The modal superposition technique has already been applied in nonlinear structural dynamics anaiysis [17-201 by following the actual tangent spectrum of the structure at every instant. However, the effectiveness of the employed techniques in practical analysis is questionable due to the excessive effort involved in solving the tangent eigenprobiem at every time step. Another approach has been used: a constant stiffness with a nonlinear internal-forces vector [21-221. The basis vectors are computed only once as the modes of the initial eigenspectrum of the structure. This procedure can be applied accurateIy in the analysis of only slightly nonlinear systems or systems with only local nonlinearities. The reduction method developed by Noor and co-workers, to be applied in nonlinear static structural problems [l&16], is recognized as well-suited for treating highly nonlinear behaviours. It is based on choosing as basis vectors those commonly used in the static path perturbation technique, that is to say a nonlinear solution and its various-order derivatives. Based on Noor’s idea, we propose a new choice of the basis vectors to be used in nonlinear structural dynamics, which we have called the ‘modal derivatives method’ [23]. It combines some of the main features of the above-mentioned papers. In it the modes of the actual spectrum are employed together with some modal derivatives that indicate the way in which the spectrum is changing. By doing so, a more complete knowledge of the nonlinear solution path is gained and the number of required basis updatings is highly reduced. When a basis change is performed the incompatibility between the old basis and the new one obliges the projection of the computed velocities and accelerations onto the new spanned subspace, thus resulting in a truncation error. We emphasize the importance of this truncation error as it introduces a continuously growing lack of equilibrium that finally deteriorates the solution. This spurious effect is illustrated by Example 5.1. The basis proposed in this paper, does not require such frequent updating due to its power in following the solution path, thereby considerably reducing the associated truncation error. A strategy for controlling the basis vectors updating process is proposed based on maintaining the residual norm of the complete system within a prescribed tolerance all along the time history. The basis updating has been implemented so as to diminish the incompat~biljty between both bases, thus reducing the truncation error. The reduced system of equations is integrated by means of Newmark’s implicit algorithm. Nevertheless, we think that more elaborate schemes, which incorporate some features such as automatic step-size control, can be economically used in view of the reduced size of the final system of equations. Several examples dealing with geometrically nonlinear problems are shown. No difficulty exists in extending the concept to material nonlinearity. 2. Equations of motion Let Ui denote the three components of a displacement field that satisfy ‘a priori’ the compatibility equations, related to the stress and strain fields, by the constitutive equations. The structural dynamics equilibrium equations are then: trij,j+fi = piij
in V,
O;,?Zj= ti
on r,,
(2.1)
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
25.5
where uii is the stress tensor, fi the body forces, p the density, ti the imposed surface tractions and nj the outward pointing unit normal component to the boundary surface r,. Dots are used to denote differentiation with respect to time. Damping effects were not taken into account, but the generalization for their introduction is obvious. ‘First-stage’ finite element discretization may be written in general as follows: NN
(2.2) where Lii represents the discrete approximation to the displacement field Ui, h;(xj) are the interpolation functions adopted to represent the displacement field, a,(t) characterizes the amplitude assumed in time by the interpolation functions and NN is the number of nodes of the discretization. By making use of Hamilton’s variational principle or the weighted-residual technique equation (2.1) may be expressed in a discrete manner by the N-dimensional system of differential equations: G(a)+ Mti = F,
(2.3)
where G(a) is the internal-forces vector, which is a nonlinear function of the unknown parameter a, M is the mass matrix and F is the applied-forces vector. With given initial conditions at t = to, (2.3) can be integrated to produce the time-history response of the structure. An incremental form of the discretized equilibrium equations with respect to a deformed position a0 can be obtained by defining the tangent stiffness matrix as:
Then the discretized incremental K(uo)Au + Mii = F-
and higher-order
equilibrium
G(uo),
where
equation reads: Au = a - a,
terms in Au were neglected.
3. Second-stage discretization A Rayleigh-Ritz technique is used to replace (2.3) by a reduced system of equations. This can be done by approximating Au by the linear combination of R linearly independent vectors: Au=Aa=+y,
(3.1)
256
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
where $ is the (N x R) matrix formed with R basis vectors and y is the (R unknowns. The system of differential equations to be solved now reads:
G(ao+ J/y)+ATj
=
x 1)
vector of new
F,
(3.2)
where (liy) = W@o+
&o+
$y),
ii? = $‘Mt,b,
F = t,VF.
(3.3)
The success in the performed approximation depends, to a great extent, on the proper choice of the basis vectors. These vectors must satisfy the following basic requirements [15]: (a) linear independence; (b) sufficient characterization of the nonlinear dynamic response of the structure in the neighbourhood of the considered point; (c) good approximation properties, in the sense that they can approximate the solution in a significant time interval. The first idea that appears, in order to generate a convenient basis, is to consider the incremental equilibrium equations (2.5). By assuming small increments about the position a0 we can accept that the system behavior will not differ significantly from that of a linear system. We can then conveniently express the displacement increments in terms of the homogeneous incremental linearized problem modes [17]. Then, by solving (K(ao) - w;M)#’
= 0
and assuming that the R lower frequencies response, we can express
(3.4) and their corresponding
modes
govern
the
R
Aa=C4’z,=cCly,
(3.5)
r=l
where z, are the modal participation factors which are directly identified with the generalized displacements y. The displacement increments should be small enough to maintain the validity of the linearization hypothesis. Numerical experiences show that this fact forces us to update the basis too frequently (specially when treating geometrical nonlinearities), and so the method becomes expensive. In order to decrease the number of updatings we will generate some new vectors that complete the previously used basis. By considering equation (3.4) this one gives a function dependence of the eigenmodes 4’ on the instantaneous deformed position a. The displacement vector a is changing continuously and, taking this fact into consideration, equation (3.5) may be rewritten as Aa =
5 qb’(a)~.
,=I
This equation expresses the displacements
(3.6) a as nonlinear function of the modal amplitudes
257
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
z,. If we now assume that Au can be developed position ao, we get
Aai = Aai(z,)=
[210z,+
[s],T+
into a Taylor series around the reference
(3.7)
* * * .
r
It is readily seen in (3.6) that the first term in the Taylor expansion is simply the eigenmode, and we realize that in (3.5) only these first derivatives were taken into account. Nevertheless, in order to obtain a better basis we will also include the higher-order terms. For the purpose of computing the second and higher derivatives, a first idea could be to differentiate (3.4). This procedure leads to a linear system of equations in which the unknowns are the desired derivatives; unfortunately, the coefficients matrix is singular, and special questions arise as to whether the independent term belongs to the image or not. So, instead of doing that, we propose an approximate computation as follows. We will assume that the accelerations can be expanded around the reference acceleration term ii0 in terms of the same parameters z,. We also assume that in this expansion the secondand higher-order terms can be neglected. Thus: ..
ii = iOi +
[aai az,Ioz’*
(3.8)
This proposal is in accordance with the well-known fact that the accelerations approximated to a lesser degree than the displacements. By differentiating the equilibrium equation (2.3) with respect to z,, we obtain:
can be
(3.9) This equation, computed
at instant to, gives the following expression:
= -CO@;.
(3.10)
Now, by differentiating (3.9) and noting that the second-order terms vanish, (3.8), we obtain aGi a2a. --+ &Zj
82,
J’Gi
a.2,aaj
ask
&j auk --=
a.2, a.2,
0.
derivatives
of the inertia
(3.11)
An expression which, evaluated at to, reads: (3.12)
258
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
In geometrically nonlinear problems the internal-forces vector is a cubic function of the displacements and, therefore, its derivatives of higher order than three vanish identically. Then, by consecutively differentiating equation (3.11) we can obtain recurrence relations as:
(3.13) Then, the displacement increment will be written as the linear combination number of modal displacement derivatives with respect to the parameters 2,:
of a given
(3.14) where y is the R-dimensional vector of generalized displacements in the new basis (R G N). It may be seen that this basis is made up of the tangent eigenmodes at to and a certain number of derivatives of these modes with respect to the modal amplitude parameters. The assumption that the first modes govern the response still holds. Therefore, these will exhibit a largely nonlinear behavior and so will require the computation of derivatives of higher order. This criterion helps the user to choose the number of basis vectors to be computed and how to calculate them. Finally, it should be mentioned that by slightly modifying the basis vector’s generation, (3.13), the scheme employed in [14-161 for the basis vector’s generation when solving static non-linear problems can be obtained. Then, both techniques can be easily programmed in a single computer code.
4. Evaluation of the internal-forces
vector derivatives
In order to calculate the modal derivatives, second and third derivatives of the internalforces vector have to be computed, (3.13). These are normally not evaluated in a standard finite element code and so, special routines for the element’s computations should be provided. Next we will summarize the main computations to be carried out when using a displacement-based finite element in which geometrical nonlinearities are considered. Although some advantages were reported when employing mixed models in conjunction with reduction methods (in static nonlinear problems), we have used here a displacement’s one for the sake of simplicity. The vector of strains within the element is expressed, as usual, as Ei =
Bjj
(Uk)Uj
(4.1)
5
where Ei is the strains vector, Uj is the vector of elemental
nodal displacements,
and Bij is a
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
259
matrix which relates both of them being a linear function of the displacements. The variation of strains is denoted as: S&i = Bi ( Uk)stJ
(4.2)
where the matrix B: relates variation of strains with variation of displacements. function of the displacements, too. The elemental internal-forces vector Gy is then computed as BXiuk dV,
G;=
It is a linear
(4.3)
where V, is the volume of the element and ok is the stress vector. By differentiating (4.3) the elemental tangent stiffness matrix Kfj is obtained: BtiGk dV.
(4.4)
If we take into account that s
(4.5)
= D,,B; , I
where Dkl is the Hookean matrix, and defining the three-dimensional
constant array B:, as (4.6)
then matrix K; is computed K&z
as
I
(BziDk,BE + BiijCTk)dV. VC
(4.7)
The elemental array of second derivatives of the internal forces Htk is computed by differentiating the last equation. By noting that: aB:,/aU,,, = 0
afterwards
(4.8)
we obtain (4.9) If we differentiate internal forces F;kl:
again, we can get the four-dimensional
array of third derivatives of the
S.R. Idelsohn, A. Cardona, No&near
260
structural dynamic analysis
It should be pointed out that the arrays H& and F&, are never explicitly computed. Directly, the code evaluates the contribution of each element to the right-hand side of (3.13). So, the involved computational cost is of the same order as in an evaluation of the internal-forces vector times the number of terms of the right-hand side. For instance, for the evaluation of f~3ai/&:]0 the following system of equations has to be solved: (4.11) Each term between brackets in the right-hand side is directly evaluated at the element level and then assembled.
5. Numerical time integration of the reduced equations of motion Either implicit or explicit time-integration methods may be conveniently employed for the integration of the reduced equations of motion. Since the reduction process increases the least period of the discretized model, the time step is not limited seriously by stability considerations when using an explicit scheme. But also, as we are integrating a reduced system, we are able to choose more elaborate schemes without having problems with excessive storage or CPU time consumption. In this paper, we have worked only with the Newmark’s time-integration scheme. This was implemented as follows: Step 1. Computations
associated with the basis change: ji, = J/‘Mii .
jn = @‘MA,
Step 2. Initial time-step computations: “1 Yti+i= -(l/j3 At)& + (I$+I= 1
y,+1=
1/2P)j;, ,
jn + (1- y)Atji + y Atjt+, , 0,
1 a,+~
=
a,,
pn+,
=
t+VF,+,- A?&;+, .
Step 3. Iterative cycle: G- i,,+, = JI’G;+&z:+J,
&+I = F,+j - G’,+, - (l/p At2)n;iy;+~,
If:
@t+,lj< TOL, +
end of the iterations,
If:
llRa+,ll< TOL, j
the stiffness matrix is not recomputed ,
(5.1)
S.R. Idelsohn, A. Cardona, Nonlinear structuraldynamic analysis
si = (,+&pqi, y
i+l_ n+l
-
y
i n+l
Ay’,
+
Ayi =
261
(S’)-lRi+l,
azl, = a;+,+ q!Jby’.
Step 4. End of the iterative cycle:
jn+1= j;+, + (r/P Af)yn,t,
jL+, = $+I + (l/P At2)yn+, .
If the velocities and accelerations (in+, = $j,+,,
are to be computed:
H,+* = $jL+1 *
As it may be seen in the described algorithm, the incompatibility between the old and the new basis introduces a truncation error in velocities and accelerations every time a basis change is made. This error produces a continuously growing lack of equilibrium in the complete system that affects the final results. The following example illustrates this effect. EXAMPLE 5.1. A clamped beam with a sinusoidal load on its end, was modelled by using five equally spaced beam elements (Fig. 1). The reduced system was built by using as shape functions the eight first instantaneous eigenfunctions of the system, and these were changed at every time step. The vertical displacement obtained for the end versus time is shown in Fig. 2 together with the correct answer given by integrating the complete system. Also shown (Fig. 3) is the weighted-residual norm evaluated as:
&= where
7i
IIRII
(5.2)
(pq + IIMtill) ’
R=F-G(a)-ML
(5.3)
r
I
2
A= 70.
L = 100.
I * 585.
p4.0070
E = 16000.
F 8 2500. w= 9.5
Fig. 1. Cantilever
beam. Problem description.
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
262
s x
Fig. 2. Vertical appears.
displacement domplete
at node 6. When the basis is changed at every reduced model (8 modes). system; --
time
step
I
I
a divergent
behaviour
8d
Sd
%m Fig. 3. Error measure. norm (8 modes).
1
I
0.40
The divergent
1
I
0.80
I
I
I
1.20
I
1.80
I
2.00
I
I
2.40
I
2.80
T behaviour
is clearly
pointed
out by the error
measure.
-
-
Residual
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
263
This residual norm exhibits a clearly divergence-type behavior that we attribute mainly to the truncation error committed at each time step. The same phenomenon occurs when employing some modal derivatives as shape functions and changing them at every time step. It will be shown later in the examples, that this problem can be solved accurately by using three modes and three modal derivatives, without changing the basis during the whole time period considered (the first 4sec. were solved). Nevertheless, once the system’s characteristics have changed sufficiently, a new basis needs to be computed. Therefore, some means for eliminating this error are needed. A first idea is to perform the velocities and accelerations computations over the complete system. While doing so, the displacement increment is computed as the linear combination of some pre-established basis vectors. This procedure leads to the formulation of a reduced system of nonlinear equations to be solved at every time step. Some details of such an integrator are given in Appendix A. Unfortunately, the analysis made in the appendix shows that this scheme is unstable in velocities and accelerations. As a result of this analysis, it is also concluded that any particular expression for the velocities and accelerations after performing the change of basis should be contained exactly in the new basis to avoid any exponential growth, and the truncation seems to be the most adequate form of computing them. Bearing in mind the idea of reducing the truncation error associated with the change of basis, we propose the general strategy illustrated in Fig. 4. In this scheme, the initial computed basis is retained during the whole integration. When an error measure exceeds a prescribed value, the last m steps are rejected and then, some new vectors are computed, added to the initial
t GENERATE
FIRST BASIS
+ ADVANCE
I STEP I
COMPUTE
ERROR MEASURE E I
1 REJECT
LAST
’
m STEPS
I
1
Fig. 4. Selected strategy for the basis updating process.
264
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
computed basis, and orthonormalized. The user has to choose the vectors in the initial basis, the vectors to be added, the number of steps to reject, and the error tolerance. The whole strategy is based on obtaining a computed response for which the error measure remains beneath the prescribed tolerance for the complete time history. The adopted error measure is the weighted-residual norm defined in equation (5.3).
6. Some comments about the selection of the basis vectors The strategy for the computation of a suitable basis, where the nonlinear response is to be expanded, can be seen as formed by two different phases: (1) The computation of a set of basis vectors capable of solving the linearized dynamic problem, i.e. the tangent eigenmodes. (2) The computation of a collection of vectors giving the modes of the new tangent eigenproblems that appear when the system deviates from the current position, i.e. the modal derivatives. Concerning the determination of which vectors have to be included in the basis the criterion for their selection should consider, side by side, these previously mentioned aspects. Thus, the two following criteria have to be verified: (1) Convergence from the point of view of a dynamic problem. (2) Convergence from the point of view of a nonlinear problem. The convergence of the dynamic problem can be assured by following similar criteria to those commonly employed for linear dynamic problems. For this reason, a recommended practice is to verify that the projection of the load onto the orthogonal subspace to the span of the basis does not contain any important component [24]. After having defined which modes are to be included, we should analyze which modal derivatives shall be computed and added to the basis. The most desirable formulation would be one in which derivatives of every selected mode are included. However, this criterion could lead to an extremely large number of degrees of freedom in the reduced system. It is not difficult to recognize that only those modes with large enough participation factors will develop large values in the associated derivative terms of the Taylor’s expansion. Thus, we can select derivatives of the modes of higher participation, i.e. the modes of smallest associated eigenfrequency. The enunciated general criteria serve only as a guide for the selection of the basis vectors. Once a particular basis has been chosen, the system is integrated by following the standard Rayleigh-Ritz procedure. However, if the basis has been selected badly, we will not be able to maintain the error norm below an acceptable tolerance, although the basis is recomputed too frequently. Thus, we can see that the process itself indicates when additional vectors must be included in the basis in order to achieve a good approximation to the true response.
7. Numerical examples Three geometrically nonlinear examples were tested: a clamped beam and two arches. All of them were modelled by using a 8-DOF (6 connected, 2 internal) straight-beam element with
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
265
cubic interpolation for both the axial and transversal displacement fields f25]. The reduced solutions were compared with the solution obtained by integrating the complete system of nonlinear equations with the Newmark’s time integrator. Both integrations were performed using the same time step. It is not our intention to discuss whether the results are physically meaningful, but rather to compare both numerical solutions, taking as exact the complete system solution. 7.1. Clamped beam with a sinusoidal varying load In Example 5.1, a clamped beam was solved by using a basis formed by the eight first tangent modes. This example showed that although these modes were updated at each time step, the system’s nonlinear response could not be predicted accurately due to the divergent behaviour associated with the change of basis. Now, the same beam was solved by using the six following basis vectors:
that is to say, the three first modes, the first derivative of the first mode with respect to its own amplitude, the first derivative of the first mode with respect to the second mode’s amplitude and the first derivative of the second mode with respect to its own amplitude. This basis was computed only once, at the initial time instant, and no changes of basis were required during the whole integration. The nonlinear system was integrated up to t = 4 using a time step of 0.025 (160 steps). Figs. 5 to 8 show the computed vertical and horizontal displacement in time of nodes 6 and 3 compared with the exact solution obtained by integrating the complete system of equations. We can see that, simply by including three modal derivatives in the basis, an almost complete agreement with the full system solution was obtained without recomputing the basis, Previously, when eight modes were taken as generalized degrees of freedom, the results diverged even though these modes were recomputed at each time step (Fig. 2). Fig. 9 displays the evolution in time of the error measure, its maximum value being under 5.7 per cent. When eight modes were employed, the error measure reached peaks of 120 per cent, as indicated in Fig. 3. We should remark that not only the vertical but also the horizontal displacement components which have a relatively high frequency content have been approximated correctly. 7.2. Clamped arch with a suddenly appi~ed load One half of the clamped arch shown in Fig. 10 was modelled by 14 equally spaced beam elements leading to 73 degrees of freedom. The example in Section 7.1 displayed a case for which the inclusion of the modal derivatives in the first computed basis allowed to predict the whole nonlinear response. The present example is intended to show that, for most structures, the basis also has to be updated periodically to conform the system’s local behaviour. We also show that the time instant in which the basis needs to be updated is clearly indicated by the error measure. In order to see the need of performing the basis updating, the arch was first solved by
266
8,
‘0.00
t
I
I
I
,
1.20
0.60
I
I
I .80
t
i
2.40
3.00
i
I
3.60
I
4.20
T
Fig. 5. Cantilever beam. Vertical displacement at node 6. The three first modes together with three modal derivatives were employed in the basis. The same set was used throughout the response analysis. ~ reduced model (3 modes, 3 derivatives, no basis change). Complete system; -8 *;
8
6
0.00
’
1
0.60
Fig. 6. Horizontal displacement derivatives, no basis change).
I
I
1.20
at node 6.
I
I
1.80
I
2.40
I
T Complete
I
1
3.00
system; --
I
3.0
I
,
I
4.20
reduced
model (3 modes, 3
267
gl] 0.00
0.60
1.80
1.20
2.40
3.00
s-60
4.20
T Fig. 7. Vertical displacement derivatives, no basis change).
I
Complete
at node 3.
I 0.00
1
I
1.20
1
1
1.80
,
I
2.40
system; --
1
I
3.00
reduced
1
I
3.60
I
model (3 modes, 3
1 4.20
T Fig. 8. Horizontal displacement derivatives, no basis change).
at node 3.
Complete
system; --
reduced
model (3 modes, 3
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
268
“I
00.00
I
I
I
0.60
1
1.20
I
I
I
I.80
I
I
2.40
I
3.00
1
1
I
I
3.60
4.20
T
Fig. 9. The error measure remained derivatives, no basis change).
t
below 5.7% during
the whole
E. 6.895
2
I = 1.0925
I
analysis.
x IO”
h= 5.08
x lO-3
Fig. 10. Clamped
arch. Problem
Residual
x IO- 2
8 =0.245 fib 2.607x
R * 2.54 PO’ 4.7034
--
x IO’
description.
rod IO3
norm (3 modes,
3
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
269
leaving the first computed basis intentionally fixed. That is, no tolerance was set to the error measure and the basis remained unchanged although the error measure reached high values. The first mode and its first derivative were used as basis vectors and a total of 50 time steps of 0.0002 were computed. Fig. 11 presents the response obtained at the apex vertical displacement compared with the complete solution. As it was expected, reduced solution does not agree with the complete one, and the error measure indicates that the results were not being computed correctly (Fig. 12). The same structure was again solved by following the next procedure (Fig. 4): each time the error measure exceeded 0.05, the last two steps computed were rejected and the first tangent mode and its first derivative, were calculated and added to the initial basis made of the first mode and its first derivative. Only two changes of basis were made at times 0.0016 and 0.0034 and 50 time steps of 0.0002 were solved. The vertical displacement of nodes 15 and 8 and the horizontal displacement of node 8 are plotted together with the complete system responses in Figs. 13 to 15. Fig. 16 shows the evolution in time of the residual norm. We can appreciate that now the 4-DOF reduced solution and the complete system solution have a close agreement. 7.3. Clamped arch with an unsymmetrical suddenly applied load Here the same arch as in Section 7.2 was loaded with an unsymmetrical load. The entire was now modelled with 14 elements, leading to 73 degrees of freedom (Fig. 17).
Fig. 11. Clamped arch. Vertical displacement at node 15. The basis formed by the first mode and its derivative was left fixed intentionally. Thus, the results were highly inaccurate. Complete system; -reduced model (1 mode, 1 derivative, no basis change).
Fig. 12. The error measure indicated clearly that the basis should have been recalculated. (1 mode, 1 derivative, no basis change).
-
-
Residual norm
t 3
8 d ::
9 L*
Cn?
;;P 0
4 t 9 c
90.00
0.16
0.32
0*46
0.64
T
0.60
0.96
1.12
*10-2
Fig. 13. Clamped arch. Vertical displacement at node 15. The first tangent mode and its derivative were added to reduced model (1 Complete system; -the first step computed basis at times 0.0016 and 0.0034. mode, 1 derivative + 1 mode, 1 derivative). 270
?
P
Fig.
14. Vertical
derivative
displacement
+ 1 mode,
at node
8. -
Complete
system;
--
reduced
model
(1 mode,
1
1 derivative).
Fig. 15. Horizontal displacement derivative + 1 mode, 1 derivative).
at node
8. -
Complete
271
system;
--
reduced
model
(1 mode,
1
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
272
8 d
8d
Sd ‘;
0
?\ f\ I \I \ I \
A II II
-
a-/
/1
I
I
\
I_ir;pi
f
I
&~I, / i I-,I”1 I 1 ,i’” I (1
\/,I
E-1 4,
,
I
Oo.00 0.16
I
I
I
0.32
I
I
0.46 1
Fig. 16, When the error measure -Residual norm (1 mode,
I
0.64
I
llo-2
I
0.60
I
I
0.96
exceeded 0.05, a basis updating was performed 1 derivative + 1 mode, 1 derivative).
I
(at t = 0.0016 and 0.0034).
I E=6.895
x IOK)
h=S.O8xIO+
I = 1.0925 x 10-s
8 ~0.245
R=2.54
fJ -2.607x
rad IO’
PO’ lx107
Fig. 17. Unsymmetrically
loaded
clamped
arch. Problem
I
1.12
description.
273
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
Various reduced solutions were computed and compared with the exact complete one. The example was intended to show how the response is more closely approximated, when the number of generalized degrees of freedom is increased. All integrations were performed up to t = 0.1, with a time step of 0.0002 (50 steps). First of all, a reduced model made up by three modes and the first three derivatives of the two first modes was solved without performing, intentionally, any basis updating during the whole integration. Figs. 18 and 19 show the computed responses for the vertical component of displacements at nodes 5 and 8. We clearly see that the reduced solution is badly computed, and that the error norm (Fig. 20) indicates the need of changing the basis. Then, an error tolerance of 0.045 was set. The first mode and its first derivative were added to the previously mentioned basis, and the last three time steps were rejected, every time a basis change was performed. A total of two changes of basis were made at times 0.0016 and 0.0024. An evident improvement in the results can be noted (Figs. 18 and 19). At last, the same error tolerance of 0.045 was imposed while the first two modes and their first three derivatives were added and the last three steps were rejected at each basis change. Only one change of basis was required at time 0.0016. A high degree of accuracy was obtained evidencing a clear asymptotic behaviour, when the number of degrees of freedom in the reduced model is increased. Finally, we should comment that the horizontal component of displacements, which in this case has an extremely high frequency content, was estimated within an acceptable margin of
tit .
.
0.32
0.42
T
0.64
llo-2
0.20
0.26
1.12
Fig. 18. Unsymmetrically loaded clamped arch. Vertical displacement at node 8. Complete system; -reduced model (3 modes, 3 derivatives, no basis change); ----reduced model (3 modes, 3 derivatives + 1 mode, 1 derivative); --reduced model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
8 d
1
$1 I 0.00
1
0.16I
I
0.32I
1
0.48 I
I
1
0.64 1
1
I
I
0.80
1
0.99
1
1
1.12
a10-2
Fig. 19. Vertical displacement at node 5. ~ Complete system; derivatives, no basis change); ----reduced model (3 modes, 3 derivatives model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
-+ 1 mode,
reduced model (3 modes, 3 1 derivative); --reduced
(Y c;
0 d
? I I 00.000.16
I
0.22
I
1 0.481 I 1
0.94
I
llo-2
1 0.80r I
0.99
I
I
1
1.12
Fig. 20. Error measure. -Residual norm (3 modes, 3 derivatives, no basis change); ----residual norm (3 residual norm (3 modes, 3 derivatives + 2 modes, 3 derivatives). modes, 3 derivatives + 1 mode, 1 derivative); ---
274
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
275
error (Figs. 21 and 22). It must be remarked that these displacements have an amplitude several orders of magnitude below the vertical displacements amplitude, thus their influence in predicting the overall structure response is minimal. The same behavior may be found in similar linear problems solved by the standard mode superposition method.
8. Concluding
remarks
A computational algorithm for predicting the nonlinear dynamic response of a structure by means of a reduction scheme has been described. In it, the nonlinear system of ordinary differential equations obtained from the finite element discretization, is reduced by employing a Rayleigh-Ritz analysis. Some new aspects have been discussed, i.e. the employment of the modal derivatives as basis vectors, the truncation error introduced at each change of basis, a specially oriented implementation of the Newmark’s time-integration scheme, and a reductionerror measure. The use of the tangent modes as basis vectors and the need for performing a periodic basis updating in nonlinear dynamics have been pointed out by several authors. However, the cost of such method seriously questions its applicability, specially when treating severe nonlinearities, due to the large number of required basis updatings. Besides, in this paper, we have
R d
a d
h
c: 0
d
‘g
d i% ;; 6 d 8 + n
40.00
0.16
0.32
0.40
0.64
1
llo-2
0.00
0.m
1.12
Fig. 21. Horizontal displacement at node 8. Complete system; -reduced model (3 modes, 3 derivatives, no basis change); ----reduced model (3 modes, 3 derivatives + 1 mode, 1 derivative); --reduced model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
276
S.R. Idelsohn, A. Cardona, Nonlineur structural dynamic analysis
Fig. 22. Horizontal displacement at node 5. reduced model (3 modes, 3 Complete system; -derivatives, no basis change); ----reduced model (3 modes, 3 derivatives f 1 mode, 1 derivative); --reduced model (3 modes, 3 derivatives + 2 modes, 3 derivatives).
shown the important role played by the truncation error introduced at each change of basis. We have pointed out that it produces a divergence in the solution when too many updatings are performed. We think that this error can be more drastic in those schemes where not only velocities and accelerations but also displacements are re-expanded at each change of basis. We are now proposing to employ some modal derivatives in the basis together with the tangent modes. These derivatives give to the basis a valuable knowledge on how the system nonlinearity is being developed. By employing them, the number of required basis updatings is highly reduced and the nonlinear dynamic solution path can be followed using only a few basis vectors as generalized degrees of freedom. The proposed method of time integration and basis updating was developed, bearing in mind the purpose of preventing the separation of the reduced solution from the true solution path at any time. Their separation is efficiently detected by the error measure and whenever it occurs the last computed steps are rejected so that the solution path is reached again. Then, the basis is changed and the time integration is continued. The additional computational effort the method requires is negligible compared with the computer time gained during the reduced systems time integration. This is particularly true due to the limited number of basis updatings needed in the proposed formulation. Results concerning the application of this technique to some nonlinear undamped dynamic problems were presented, which showed the power of the method.
S.R. Idelsohn, A. Cardona, Nonlinear structuraf dynamic analysis
277
Appendix A An alternative integration scheme to that presented can be thought of as follows:
in (5.1) in which no truncation
is made
Step 1. Initial time-step computations: 1
art+1 = 48 t 1 = - (l/p At)& + (1 - l/2/3)&, a*aPI+1
‘1 a,+1 = h, +
1 yn+1=0,
(1- y)At& + y Atii:,,
,
(A4
F,,, = Cft’F,,,, - $‘iWi:+~.
Step 2. Iterative cycle: The same as in the scheme (5.1). Step 3. End of the iterative cycle:
This integrator would lead to a continuously growing error in velocities and accelerations as it can be seen from the following analysis. The velocities and accelerations at time t,, can be written as hi, = (1 - #)Li,-I+ &, =
(-l/p
At{1 - 7/2/3)ii,-, + (r/p At)Aa, ,
At)&
64.2)
+ (1- 1/2@&_1+ (l/p At2)Aa,.
By applying the iatter repeatedly
and writing it in matrix form:
+ 5 A’-‘L Aai, It] =Af;:l i=l
(A-3)
where
and A” denotes the nth power of A. Suppose that at instant to a basis change was performed. read:
So, the velocities and accelerations
(A-4)
278
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
The error terms appear because the velocities and accelerations can not be exactly re-expanded in the new basis. By replacing in equation (A.3) we obtain:
44 e (4
].
(A.9
The last term introduces an error that grows exponentially with II. It cannot be balanced out by the first terms because they lie in a subspace orthogonal to the space spanned by e(&) and e(&). So, it is necessary to make the corrections on velocities and accelerations indicated in Step 1 in (5.1).
References [II R.W. Clough and E.L. Wilson, Dynamic
analysis of large structural systems with local nonlinearities, Comput. Meths. Appl. Mech. Engrg. 17/18 (1979) 107-129. in linear and nonlinear analysis, Internat. J. Numer. Meths. PI R.H. Dodds and L.A. Lopez, Substructuring Engrg. 15 (1980) 583-597. techniques in material nonlinear analysis, Internat. J. [31 D.R.J. Owen and O.I.A. Goncalvez, Substructuring Comput. & Structures 15 (1982) 205-213. [41 L. Meirovitch and A.L. Hale, On the substructure synthesis method, AIAA J. 19 (1981) 940-947. in: K.J. Bathe, J.T. Oden and W. Wunderlich, eds., 151 E.L. Wachspress, Two-level finite element computation, Formulations and Computational Algorithms in Finite Element Analysis (MIT, Cambridge, MA, 1977) 877-913. scheme for the solution of systems of equations resulting bl L.C. Wellford Jr. and B. Vahdani, A block-iteration from linear and nonlinear finite element models, Comput. Meths. Appl. Mech. Engrg. 26 (1981) 33-52. An iteration procedure for reducing the expenses of static, elastoplastic and eigenvalue [71 W. Dirschmid, problems in finite element analysis, Comput. Meths. Appt. Mech. Engrg. 35 (1982) 15-33. in structural analysis, Internat. J. Numer. Meths. PI Y.T. Leung, An accurate method of dynamic condensation Engrg. 12 (1978) 1705-1715. 191 R. Guyan, Reduction of stiffness and mass matrices, AIAA J. 3 (1965) 380. [lOI B.M. Irons, Structural eigenvalue problem: elimination of unwanted variables, AIAA J. 3 (1965) 961-962. of Ritz vectors, [Ill E.L. Wilson, M.W. Yuan and Y.M. Dickens, Dynamic analysis by direct superposition Earthquake Engrg. Struct. Dyn. 10 (1982) 813-821. nonlinear finite element behavior using buckling mode superposition, [121 D.A. Nagy and M. Konig, Geometrically Comput. Meths. Appl. Mech. Engrg. 19 (1979) 447-484. [13] E.G. Carnoy, Asymptotic study of the elastic post-buckling behavior of structures by the finite element method, Comput. Meths. Appl. Mech. Engrg. 29 (1981) 147-173. [14] A.K. Noor and J.M. Peters, Reduced basis technique for nonlinear analysis of structures, AIAA J. 18 (1980) 455-462. [15] A.K. Noor, Recent advances in reduction methods for nonlinear problems, Internat. J. Comput. & Structures 13 (1981) 314. [16] A.K. Noor, On making large nonlinear problems small, Comput. Meths. Appl. Mech. Engrg. 34 (1982) 955-985. [17] R.E. Nickell, Nonlinear dynamics by mode superposition, Comput. Meths. Appl. Mech. Engrg. 7 (1976) 107-129. [18] N.F. Morris, The use of modal superposition in nonlinear dynamics, Internat. J. Comput. & Structures 7 (1977) 65-72. [19] S.N. Remseth, Nonlinear static and dynamic analysis of framed structures, Internat. J. Comput. & Structures 10 (1979) 879-897.
S.R. Idelsohn, A. Cardona, Nonlinear structural dynamic analysis
279
[20] L. Landau, N.F. Favilla Ebecken and E.C. Prates de Lima, Comportamiento dinamico nHo-linear de estruturas pelo metodo de superposicao modal, Proc. I Cong. Lat. Amer. Met. Comp. Engng., Porto Alegre, Brasil (1979) 129-151. [21] V.N. Shah, G.I. Bohm and A.N. Nahavandi, Modal superposition method for computationally economical nonlinear structural analysis, J. Press. Vessel Tech. 101 (1979) 134-141. [22] K.J. Bathe and S. Gracewski, On nonlinear dynamic analysis using substructuring and mode superposition, Internat. J. Comput. & Structures 13 (1981) 699-707. [23] A. Cardona, S. Idelsohn and G. Sander, Two stage finite element discretization in structural dynamics, Proc. Third Internat. Symp. Numerical Methods Engineering, March 1983, Paris (1983) 437-445. [24] A. Huck, Application de l’analyse modale par elements finis au calcul de la reponse transitoire et forcee des structures elastiques, Rapport VF-14, Laboratoire de techniques aeronautiques et spatiales, Universite de Liege, Octobre 1972. [25] H.R. Mimer, Accurate finite element analysis of large displacements in skeletal frames, Internat. J. Comput. & Structures 14 (1981) 205-210.