004s7949/M53.00 + 0.00 PcrpmonJoutnais Ltd.
conprte,sd sm,c~wcr Vol.24.No. 4. pp. 559-W. 1986 tinted ia Grat Britin.
AN EXPERIMENTAL-MODAL COMPUTATIONAL TECHNIQUE TO FIND DYNAMIC STRESSES MICHAEL A. TUCCHIOt U.S. Naval Underwater Systems Center, New London, CT 06320, U.S.A. HOWARD I. EPSTEIN~ Department of Civil Engineering, University of Connecticut,
Storrs, CT 06268, U.S.A.
and JOHN
F.
CARNEY 111s
Department of Civil and Environmental Engineering, Vanderbilt University, Nashville, TN 37235, U.S.A. (Received 4 December
1985)
Abatraet-It is shown how the dynamic stresses in a structure can be obtained from measured acceleration data in conjunction with computer calculations of the normal modes of the structure. The proposed technique is demonstrated and the treatment of experimental results is discussed.
NOTATION
B C E, K M P Q
q X, x, i, f Y y z A 4 $ u
mass and damping coefficient matrix damping matrix Young’s modulus stiffness matrix mass matrix forcing function vector stiffness and mass cc&cient ma&x generalized coordinate vector normal function displacement, velocity and acceleration vectors expanded forcing function vector expanded displacement/velocity vector expanded generalized coordinate expanded eigenvector matrix eigenvector matrix stress matrix stress vector INTRODUCTION
This paper shows how the dynamic stresses in a structure can be obtained from measured acceleration data in conjunction with the calculation of the normal modes of the structure. The driving force in the performance of this analysis has been the absence of a method that can be used to determine the dynamic stresses in an entire structure without resorting to hundreds of strain gages. This can be a significant problem when it is not possible to pass the many wires for the strain gages through such structures as the pressure hull of a submarine, the fuselage of an aircraft,
or the body of a missile
t Team Leader. $ Associate Professor. QProfessor.
silo. Penetrations
through structures such as these are usually kept to a minimum because of space limitations and safety reasons. It is shown in this paper that the proposed technique works well, but is subject to problems associated with the current state-of-the-art in gathering displacement data from accelerometers. It is possible to obtain displacement data fairly accurately by other means such as photographic techniques and direct displacement measuring devices (light probes and eddy current meters are but two of the commonlyused laboratory techniques). For the previously cited field problems the accelerometer provides the widest range of application. Techniques exist that can be used to determine stresses from measured displacement data but they must use some form of the strain displacement relationship. These techniques are subject to the error involved in taking spatial derivatives or experimental data and can be applied only to small areas within the field of measurement. Development of the technique to use accelerometer information to determine stresses requires the organization of a modal analysis in a form to accept measured data. That is, the equations of motion of a general structure are organized in such a manner that practical assumptions can be made regarding damping, measurement
limitations
and removal
of rigid
body modes. In the solution of the equations of motion, the application of a generalized inverse is shown to increase the accuracy without increasing the number of measurements. However, the positioning of measurement points is demonstrated to be extremely important.
560
MICHAEL
A.
FORTRAN IV programs have been written to process the data. The method is checked out both analytically and experimentally. A closed form solution for a simply supported beam with suddenly applied end moment and a finite element solution for a single-bay frame with suddenly applied concentrated load are used to provide the analytical support. Base excitation for a cantilever beam and a multibay frame are used in the experimental verification. In all of the illustrative problems, spatial distribution of measurement points, number of measurement points and number of modes considered in the solution substantially affect the accuracy of the calculated stresses. Tucchio et al. [l] have shown that the technique works well for a simple one-dimensional structure (i.e. beam). This paper expands the technique by including a frame for both analytical and experimental verification. It is demonstrated that the proposed technique works well but is subject to problems associated with the current state-of-the-art in gathering displacement data from accelerometers. Note that this work does not reflect the employment of the expensive gyrostabilized platforms for mounting accelerometers. Other techniques exist that can be used to determine stresses from measured displacement data but they must use some form of the strain displacement relationship. These techniques are subject to the error involved in taking spatial derivatives of experimental data and can be applied only to small areas within the field of measurement. Traditionally, measured accelerations have not been used in the determination of dynamic stresses. There are basically two reasons for this, both of an experimental nature. The first reason is the difficulty arising from methods used in the past for removal of rigid body modes from a structure. This was treated recently by Padagaonkar [2]. wherein he avoids the numerical problem associated with accumulative measurement error. The second reason is the difficulty of determining displacements from acceleration data. Several different procedures for this have been proposed over the years. Trifunac [3] and Pecknold [4] attempt to produce corrected accelerograms by giving estimates of initial ground velocity, ground displacement and ground acceleration. According to Pecknold (41, “Unfortunately, there is no completely satisfactory way to avoid the problem at present.” In the present paper. the parabolic type correction, as established by Trifunac [3], is used with partial success. The theorv of the method is developed from a modal solution of the general equations of motion of an elastic structure in matrix form. A generalized inverse is used to solve for the stresses based on a rectangular eigenvector matrix. This allows the number of accelerometers to be less than the number of modes. The generalized inverse provides a least squares approximation to the overspecified set of equations. A closed form solution for a simply sup-
TUCCHIO et al.
ported beam with suddenly applied end moment is used to provide the analytical support. Considerations of experimental limitations are investigated and discussed in the demonstration problems. The accuracy of measurements as related to the spatial distribution of accelerometers and number of modes is shown to be of major importance. In the solution of the equations of motion, the application of a generalized inverse is shown to increase the accuracy without increasing the number of measurements. However, the positioning of measurement points is demonstrated to be extremely important. The method gives excellent results within the limits of application delineated in the paper.
MODAL
ANALYSIS
The required modal analysis is based on the method developed by Hurty and Rubinstein [S] and Beliveau [6]. It will be assumed that the equations of motion can be represented by
WI{*) +[Cl{i) + WI(x) = {P(f)),
(1)
where M, [C] and [K] are the mass, damping and stiffness matrices, respectively, the force, {P(r)}, is a function of time, and {x} is the displacement vector. The damping force will be assumed to be linear. For large amounts of damping, [Cl, it is not possible to uncouple eqn (1) directly by the separation of variables technique, {x} = [+]{q}, where [4] is the eigenvector matrix, and {q} is the generalized coordinate. The damping need not be proportional, e.g. for any component of (4’) for an undamped system, or for one in which damping is proportional, {4r) is distinguished from other components by amplitude only, the phases being equal or 180” apart as determined by sign. By utilizing the identity [M]{i}-[M](i}=O,
(2)
the n eqns (1) can be expanded to 2n eqns (3): [~l~]{r}+[~l~]~~~={bj,
(3)
or, if eqn (3) is written more compactly
PI{);) +[Qlb1=
{YJ-
(4)
where
PI = [;I$]:
{Q) = [#+
and
(5)
Experimental-modal
computational
In order to uncouple eqn (4) and use the separation of variables technique, y must be of the form: {Y} =
[Al{zI.
technique to find dynamic stresses
561
matrix. If the earlier complex eigenvector matrix were retained, it would still be possible to determine a pseudoinverse, as shown by Goldstein [8].
(6) Removal of rigid body modes
The transformation
matrix [A] is constructed column by column using the 2n eigenvectors {y(h), {u -(I)}, {y(z)}, {.v-‘2’},..., and {y ‘“‘} and { y -(“)}, and the bar signifies complex conjugate. The generalized coordinate, z, is also complex. Investigation of the availability of measurement techniques to determine {JJ}, where it is necessary to determine phase as well as amplitude, leads to the conclusion that accurate phase information cannot be obtained with today’s technology unless the application is limited to narrow frequency bands where special application equipment for a particular frequency response can be designed. Because of this and the inaccuracy in the determination of displacement from accelerometer measurements, this study is limited to cases where the C matrix is small. The modal analysis, therefore, is unaffected by the presence of damping. This has application in many areas including steel structures, where it is common to experience less than 6% critical damping. With this restriction, the solution of eqn (1) is simply expressed as {xl=
Wlb71-
(7)
The acceleration amplitude is the only measurement needed to calculate the displacement {x} of eqn (7). The eigenvector can be calculated using any suitable technique. Generalized inverse
In eqn (7), it is necessary to assume the number of modes to be considered. This is common in all modal analyses. In many structural applications, only the first few modes are of interest. This paper will be restricted to consideration of such problems. In anticipation of problems associated with the accuracy of {x), let it be assumed that fewer measurement points are available than the number of modes to be considered. That is, let eqn (7) be written as iXJN~c=[91~~c.~hi0DEs
14l~h40Dm3
(8)
in which NAC = number of accelerometers and NMODES = the number of modes. Now [4] is rectangular and, in effect, {x} is overspecified. To maximize the effectiveness of the number of measurements, so that as many modes as possible can be considered, it is assumed that the number of modes, NMODES, is equal to or greater than the number of measurement points or accelerometers, NAC. The determination of the generalized coordinate, {q}, from eqn (8) involves the pseudoinverse of [#I. This can be found by using a least squares polynomial fit [7l to an NAC x NMODES
Removal of the rigid body modes is necessary so that [K] will not be singular. This is done in a manner that leads to a rigid body solution in which cumulative numerical integration errors are eliminated. Padgoankar [2] developed a nine accelerometer scheme which does not use information from a previous time step. For planar motion, only four accelerometers are required. In the present paper, Padgoankar’s method is used directly. Normalized stress matrix
At this point in the development, the generalized coordinate has been determined from eqn (8). A stress matrix is now constructed that is compatible with the form of the generalized coordinate. The stress analysis is performed on a modal basis so that a stress matrix is formed with columns that represent the stress in the structure due to a displacement equal to the eigenvector for that mode. Theoretically, the number of positions considered is not limited. For example, let X represent a normal function for a beam; the flexural stresses can be obtained from eqn (9). e,=
-E
d2X,
-,
’ dx2
where u” represents the stresses for the n”’ mode. A stress matrix, [#I, for a number of modes equal to NMODES and for a number of positions equal to N, can be constructed in the form rg,,
u12 -----
L 1--------
OlN
1
- uNMODE.N
The stress matrix can be used in conjunction with the generalized coordinates to determine the required stress distribution,
(11) In eqn (1 I), {q} is the generalized coordinate experimentally determined by eqn (8) and the unsubscripted {c} is the desired stress. It should be noted that unless one knows a priori the direction of the input displacement, the corresponding calculated stresses can be in phase or 180” out of phase with the measured results. Since the [+] matrix is constructed based on an assumed eigenvector positive sign convention, the actual displacement might very well be of opposite sign.
562
MICHAEL
A. TUCCHIOer al.
I A
123456 -I
l-&x
7
J-08
FP)
I
”
L
Ll4
L/2
-9 t
+l
8
////////////J L = 10.0 ” (254 mm)
Member m is Rigid Point 4 is Midpoint of Horizontal Member Ax = ~18
Resonant Frequencies Timoshenko
[iol
Finite Element
[9 ]
Fngurncy in Hwtz
230.5
230.3
295.4 431.6 921.5 1050.4 1274.0
295.1 430.1 919.6 1047.7 1277.7
Fig. 1. Plane frame with a step function transverse load.
ANALYTICAL
SOLUTION
To test the method, the stresses in a plane frame are determined. The loading consists of an instantaneously applied concentrated force. A finite element program, NASTRAN [9], is used to obtain the solution to the problem of an impulsive force applied to the plane frame depicted in Fig. 1. The displacements of certain points, as determined from NASTRAN, are referred to as ‘measured’ displacements. The stresses are found by processing those data as an ‘experimental’ solution via eqns (8), (10) and (11) and comparing the result with the NASTRAN stress output. The finite element model was verified by comparing the first six free vibration frequencies with the Timoshenko solution. Results are shown in Fig. 1. The duration of the impulsive force, i, is 0.0001 set which corresponds to a circular frequency of 1592 Hz. Note that this is slightly higher than the sixth resonant frequency of 1278 Hz. The stresses at the midpoint of the horizontal member as determined from the finite element solution are plotted in Fig. 2. Also included in the figure are various ‘experimental’ solutions. As can be seen from Fig. 2. when six ‘measurement points’ were used
(NAC = 6), together with six modes (NMODES = 6), the resulting 6 x 6 [I+%] matrix gives excellent results for this case. Generally it is desirable to reduce the number of measurement points in order to reduce the number of accelerometers, wires, amplifiers, signal conditioners and recorder channels. Reduction of the size of the eigenvector matrix to a 3 x 3, where ‘measurements’ were taken at points 1, 2 and 3 (see Fig. 1) for NAC = 3, does not approximate the true stress very well as can be seen in Fig. 2. However, the dashed curve in the figure shows that excellent results are obtained when the eigenvector matrix is expanded to a 3 x 6. (Six modes are used with only three measurement positions.) Therefore, the three measurement positions seem sufficient for this example. As discussed by Tucchio el al. [l] the accuracy of measurement of actual displacements is usually not high, ahd when the spatial separation is close relative displacements may be meaningless. Therefore the decision was made to investigate a coarse spatial distribution of measurement positions and to use this data to determine stresses. The measurement points were placed close to the antinodes of the first node of the structure (points 4, 7 and 8 on Fig. 1) rather than using the closely spaced positions 1, 2 and 3. When
Experimental-modal computational technique to find dynamic stresses
563
this is done, eqn (8) becomes:
rx41
r
1.oOOO
1OOOO 0.002471 0.9970
l.OOOO
-0.9982
0.9992
-0.5024
0.003166 -0.003211
0.5003
-0.002727
Note that the numbers in the fourth column of the eigenvector in eqn (12) are very small. This is because points 4, 7 and 8 are close to the nodes of the fourth mode of this structure. Stresses determined from the {q} matrix calculated with the eigenvectors of eqn (12) give erroneous results. The cause of this is the numerical problems associated with the inverse of a matrix that is nearly singular. Satisfactory results are obtained by substituting position 9 (Fig. 1) for position 8. Equation (8) then yields: l.OOOO
l.OOOO 0.002471 0.9970
l.OOOO
-0.7074
-0.5024
SWOSB
-0.03837
0.2153
0.3210
-0.2145
(12)
0.3125 J
146
-0.003211 1.oOOO
-0.4319 0.2153 0.9586
r q, 1
-0.03837 1
-0.5414
L
I
.46 -
A\ I
.30 -
A
.15 -
O-
\
I
2 -.15 -
A4
\
Spamoly Swcad-13 x 61
I
I
I
I
I 10
1 12
6
A\
i 14
Tima
Sot x lo-4
Cloroly Spacod p Y 31 -.30 -
-.4s -
-.w
-
-.75 -
16x 61 -30
-
Fig. 2. Comparison of stresses at midpoint of horizontal member of plane frame. C.&S.
2414-O
I (13)
0.3210
.n
MPa 30 -
J
Therefore eqn (13) is used to determine {q} and the stresses are calculated again at the midpoint of the horizontal member of the frame. The triangles of Fig. 2 show the results of this analysis. These stresses, although not as accurate as those determined for the closely spaced measurement points, are good approximations to the finite element solution. However, the coarse spatial positions 4, 7 and 9 have the
0.003166
0.4824
0.7645
-0.4319
MICHAEL
564
O-
Accelerometer
o-
Slrain Gage
A.
TIJCCHIO et al.
H = 1.0”(25.4mm)
15H
40H
#l( re
t 1SH
1
113
c
u141 I -
-
15H ~---j L(
30H a I1
Colocated Channels 3 6 4 and
5 16
(ON BASE) are used for redundancy.
Fig. 3. Multibay frame-plan of being suitable for actual measurements. For this example, the use of the generalized inverse provides an excellent method of minimizing the number of measurement points. At the same time, caution must be exercised in the spatial distribution of these points. advantage
EXPERIMENTAL
RESULTS
The test specimen is the multibay fraine shown in Fig. 3. The frame is of welded construction and is coated with metallic copper tape to provide damping. The eigenvector solution for this structure was obtained with NASTRAN [9] using 425 degrees of freedom. Base motion is imparted to the frame with an L.A.B. Corp., Model SPA-400 shock machine. The digitized raw data for typical accelerations are shown in Fig. 4. Integration of the accelerations.to obtain displacement is done using Simpson’s rule. However, this does not give satisfactory results for the purposes of this investigation. A review of the recent literature [3,4] and discussions with the measurement community revealed that this is a common problem. Noise and some low frequency initial conditions can seriously affect the measured response. In this work, it was decided to attempt to solve the integration
6H x H/2
view.
problem by using a simple correction technique. A constant initial acceleration is assumed and the velocity or displacement results are corrected based on this initial condition, Physically, the shock machine carriage is driven down an arbitrary amount and is decelerated by striking lead pads between the moving carriage and the stationary anvil. The velocity of the carriage should be zero after a short period of time. This zero velocity condition is imposed on the first integration of acceleration as a boundary value. This results in a corrected velocity which, upon integration, gives the proper displacement. Figure 5 shows the output from four typical channels: I, 3, 6 and 7. The relative displacement of the frame with respect to the carriage must be zero after a period of time. The drift, or accumulated error, results from an initial small value of acceleration and is assumed to be constant. The velocity correction is linear and the displacement correction is parabolic. Corrected displacements are shown in Fig. 6. It is seen that this correction technique works very well for channels 1 and 7, fairly well for channel 3 and not very well for channel 6. This observation is based on the attainment of steady state. Unfortunately, it is difficult to select the proper displacement formulation for correction and only by engineering inspection, considering the physical realities, can the quality of the result be judged.
Experimental-modal
computational
technique
MIS2
MIS2 CHANNEL
50-
565
find dynamic stresses
to
1
=5r
CHANNEL 3
I
0
-5o-
-lOO-
-150
-
0
-2000m
1
.20
SECONDS
A0
.60 I
SECONDS 125
125
25
25
-25
0 SECONDS
.20 .40 SECONDS
.50
Fig. 4. Multibay frame-accelerations.
Experimental data for comparison are given by channels 10, 1I and 12 of Fig. 7. Processing of the stress data uses the normalized stress matrix and eigenvector matrices determined from the NASTRAN [9] eigenvector analysis. These matrices are
The arrangement of the elements in eqn (15) depends upon the sequence of the terms in the displacement vector. The stresses calculated using the processing scheme such that the sequence for the displacement vector is 7-3-6 are shown in Fig. 8. For
- 3883.
4855.
26988.
- 20902.
20602.
85166.
- 3676.
501.
- 28594.
10540.
-493380.
41460.
-26W
- 17917.
8224.
- 77850.
11491.
- 0.74846
0.55155
0.98259
-3172.
(14)
and 1.0
0.19511 0.20688
0.13420 -0.24302
-0.56086
-0.42875
-0.78317
-0.01702
-0.63006 -0.73593 0.59193
-0.076716 0.061975 0.913416
.
(15)
566
MICHAEL
A.
%CCHIO
el al.
I
I
I
Experimental-modal computational technique to find dynamic stresses
567
CHANNEL 10
1
0
4
1
20
.40
30
‘SECONDS
CHANNEL 11
I 0
-20 SECONDS
,
17r
0
Fig. 7. Multibay frame-measured
one measurement point (channel 7) the comparison of results is very good for channel 10, considering one mode, and also for channels 11 and 12, considering two modes, as in Figs 8a and 8b. Inclusion of higher order modes give lower values than measurements. Similarly, when two measurement points (channels 7 and 3) are used, the results compare favorably if the lower modes only are used (see Fig. 8~). Expanding the matrix to a 3 x 3 by including channel 6 gives poor results because of the inaccuracy of the corrected displacement from this channel. The displacement vector sequenced 3-7-6 also gives good stress results as compared to those of channel IO in Fig. 7 for one measurement position (channel 3) and one mode. Additionally, good agreement is obtained when two measurement positions and two modes are included. Again, expanding the measurements to three gives results similar to those encountered above with a displacement vector 7-3-6. When channel 6 is processed first, the results are not in agreement with measurements. This is expected since the quality of the displacement is definitely not as good as for either channel 3 or 7. The curves have not been included here. From the structure of eqn (15) it can be seen that if the first measurement,
CHANNEL 12
.20 .40 SECONDS
.60
stresses.
x6, is inaccurate the solution for {q 1 is inaccurate. Subsequently, the solution for {u) is also inaccurate.
CONCLUSIONS
The purpose of this investigation was to develop a method to determine the dynamic stresses in the structure using processed accelerometer data together with a modal analysis of a structure. In the experiment, several of the comparisons were very good. The analytical demonstration showed excellent agreement. Unfortunately, the use of displacement data derived from accelerometers is not of sufficient quality to reliably predict the stresses in a structure. The state-of-the-art in deriving displacements from integrated accelerometer data has not progressed to the point where the scheme derived in this paper can be applied consistently. This is borne out by recent literature as well as the work presented here. However, the scheme itself is viable as has been shown by the sample analytical problem. The placement of the accelerometers on the structure must be chosen carefully and should be based on
568
MICHAEL A. TIJCCHIO et al.
CHANNEL 10
NMODES = 1 NAC = 1
-34-
I
0
.20 .40 SECONDS
A0
17-
CHANNEL ii
O-
0
83
%
Z-17 t;;
+I -17-
G -34
t
w
0
SO
0-
.60
SECONDS
SECONDS
CHANNEL
10
(b)
0-o SECONDS
CHANNU 11
17
I
0
17
I
20
.40
.60
CHANna
0NO
. SECONDS Fig. 8 (continued
SECONDS on following
page).
12
Experimental-modal computational technique to find dynamic stresses
569
CNANNEL 10
(cl NWOftB-2 NACI2
I 0
.20
.40
.sb
SECONDS
-2Te-G
0
CNANNU
l.2
-34L.0
.20
.40
.60
SECONDS
SECONDS
Fig. 8a-c. Multibay frame calculated stresses, processing order 7-3-6. Note: compare to measured values, Fig. 7.
the eigenvector solution. It is possible to locate the measurement positions in a manner such that a nearly singular eigenvector matrix is obtained. The resulting column of very small numbers in the eigenvector matrix significantly affects the quality of the inverse, and the magnitude of the resulting generalized coordinates may be very inaccurate. This problem arises because of the number of significant figures obtainable from the accelerometer integration in conjunction with the small numbers in the eigenvector column. It would at first appear that placing additional accelerometers at antinodes to keep at least one member in the eigenvector column large could provide a solution. Analytically, this is proven to be a very good technique. In the experiments, however, the small number of significant figures obtainable with the integrated accelerations precluded this strategy. The futility of accurately measuring relative displacements between a large number of positions is the basic problem. The processing order of the displacements was also shown to affect the accuracy of the calculated stresses. This occurs when two or more accelerometers are used to obtain the stresses and the quality of the integrated acceleration data is not the same. Processing the least accurate displacement first gives erroneous calculated stresses.
REFERENCES
M. A. Tucchio, J. F. Carney III and H. I. Epstein, Dynamic stresses from experimental and modal analysis. ASCE/EMD Specialty Conference on the Dynamic Response of Stru&res, Atlanta, GA (1981). _ 2. A. J. Padaaoankar. K. W. Krieaer and A. I. Kina. MeasuremeZ of angular accel&ation of a rigt% body using linear accelerometers. J. uppl. Mech. 42 (1975). 3. M. D. Trifunac, Low Frequency Digitalion Errors and a New Method of Zero Baseline Correction of StrongMotion Accelerograms. Earthquake Engineering Research Laboratory, EERL 70-07: California-brstitut~ of Technoloav. __ Pasadena. CA f1970). 4. D. A. Pecknold and R. Riddeh, E&t of initial base motion on spectra. J. J%gng Me& Diu. ASCE 104, 1.
485-489 (I 978). 5. W. C. Hurty and M. F. Rubinstein. Dynamics of
Structures. Prentice-Hall, Englewood Cliffs, NJ (1964). 6. J. Reliveau, Identification of viscous damping in structures from modal information. J. appl. hfech. 43 (1976). 7. M. J. Goldstein, Linear least squares estimation using the generalized inverse. UNIVAC Scientific Exchange Meeting, Atlanta, GA, IS October 1970. 8. M. J. Goldstein, Reduction of the pseudo-inverse of a Hermitian persymmetric matrix. J. Math. Compur. 28 (1974). 9. NASTRAN-NASA
Structural Analysis-NASA SP. 221 (OS) National Aeronautics and Space Administration, Washinaton. DC. December 1978. 10. S. Timoshenko, .Vibrotion Problems in Engineering. D. Van Nostrand Company, New York (1955).