PULSE T E S T I N G FOR D Y N A M I C A N A L Y S I S A n Investigation of Computational Methods and Dzficulties W I L L I A M C. CLEMENTS, JR.,
AND KARL B. SCHNELLE, JR.
Vanderbilt University, Nashville, Tenn.
Pulse testing has been suggested as a simple method of obtaining system frequency response. This paper presents details of certain theoretical and practical aspects. It is simple in principle and in application. The system input is given an arbitrary closed pulse and the resultant output recorded. Under proper conditions the entire useful portion of the frequency response curve can be calculated from these data. Working equations for reduction of pulse data to frequency response are presented, using various approximation techniques. Problems concerned with the transcription of data from experimental curves and application of various data-reduction procedures were examined, using pulse tests on an analog computer and analytically computed pulse response. response data.
Specific recommendations are given to maximize the recovery of frequency
of frequency response methods for the dynamic analysis of physical systems are well known. Some of the most fruitful techniques for process control and analysis are based on frequency response methods ( 7 , 2, 3). One especially interesting use is the comparison of experimental frequency response data for a system with the response calculated for a postulated model, to obtain a quantitative estimate of the range of validity of the model in describing the system. This approach is useful even for distributed parameter systems, where it is often difficult to obtain complete time domain solutions to the differential equations (4, 5, 77). The classical method of obtaining a frequency response experimentally is based on its definition. The system input, x ( t ) , is varied in a sinusoidal manner a t frequency w ; the system output, ~ ( t )is, measured simultaneously. If the inputoutput relationship is described by a linear differential equation with constant coefficients, then y ( t ) will be sinusoidal and of the same frequency as x ( t ) . T o describe the effect of the system mechanism on the input signal, it is convenient to form the ratio between the amplitudes of the output and the input; this is called the magnitude ratio (M.R.). The phase difference, 4, between y ( t ) and x ( t ) is also of interest. T o obtain enough points to define the frequency response curve adequately, it is necessary to make a number of sinusoidal tests at various frequencies. This is a significant disadvantage in many cases. Not only is it awkward to generate reasonably pure sine waves in the input variable, which may be temperature, pressure, concentration, etc., but the sinusoidal upsets thus introduced cause a corresponding disturbance in the normal operation of the process, and this disturbance must be maintained until the sinusoidal steady state is reached. This is obviously very undesirable in the testing of process equipment in actual plant use. Furthermore, sinusoidal methods are incompatible with many systems-for instance, one involving a fast chemical reaction. The pulse testing method was developed to circumvent difficulties of these types. The saving in testing time and money, in equipment needed, in manpower, and in upset of normal operation can be very great. Pulse testing has been used in a number of studies reported in the literature. The earliest of such test programs were concerned with analysis of aircraft dynamics during actual HE ADVANTAGES
94
I & E C PROCESS DESIGN A N D DEVELOPMENT
flight (78, 79, 20). Other studies have been made on servomechanism and control apparatus (6, 75). Several papers of particular interest to the process industry report the use of pulse testing in investigating the dynamics of certain chemical engineering processes (77, 73, 27, 22). Pulse methods, however, have not been applied nearly so widely as one would expect, when their many advantages are considered. This unfortunate situation is probably due to the lack of published material aimed at familiarizing the practicing engineer with the theoretical and practical aspects of the calculations involved. There needs to be made generally available a treatment of pulse reduction techniques that carefully explains the physical situation involved, in such a manner that it will be immediately of use to those who must perform the pulse tests and program the calculations on a digital computer. This paper discusses the topic of pulse testing from these viewpoints. I t is the authors’ hope that it will be of value in aiding instrument people, practicing engineers, and other workers involved with process and instrument dynamics, in applying pulse testing techniques in their own work. Mathematical Background
Frequency Content and Pulse Reduction Equations. Draper, McKay, and Lees (7) describe in detail the general mathematical background for certain methods of numerical calculation of Fourier transforms and, through these, the reduction of pulse data. Lees (72) describes the general area of dynamic measurement. Smith and Triplett (20), Hougen and Walsh (77), and Dreifke (8) discuss pulse frequency content and the effect of the frequency content of the input pulse on frequency response information recovered from pulse tests. Additional discussions of Fourier transforms and the Fourier integral are widely available ( 7 , 70, 74, 77). I t can be shown (74) that the frequency content, S(w), of an aperiodic function is given by its Fourier transform
From Equation 1 it follows that S(w) is a continuous function of w, under certain mathematical restrictions not discussed here. I t is convenient, in discussing the frequency content of
RECTANGULAR HALF-SINE
PULSE
A
RAMP PULSE SQUARED
n
PULSE
RAMP PULSE
TRIANGULAR
A
PULSE
DISPLACED COSINE PULSE
A
WEIGHTED DISPLACED COSINE PULSE
A
SQUARED TRIANGULAR
PULSE
CUBED TRIANGULAR PULSE DIRAC
IMPULSE
FUNCTION
A
A
t
0.6
0.5
0A
0.2
0.2
0.1
0 2
0
6
4
8
12
IO
14
16
20
* (FREQUENCY
Finure 1 .
X
PULSE WIDTH)
Freauencv content of some mathematically exwessible Dulses
different pulses, to normalize S(w) by dividing it by the content a t zero frequency ( 7 7,20) : (3) S(W)” = S ( W ) =
S(0) ~
S(w)
[total area under x ( t ) ]
This quantity is extremely useful in comparing the relative frequency contents of various pulse shapes. A comparative plot of the normalized frequency content magnitude, ~ ( c J ) , , for a number of pulse shapes is shown in Figure 1. T o reduce clutter of the plot a t high frequencies, the curves are not plotted past their first zeros. This paper is concerned only with pulse functions-functions that are zero a t time zero and assume definite values for a finite length of time (called the “pulse width,” T J , after which they return to zero and remain there indefinitely. From these statements it follows that
T h e magnitude ratio is the ratio of output frequency content to input frequency content, and the phase angle is the phase where the subscripts indicate angle between S(w), and S(w),, frequency content of system output and input, respectively. Thus
(4)
and VOL. 2
NO. 2
APRIL
1963
95
9 = angle [S(w),] - angle [S(w),]
= uz. - us
(5)
where S(w) =
s(w)e1u(w)
(6)
Applying the identity e--IWt = cos wt - j sin wt to both integrals in Equation 4 and subsequent algebraic reduction give the basic equations used for obtaining frequency response data from pulse tests: M.R.
(w) =
d R e 2( w )
$, = tan-’
[“I
+ Irnz ( w )
(7)
Re(w) (9)
where
Irn(w)
B =
c
~
so”
y(t)
= JoTz
D =
A D - BC C2 D2
=
+
sin
ut dt
x ( t ) cos wt dt
soTz
x ( i ) sin wt dt
The computational problem is the evaluation of A , B, C, and D from experimental data. The most direct method is the application of a quadrature formula, such as the trapezoidal rule
etc., found in Equations 12 through 15 oscillate faster and faster. Thus the approximation of the curves, no matter which quadrature formula is used, rapidly deteriorates so that with a given accuracy in x ( t ) and y ( t ) a value of frequency will eventually be reached where the quadrature formula becomes completely useless. This high frequency restriction is common to all trigonometric integrals of these types, and it can be removed by several special approximation methods. These methods are based on approximating a small portion of the pulse curves, x ( t ) and y ( t ) , followed by integration of the approximating expression multiplied by sin ut or cos ut. The resulting equation gives the area under one segment of the product curve. The total area is obtained by summation of all such segments contained between the limits of integration. Since the trigonometric integration is done analytically, large values of w do not introduce error into the quadrature. The only factor left which can influence the accuracy is the degree of approximation of the pulse curve by the chosen approximating function. Analytical integration using the approximate expression for x ( t ) can also be applied to the entire Fourier transform at once; however, this result for the more sophisticated approximations may be too complicated to be convenient for computation. In such cases it is best to approximate the individual product curve integrals and then combine them to get M.R., 4: etc., as has been done throughout the study reported in this paper. The simplest approximation to a pulse curve which has been used is the stepped curve, or rectangular-area, method (7). This involves the approximation of the function of time by the tops of a series of rectangles of equal width, and is the least accurate approximation that could be expected to have any practical use. The stepped-curve method is simple enough to be applied to the entire Fourier transform, giving for the ratio of transforms expressing frequency response :
i =1
directly to the trigonometric integrals in Equations 12 to 15, giving n- 1
A
=
y ( k A t ) cos wkAt
At,
A second widely used method of numerical calculation for these transforms is the approximation of the pulse curve by the tops of a number of trapezoids of equal width (7). This approximation can also be applied to the entire Fourier transform, giving
k=l
n- 1
B = at,
C y ( k ~ tsin) w k ~ t
k=l
c
n-1
C
=
Atz
x ( i A t ) cos wiAt
i = l m- 1
D = Atz
x ( i A t ) sin wiAt
(20)
i = l
where (excluding zero end points), n- 1 points are read on the points are read on the x ( t ) curve at intervals At,. Substitution of a particular value of w in Equations 17 to 20 and performing the indicated summations using experimental x ( t ) and y ( t ) data will yield A , B, C, and D. M.R., 9, and s ( w ) ~then follow from Equations 7 to 11. The accuracy of the evaluation of the integrals should be improved by using Simpson’s rule for the quadrature. However, as frequency increases, the product curves x ( t ) sin ut,
y ( t ) curve a t intervals Atu and rn-I
96
l & E C PROCESS D E S I G N A N D D E V E L O P M E N T
a result which differs from the previous approximation to G(j w ) only in that the trigonometric functions multiplying the summations are squared. As Hougen and Walsh ( 7 7 ) have noted, Equation 22 gives the exact transform when applied to a pulse composed of straight-line segments, if segment boundaries coincide with 4t boundaries. Equations 21 and 22 may be compared to the previously discussed method in which the product curves themselves are approximated by trapezoidal segments, giving
Y
G(jw) =
k=l rn Y
i=1
Programming the Calculation
Thus the trigonometric functions
act as correction factors by which the summations themselves are multiplied to account for the variation in w, under the stepped-curve and trapezoidal area assumptions, respectively (72). All three approximation techniques produce the identical result for G ( j w ) if Atz and At, are equal, for in this special case the correction factors as well as the At values in the numerator and denominator of Equations 21: 22, and 23 cancel each other. Since practically all pulses that would be encountered in actual plant tests will be smooth curves, a suitable curvedsegment approximation to x ( t ) and y ( t ) should allow a n adequate approximation with less curve reading than the methods mentioned previously. T o circumvent difficulties encountered when w in the product curves becomes large and to provide a curved segmentation of the time functions, Filon ( 9 )developed
T o facilitate the computations required by Equations 7 through 15, the authors have written programs for the IBhf-650 computer which carry out the pulse evaluation using both the trapezoidal rule (Equations 17 to 20) and Filon’s formula (Equations 24 through 31), assuming the input and output pulses for a particular test are given by Figure 2.
t
a quadrature formula for f # ( t ) sin wtdt and Ib#(t) cos wtdt, based on approximation by parabolic segments as in Simpson’s rule. However, in Filon’s method the 2/3, and ‘/3 coefficients of Simpson’s rule are replaced by trigonometric functions of w1.t. To apply Filon’s formula? the curve is divided into a n odd number, 2n 1, of points a t intervals At. Denote these points by I O ,I,, I z , . . . I?, where IO = $(a) sin wa and I*,, = $ ( b ) sin ab. Let
+
SZ,, =
1
lo
+ 12 + . . . +
S?n- I = J I
2
+
+ Ja + . . . + Z?ne
1 2
- z2,z
(24)
I
(25)
(26)
= PI
Figure 2. dissection
a, p, and y must be computed from Taylor’s series expansions of Equations 27 through 29 when 8 is small, to prevent loss of
significant digits (9). Then a [ $ ( a ) cos wa
- $ ( b ) cos w b ] + P S2,
+
Y
(30)
&“-I}
If Io is now denoted by $ ( a ) cos wa and Izn by $ ( b ) cos wb
*/2)1
+ P + Y &$2,
11
(31)
Filon’s formulas give the parabolic approximation to the pulse curves with an accuracy that is independent of the value of w , and reduce to Simpson’s rule as w 0. Like Equation 22, Filon’s method gives the exact Fourier transform for pulses composed of straight-line segments, and also computes exactly for pulses composed of parabolic segments, when segment boundaries and At boundaries coincide.
-
Typical input and output pulses and their
In many cases, as in this illustration, one portion of the curve changes more rapidly than another portion. In such cases it is convenient to read a large number of points in the portion of fastest change, and less points when the change is more gradual. The programs of the authors are written to accept two different increments of time on both the input and output pulses. Thus the program input consists of the input data points, xo through x n ; the output data points, yo through y n ; the four At values, two for x ( t ) and two for y ( t ) ; and the number of points read using each At. These quantities are all indicated on Figure 2. The program then calculates numerically the area under the product curves as specified by Equations 12 through 15. These areas are calculated separately for each At portion and then summed to yield the total values of A , B, C, and D. Magnitude ratio and phase angles are then computed using Equations 7 and 8, beginning at a specified value of w and incrementing it by a specified amount for each additional point until the highest value desired has been reached. Also calculated and punched simultaneously is S ( W ) ~ , the frequency content of the input pulse at the w value in question, as well as the real and imaginary parts of G(j w ) . Detailed descriptions of both the trapezoidal rule and Filon’s formula programs are available from the authors, VOL. 2
NO. 2
APRIL
1963
97
Practical Limits of Pulse Evaluation Calculations. When a pulse input is applied to a system, the system will be excited a t all frequencies except those frequencies, wo, at which the frequency content, $(a),,becomes zero. At these points the summation of the positive area under the product curves x ( t ) sin wt,x ( t ) cos at, y ( t ) sin at, and y ( t ) cos ut just equals the summation of negative areas, and their integrals are simultaneously zero, reducing the expressions for magnitude ratio and phase angle to the indeterminate form
0
T.
In practical computations, small inaccuracies arising from errors in reading points from the experimental curves and truncation errors inherent in the numerical approximation of the areas under the product curves will usually cause the computed areas to assume some meaningless small values at the wo points. Thus highly erroneous values of M.R. and 4 will be found if calculation is attempted a t these points. Furthermore, as these integrals approach zero, some frequency will be reached at which enough significant digits will be lost, because of the subtraction of two nearly equal numbers in the product-curve area evaluations, that a deviation from the true frequency response curve will occur a t frequency values less than wo. One might suppose that such deviation would be readily apparent, in that the calculated points would start to scatter noticeably as soon as enough significant digits were lost in the area approximations. This scatter does in fact occur, but only well after the onset of error. The computed curves almost invariably remain smooth even as they depart from the locus of the true system frequency response (see Figures 4, 5 , 6, 7, and 8 for examples). This effect would lead the uninformed investigator to trust his calculated points up to a n excessively high value of frequency, and might result in misinterpretation of system dynamics a t higher values of w. Since very little information has been published concerning the effects of using various data reduction procedures, the number of significant digits retained in the pulse-curve transcription, and other practical considerations upon the location of the point of first deviation, it was necessary to undertake a study to clarify some of these points before maximum benefit from the pulse method could be expected. Experimental Procedure
I t was first decided to pulse-test simple first- and secondorder systems, using an analog computer to simulate the transfer functions. Rectangular, half-sine, and triangular pulses were produced by feeding the output of a Hewlett-Packard function generator into a half-wave diode rectifier. Squared-
+;%;% FUNCTION
6AL5 D I O D E
GENERATOR
\
1
(SWITCH CLOSED 8EFORE POSITIVE P U L S E AND OPENED IMMEDIATELY AFTER)
T O VISICORDER CHANNEL I
triangular pulses were produced by feeding the triangular waveform from the function generator into a Donner multiplier unit. Randomly shaped pulses were generated when desired by merely turning by hand a potentiometer connected across a voltage source. Both the input pulses and the output pulses were recorded simultaneously on a Honeywell Visicorder oscillographic recorder. The general experimental arrangement is diagrammed in Figure 3. A large number of runs were made in this manner, testing first-order transfer functions whose time constants ranged from 0.1 to 5 seconds, and several second-order systems with time constants in the above range and whose damping ratio ranged from values giving a lightly damped response to values giving nearly critical damping. Frequency response curves calculated by the two computer programs were then compared with each other and with the exact frequency response as computed theoretically from the system differential equation. After many runs had been made on the analog computer, it became apparent that a control against which to compare the results obtained would be highly desirable. Testing the analog and using oscillographic recorder readout gave results which, although somewhat approximating conditions under which an actual plant test would be made, were naturally subject to random experimental and curve-reading errors and whose accuracy, in any case, rarely exceeded two significant digits. I t was decided therefore to compare conclusions already drawn from the analog results with results of similar pulse evaluations made using pulse test data derived by solving firstand second-order differential equations directly for the particular pulse shapes studied. I n this way one can calculate pulse data as accurately as one wishes, and determine readily the effects of rounding off data, of the value of of different calculation procedures, etc., on the computed frequency response without the interference of random experimental and curve-reading errors. Accordingly, the general linear first-order differential equation with constant coefficients was solved for both halfsine and triangular pulse forcing functions. I n addition, the general linear second-order differential equation with constant coefficients was solved for a half-sine forcing function. The y ( t ) (output pulse) data were then calculated from these results, using pulse widths T p = 1, 0.5, and 0.1 second, system time constants of 5, 1, 0.5, and 0.1 second, and, for the second-order response, damping ratios of 0.1, 0.5, and 0.8. These calculations were done by digital computer to eight decimal places. Computation was stopped wheny(t) had decayed below 1 in the eighth place. Thus having available the pulse response data computed accurately and at closely spaced time intervals, it was easy to introduce irregularities deliberately into the data and to compare the frequency response curves with each other and with the theoretical curves as done previously for the analog computer tests. Both magnitude ratio and phase angle plots were compared in these tests; however, magnitude ratio alone could have been used? since points of deviation, scatter, etc., were similar for the two quantities. Only magnitude plots are shown here.
I
Results and Conclusions FIRST ORDER CIRCUIT
I TO VISICORDER CHANNEL 2
Figure 3. Schematic diagram arrangement 98
2 A
-
TO VISICORDER CHANNEL 2
SECOND ORDER CIRCUIT
of analog
pulse test
I & E C PROCESS D E S I G N A N D D E V E L O P M E N T
A number of interesting and pertinent conclusions have been drawn as a result of the experimentation. Value of Using &),, as Calculated with Frequency Response as Indication of Point of Deviation. I t was felt to be extremely important to assess accurately whether or not
the initial approach of s(w)% to zero, as calculated numerically from the input pulse data, could be used as a guide to the value of frequency beyond which it would be unsafe to trust the computed frequency response curve. Theoretically, of course, the size of is the perfect indicator of deviation, because the magnitude ratio and phase can be calculated a t every frequency for which J.(U)= # 0. Practical calculations, necessarily involving the numerical approximation of integrals and thus including truncation and round-off errors, will always suffer enough loss of significant digits as approaches zero so that actual deviation occurs before the first zero of frequency content is reached. Thus the question remains: Can s(w)* be approximated sufficiently well near its first zero to act as a useful measure of deviation? Fortunately for the experimenter, the computation can be performed with sufficient accuracy, a t least when the input data are reasonably accurate. I t has been found through experience with many experimental and analytically computed pulse tests, that deviation usually sets in when s ( w ) % drops to between 0.3 and 0.2. These same limiting values were observed even when all eight significant digits were carried for the cases of computed pulse response, and whether Filon's formula or the trapezoidal rule quadrature was used. These bounding values agree roughly in many cases with Hougen and IYalsh's ( 7 7) observation that reliable Fourier transforms could be calculated for half-sine and displaced-cosine shapes out to 70 to 807, of the first zero of s(w),. Because of the varying shapes of pulse frequency content curves, however, the above limits on s(w),, have been found to be a more generally applicable quantitative criterion of deviation. Although it is sometimes possible to recover reliable frequency response information between zeros of S ( W ) ~ , the authors have found that this is the exception rather than the rule for most practical tests, and do not recommend using M.R. and 4 values calculated past the first zero. Effect of Various Data Errors. TRUNCATION OF DATA. Truncation of the output pulse is the sudden equating of the
\!
5OLlD L I N E
0
-
-
A C T U A L CURVE FOR S Y S T E M
0
C A L C U L A T E D MR. FOR CASE WHERE TRUNCATION HAS O C C U R R E D
I
2
4 W
Figure
4.
6 ' Q ' I b (RADIANS/SEC,)
Effect of truncation
+
, 20
of
,
I
40
.
, . . . 60
output
80
pulse
Test of (0.1s 11-l with 0.5-second half-sine pulse Trapezoidal rule used in d a t a reduction procedure
pulse to zero before it has actually reached this value. O n e might be tempted to alter the data in this manner in cases when the output pulse decays very slowly and data reading becomes tedious. T h e findings of the authors agree with those of Hougen and Walsh ( 7 7 ) . As reference to Figure 4 will show, the low frequency portion of the curve is made abnormally flat, while the high frequency part is attenuated too sharply. The steady-state gain (M.R. at w = 0) is reduced below the true value because of the elimination of part of the area under the product curvey(t) cos w t . EFFECTOF SUMBER OF SIGXIFICANT DIGITSO F x ( t ) ASD y ( t ) . Calculated magnitude and phase were compared both by visual comparison of the values and by plotting (plots could be made to three significant figures). It is reasonable to suppose that the usual engineering use of the frequency response data would require accuracy no greater than three figures. The calculated magnitude and phase values were not extremely sensitive to the significant digits carried in the program input data. Carrying three digits for x ( t ) and y ( t ) gave the same plot as when all eight digits of the analytically derived pulse response were used. IVhen the accuracy was reduced to two digits, somewhat more scatter was observed after deviation had set in, but the plots were essentially identical up to this point. Carrying only the first significant digit of the pulse curves had a more pronounced effect, moving the point of deviation back to a frequency value about 307, lower than for the runs using more precise data (see Figure 5 for this comparison). Since normal chart-recorder curves can be read to a t least two significant figures, one can conclude that recording equipment of the type generally found in industry is sufficiently accurate for pulse testing. INSUFFICIENT CURVEREADISG.Using too few data pointsthat is, reading the pulse curves a t excessively large time intervals-will naturally introduce error due to insufficient approximation of the curves by the points read. This will of course cause deviation of the calculated frequency response earlier than if the curves had been read more accurately. It is likely also to give improper emphasis to certain frequencies during the numerical data reduction, leading to some inaccuracies even before the onset of gross deviation. On the other hand, the use of more points than is really necessary will increase computation time and tedium of curve reading, often out of proportion to the small increase in accuracy obtained. I t is rarely of benefit to use more than 50 points for a normal exponential-type response curve. For cases of damped-oscillatory or other complex response forms, around 100 points are usually sufficient, although this \vi11 naturally depend upon the individual circumstances. POSSIBLEILL EFFECTSOF USINGDIFFEREXT INCREMESTAL TIMEVALUESFOR DIFFERENTPORTIONSOF SAME PULSE. As Hougen and 'Walsh ( 7 I ) have pointed out, when reading points from a slowly decaying curve, it is convenient to use a small At value on the first portion of the curve, where the ordinate is changing rapidly, and switch to a larger At on the sloivly decaying tail of the curve. Since this procedure reduces information fed into the pulse-reduction program, it was thought desirable to determine the effect of using a larger At for the tail of a pulse curve on the accuracy of the calculated results, since no possible ill effects have been mentioned in the literature. I t was found that use of two At values in this manner can reduce the computational accuracy somewhat (see Figure 6). This loss of accuracy is a function mainly of the quadrature procedure used. I t appears to be significant only in the case VOL. 2
NO. 2
A P R I L 1963
99
I DIGIT
CARRIED
no investigators actually used Filon’s formula for pulse reduction and compared it with simpler methods. Two data evaluation methods were used for all pulse tests: the unmodified trapezoidal rule method (Equations 17 to 20) and Filon’s formula method (Equations 24 to 31). The comparison of identical runs calculated by both methods gives insight into the effect of the oscillation of the product curves upon accuracy of the derived frequency response, and gives a basis for determining certain other advantages of Filon’s method. If enough data points are read so that the product curves are well approximated for all At values of interest, and if the pulse data are reasonably accurate (at least two significant digits, At intervals chosen and number of points read in accordance with previous suggestions) the trapezoidal rule method may be expected to be as good as any more sophisticated procedure. In such cases, deviation occurs from other sources, such as accumulation of truncation and roundoff errors, before appearance of deviation due to rapid fluctuation of the product curves a t large w . However, one cannot always be sure in advance how well his data satisfy the above criteria. I n addition, Filon’s formula is in many cases less sensitive to the presence of random errors in the pulse data. Filon’s method is rather insensitive to increasing At on the tail of a pulse, whereas the trapezoidal rule calculation shows some sensitivity to this. Finally, Filon’s formula invariably shows a n advantage over the trapezoidal rule when pulse data of somewhat limited accuracy are compared. Figure 7 demonstrates how Filon’s formula allows calculation closer to the theoretical curve in this latter case, even after the onset of deviation. I t may thus be concluded that it is best to use Filon’s formula for the data reduction. This method, applying a curvedsegment approximation to the pulse curves and being insensitive to large w, should prove in the long run superior to the methods other workers have used for the reduction of data. The increase in computing time over simpler methods, although noticeable on a small computer such as the IBM-650 used by the authors, would be of much less importance on larger machines. Influence of I n p u t Pulse Frequency Content on Reliability and Range of Results. This point is of great significance both in the performance of the actual test and in the proper use of the reduced data. ~
I
I
2
4
6 8 K )
20
40
6 0 8 0 )O
Figure 5. Effect of data accuracy on calculated frequency response Test of (0.1s f 1 )-I with 0.5-second half-sine pulse Trapezoidal rule used in data reduction procedure
SOLID L I N E -THEORETICAL
1
.3
CURVE
TWICE THE A t - v a ~ uAND ~ HALF A S MANY P O I N T S WERE R E A D O N T H E T A I L OF T H E y ( t ) P U L S E FOR THE TRIANGULAR POINTS AS W E R E R E A D FOR THE CIRCULAR POINTS,
4
6
W
8
10
20
40
60 80
IO
cRADlANS/SEC.l
Figure 6. Effect of doubling At on last portion of system output pulse Test of (0.1s f 1 )-I with 0.5-second half-sine pulse Trapezoidal rule used in data reduction procedure
of the trapezoidal rule approximation to the product curves themselves, for this method of calculation, representing a fixed sampling of these curves, is sensitive to increased oscillation at higher w values, whereas the trapezoidally approximated Fourier transform and Filon’s formula methods are not. The results of a study of this problem may be generalized by sa\,ing that any loss of accuracy is likely to be smallest in cases where large total numbers of points are read from the ciirves and where Filon’s formula is used for data reduction. Errors caused by using a larger At for the tail of a curve are rarely large in any case, and fear of such errors should not in general deter an investigator from taking judicious advantage of using different At values. Comparison of Filon’s Formula to Trapezoidal Rule Procedure. Schumacher (76) and Smith and Triplett (20) suggested the Filon approach, but prior to the previous study 100
I & E C P R O C E S S DESIGN A N D DEVELOPMENT
8
... u
(01s t IT’
SYSTEM
i
COMPUTED
e
4.
A
FROM P U L S E RESPONSE
- T R A P E Z O I D A L RULE USED - FILON’S FORMULA USED
SOLID
LINE
-THEORETICAL 2
4
CURVE 6
w
810
20
40
60 80 I 0
( RADIANS/SEC.)
Figure 7. Comparison of Filon’s formula and trapezoidal rule for pulse data of limited accuracy Two significant figures carried in pulse data
0
... P a w
3 c)
t 0 4 L
0.65 SEC. RECTANGULAR
PULSE
0
0.65 SEC. H A L F - S I N E
A
0.65 SEC. T R I A N G U L A R
PULSE
SOLID L I N E : T H E O R E T I C A L
CURVE
PULSE
F I R S T Z E R O E S OF
2
4
. A
J
i,!
6 8 1 0 W
Figure 8.
I
S(w):
,
20
(RADIANS/SEC.l
Comparison of pulse tests on
,
,
40
,
,
60
,
80
+ 1 I-'
(0.5s
Analog simulation Trapezoidal rule used in data reduction procedure
Reference to Figure 1 will show the comparison of frequency content, s ( w ) n, for various pulses of mathematically describable shape. Frequency content information for pulses of random shape was also available through the numerical ~ same time that the frequency recalculation of ~ ( 0 a t) the sponse was evaluated by the computer. Magnitude ratio curves derived from rectangular, half-sine, and triangular pulse inputs to the system (0.5s l)-l are compared in Figure 8. Although a large number of such comparisons were made, Figure 8 illustrates the general trend of all these results. T h e calculated points deviate first for the test using the rectangular pulse; this is to be expected, since its frequency content is lowest (frequency value a t first zero is lowest of all). Deviation of the half-sin? test occurs next, with the points for the triangular pulse test reproducing the true curve to a slightly higher value of frequency. (Reference to the frequency content curve, Figure 1, will aid in visualizing the situation.) The effect of frequency content as predicted by theory is borne out experimentally in every case studied: Frequency content of the input pulse has a great influence on the maximum value of frequency a t which the frequency response is reproduced accurately. The use of pulses having higher frequency content is thus recommended whenever possible.
+
RECTANGULAR
HALF-SINE
TRIANGULAR
DISPLACED COSINE
t WEIGHTED
SQUARED
CUBED
DISPLACED
TRIANGULAR
TRIANGULAR
IMPULSE
COSINE
Figure
9. Some standard pulse shapes
Arranged in order o f increasing frequency content
An intensive study of the effects of varying pulse shape and width, including the effects of varying the number of pulse points read and the effect of output pulse truncation, was recently reported by Dreifke (8). This study was accomplished on a n analog computer, and used an equivalent form of Equation 22 for all data reduction. Some quantitative criteria are offered in this work for maximum recovery of frequency response, and the reader desiring to use the method of Equation 22 should refer to this thesis as well as the present paper. I t is probably safe to say that an optimum pulse shape for a given process does not exist. Perhaps a more accurate statement is that such a large number of possible pulse shapes will produce the same end result that it is meaningless to bothvr with the notion of an optimum shape. Kevertheless. certain broad rules do apply. The pulse must contain a wide enough band of frequencies so that the system dynamics are excited over the entire range of frequency values in which one is interested. For a constant pulse shape, frequency content increases as pulse width, T,, decreases. But one cannot expect to apply a quick pulse to a sluggish system and obtain an output that is large enough to be accurately measurable. Attempts to raise the output above the noise level by increasing input magnitude may saturate the system or force it so far beyond its normal operating point as to cause nonlinearities to reveal themselves. Fortunately, one can control frequency content in another way. The frequency content of a pulse is profoundly affected by its shape, which ranges from the rectangular pulse a t the low extreme of frequency content to the impulse function at the high end (Figure 9). The frequency content drops off as the square-cornered, flat-in-the-middle rectangular form is approached. Going in the opposite direction toward forms increasingly rounded and then approaching a central spike, the frequency content increases greatly, until the limiting form of the pure Dirac impulse, having a constant frequency spectrum, is reached (77). The impulse function. having unit area, infinite height, and infinitesimal \vidth, is of course experimentally unattainable, But forms such as the higher degree weighted displaced cosine pulses (8) and the squared and cubed triangular pulses have great frequency content for a given pulse width and can be generated by analog methods. I n actual plant tests it is usually inconvenient to use these regular shapes, since it is much easier to generate a pulse by, say, turning a valve by hand. Although frequency contents of such random pulse shapes vary nidely as a function of the actual shape, it is usually possible to produce a pulse having a sufficient content of frequency to cover the range of interest, if the pulse Midth is approximately equal to the dominant process time constant. If this time constant is unknown. it will be necessary to make several pulse tests using pulses of different widths. I n any event the actual harmonic content that has been produced will be available for examination if, as previously recommended, S ( W ) is calculated numerically at the time the data reduction is accomplished. A small amount of experience gained either through actual tests on plant equipment or on an analog computer mill quickly lead to enough familiarity n i t h the method so that little trouble ill be experienced.
Additional Study Needed
Several topics concerned with pulse testing are of sufficient interest and importance to warrant further study. Is it possible to develop a class of functions cvhich possess VOL. 2
NO. 2 WPRIL 1 9 6 3
101
= time-dependent system input = time-dependent system output
maximum frequency content subject to certain restrictions such as existence for a finite time, T,, and suitability for physical realization in pulse tests? Such a problem might be approached using the calculus of variations, and would aid in selection of forms likely to produce the most reliable results. Pulses which consist of more than one peak, such as multiple sine pulses, often exhibit a large maximum in frequency content over a narrow range of frequencies [Figure 4 of ( Z O ) ] . I t has been suggested that several such pulses, properly chosen, could test a system over a large range of frequency with unusually reliable results (20). There needs to be a detailed study of this frequency-band-selection technique made on both simulated and actual systems to supplement what has been done ( 8 )and to establish criteria for its efficient use. A number of pulses exist whose frequency content curves drop off steadily but never reach zero-the ramp and squared ramp pulses, and the squared triangular pulse, to mention but a few (see Figure 1). Use of pulses of lower frequency content, whose s ( w ) curves exhibit zeros, was considered wisest in this study, because the frequency content of this class of pulses resembles more closely the content of random shape pulse!; produced in actual on-the-job tests. Nevertheless, use of the analytical response of simple systems to the special pulses having no zeros of s ( w ) in a study such as the one outlined might prove useful in investigating the maximum possible practical recovery of frequency response from pulse testing. No such study has yet been reported. The statistical design of experiments is finding increased application in engineering analyses. A series of experiments designed and analyzed statistically might shed new light on some of the topics discussed in this paper.
x(t) ~ ( t )
Nomenclature
n. SOC.Mech. Eng. 80,833 ( (13) Lees, S., Hougen, J. O., Znd. Eng. Chem. 48, 1064 (1956). (14) Ley! B. J., Lutz, S. G., Rehberg, C. F., “Linear Circuit Analysis, McGraw-Hill, New York, 1959. (15) Ross, C. W., Gardner, M. T., “Comparison of Sinusoidal Forcing Function and Pulse Forcing Function Techniques for Obtaining Dynamic Performance Chai acteristics of Fire Control Systems,” Instrumentation Laboratory, Mass. Inst. of Technology, Cambridge, Mass., Rept. T-30 (May 1953). (16) Schumacher, L. E., “Methods of Analyzing Transient Flight Data to Obtain Aircraft Frequency Response,” Wright Air Development Center, Memo. Rept. MCRFT-2268 (January 1950). (17) Scott, R. E., “Linear Circuits, Part 2. Frequency-Domain Analysis,” Addison-Wesley, Reading, Mass., 1960. (18) Seamans, R. C., Jr., Blasingame, B. P., Clementson, G. C., J . Aeronaut. Scz. 17, 22 (1950). (19) Seamans, R. C., Jx., Bromberg, G. G., Payne, L. E., Zbid., 15, 535 (1948). Mech. Eng. (20) Smith, G. A., Triplett, M‘. C., Trans. Am. SOC. 7 6 , 1383 (1954). (21) Stewart, W. S., Sliepcevich, C. M., Puckett, T. H., Chem. Eng. Progr. Symp. Ser. 57, No. 36, 119 (1961). (22) Watjen, J. SY., “Dynamic Behavior of a Pulsed Plate Extraction Column,” Ph.D. thesis, University of Virginia, 1962.
f(t)
F[ ] G(jw) G(s)
=
arbitrary function of time
= Fourier transformation
= system performance function = system fransfer function
arbitrary indices
,/=i I
-
magnitude ratio, units o: ~ ( t divided ) by units of x ( t ) Laplace transform complex variable sum of even integrands in Filon’s formula sum of odd integrands in Filon’s formula Fourier transform of a pulse function frequency content normalized to 1 a t w = 0 frequency content of input and output pulses, respectively time integrand in Filon’s formula imaginary part of system performance function real part of system performance function pulse width, units of time width of system input pulse, units of time width of system output pulse, units of time
102
l&EC PROCESS DESIGN A N D DEVELOPMENT
CY,
0
/3, y, e = constants in Filon’s formula
$it)
w
= phase angle, degrees = arbitrary function of time = angular frequency, radians per unit of time
Acknowledgment
The authors recognize the aid of the National Science Foundation, whose support of the Vanderbilt Computer Center through Grant No. M F - G l 0 0 8 has made possible the many computations involved in this study. literature Cited
(1) Aseltine, J. M., “Transform Method in Linear System Analysis,” McGraw-Hill, New York, 1958. (2) Caldwell, W. I., Coon, G. A,, Zoss, L. M., “Frequency Response for Process Control,” McGraw-Hill, New York, 1959. (3) Chestnut, H., Mayer, R. W., “Servomechanisms and Regulating System Design,’‘ Vol. I, Wiley, New York, 1959. (4) Clements, W. C., Jr., “Analog Simulation of a Concentric Pipe Heat Exchanger Using Frequency Response Techniques,” M.S. thesis, Vanderbilt University, 1961. (5) Cohen, W. C., Johnson, E. F., 2nd. Eng. Chem. 48, 1031 (1956). (6) Coppedge, T. N., Jr., TVolf, R. L., “Pulse Testing of Systems with Some Nonlinearities,” Instrumentation Laboratory, Mass. Inst. of Technology, Cambridge, Mass., Rept. T-68 (May 1954). (7) Draper, C. S., McKay, W., Lees, S., “Instrument Engineering,” Vol. 11, McGraw-Hill, New York, 1953. (8) . , Dreifke, G. E., “Effects of Input Pulse Shape and Width on Accuracy of Dynamic System Analysis from Experimental Pulse Data,” Sc.D. thesis, it‘ashington University, 1961. (9) Filon. L. N. G.. Proc. Rov. Soc. Edznburph 49. 38 (1928). , , Noise,” McGrawlHill
RECEIVED for review July 16, 1962 ACCEPTED December 12. 1962