Numerical simulation of an oscillating towed weight

Ocean Engng, Vol. 16, No 2, pp. 143--172, 1989. Printed in Great Britain.

0029-8018/89 $3.00 + .00 Pergamon Press pie

NUMERICAL SIMULATION OF AN OSCILLATING T O W E D WEIGHT THOMAS N. DELMER Advanced Computer Solutions, 5580 La Jolla Blvdl Suite 449, La Jolla, CA 92037, U.S.A.

and THOMAS C. STEPHENS* Horizons Technology, Inc., 7830 Clairemont Mesa Blvd, San Diego, CA 92111, U.S.A. Abstract~The TAAS (Towed Acoustic Array Simulation) computer program, previously used to calculate the large-scale, three dimensional dynamic behavior of series-connected marine structures, is used here to simulate the reaction of a cable-towed weight to small-scale transverse sinusoidal oscillations superimposed upon steady towing. The reaction shows a linear component at the driving frequency and a quadratic, nonlinear component which produces a response at twice the driving frequency as well as a systematic distortion. The amplitude and phase of the three components of the response are estimated at points distributed throughout the system. The results compare favorably with those from D.A. Chapman's linear analysis based on mathematical perturbation techniques. Our extension of Chapman's method to estimate the response at twice the driving frequency also compared well. We conclude that full simulation is an effective alternative to linear analysis; it also illuminates nonlinear behavior beyond the reach of traditional linear methods. We find that the nonlinear response of oscillating towed systems shows higher propagation speeds and less damping than the linear response and, consequently, may be significant. 1.

INTRODUCTION

THE T A A S COMPUTERp r o g r a m ( T o w e d A c o u s t i c A r r a y Simulation) has b e e n successfully used to calculate the large-scale, three dimensional d y n a m i c b e h a v i o r of o c e a n cable systems ( D e l m e r , Stephens and C o e , 1983, referred to as D S C ) . While primarily a tool for analysing towed arrays, it is useful in general for calculating the reactions of seriesc o n n e c t e d marine structures (moorings, t o w e d cables, etc.) e x p o s e d to highly dynamic situations ( D e l m e r , Stephens and Tremills, 1988, referred to as D S T ) . In this paper, we present the application o f T A A S to the p r o b l e m o f analysing oscillating towed systems. W e will begin by characterizing the systems bf interest, followed by a discussion of linear analysis techniques applicable to those systems, describe s o m e prior analyses and outline the T A A S simulation a p p r o a c h . 1.1. Oscillating towed systems T h e subject of oscillating t o w e d systems has been of interest for s o m e time in the towed array c o m m u n i t y , s p a r k e d by the observation that the long, thin, neutrally *To whom correspondence should be addressed, 143

144

T . N . DELMER and T. C. STEPHENS

buoyant arrays tend sometimes to meander in a snake-like fashion while under tow. The degree to which this behavior is due to erratic towpoint motions filtering down to the array has become a topic of investigation. Traditionally, such situations have been analysed mathematically by perturbation techniques, which characterize small-scale phenomena as a deviation from some norm. Perturbation techniques are in essence special calculations, embodying certain limitations. We wanted to see to what extent a direct simulation could duplicate, and perhaps extend, the results of perturbation calculations. The systems of interest are considered to be under steady tow (the norm), with smallscale oscillations (the deviation) superimposed on the towpoint motion. By "small" oscillation we mean small relative to the scale of the object. The systems we will discuss are on the order of a kilometer scale, and we are concerned with oscillation amplitudes (peak to trough) of 10-20 m (about 1% of the object scale size). The oscillation periods range from 1 to 10 min. The form of the oscillation, for simplicity's sake, will be sine waves--a reasonable choice since any periodic phenomena can be decomposed into sines and cosines according to familiar Fourier analysis. We are interested in the system response to the towpoint oscillation. The driving oscillation will create a disturbance that propagates in an attenuated form along the system. A fair presumption might be that the frequency response of the system would follow that of the towpoint; and that the amplitude of response would have a simple, linear relationship with the amplitude of the towpoint for a uniform system. Doubling the towpoint amplitude, for instance, should produce a doubling of the response amplitude interior to the system. Both properties are characteristics of linear systems. But we are also interested in any nonlinear effects that may turn up. These would involve, for instance, amplitude responses that scale to some power of the driver amplitude or frequency doubling effects. Also interesting would be any systematic distortions, i.e. migrations of the average system positions away from the initial steadystate. 1.2. Linear system analysis For those systems that are expected to respond linearly there exist a well-developed set of tools designed to analyse frequency response, grouped under the name "Linear System Analysis". Even though these tools are designed specifically for so-called linear systems, it is also possible in certain circumstances to apply them to nonlinear systems. This is because some nonlinear systems behave linearly in the neighbourhood of steady state solutions. The systems that we are interested in analysing, the towed systems, are actually highly nonlinear; but if we are careful to analyse oscillations about steady state solutions--steady towing--then linear analysis should apply. The two most common categories of linear system analysis are termed single-frequency response and impulse

response. The simpler of the two categories is single-frequency analysis. In essence, it consists merely of applying a continuous oscillation that drives a system at a single-frequency and at a constant amplitude. It may take some time for the system to react to the onset of oscillations, but once transient effects die out, the system will eventually reach what could be called oscillatory steady state. At that point it is relatively simple to estimate the response amplitude at any point downstream of the driving oscillation. The relationship

Simulation of an oscillating towed weight

145

between the oscillation amplitude and phase at any point on the cable to that of the driving amplitude and phase (called the transfer function) is fairly easy to discern. One can estimate the transfer function by measuring snapshots of the oscillation. Singlefrequency analysis has the advantage of being intuitive in that the relationship between the driver and the response is so simple. But it is not very comprehensive because one must perform a series of experiments on the system involving various frequencies in order to be able to draw general conclusions regarding the system response. The other method is called impulse response analysis. This involves driving the system as hard as possible with a sharp,'spike-like impulse which should be of as short duration and as large amplitude as possible. This is akin to hammering a bell and listening to the ringing response. The theory is that an ideal impulse (actually a mathematical delta function of zero time-width and infinite magnitude) is composed of all frequencies, equally and simultaneously. An impulse response analysis is comparable to performing all possible single-frequency analyses in a single shot. The problem with this method is that the system response to such an impulse may appear very complex; there may be no obvious relationship between the form of the impulse and the form of the response. Relatively sophisticated digital signal processing (DSP) techniques are required to extract the frequency contents of the system response. In reality, perfect impulses do not exist. Laboratory impulses always have a finite time duration and a finite magnitude. So the same DSP techniques have to be applied to the impulse as well as the response. (The transfer function is calculated as the Fourier transform of the response divided by the Fourier transform of the driver.) Impulse response methods are appropriate for true linear systems, but the danger when applied to nonlinear systems, even those near steady state, is that they will excite nonlinear responses and possibly jar the system off of its steady state. The method has the, advantage that it is comprehensive and economical in that one can draw general conclusions covering the entire range of frequency response of the system from the results of a single experiment. It gives, in theory, the same result as the single-frequency method for all those frequencies where the driver and sampled output yield signals above noise (either measurement noise or noise from the numerical calculation). But it is not very intuitive due to the complex responses possible and the required use of DSP techniques. 1.3. Some prior analyses Two particular examples of linear system analysis applied to towed systems are mentioned here. The first is due to Kennedy and Galletta (1981) and Kennedy and Strahlen (1981, two related papers referred to collectively as KGS) at the Naval Underwater Systems Center. KGS performed full scale ocean experiments involving a ship-towed array system. They devised several types of towing schemes to drive an oscillating response. The schemes have to be viewed as hybrids of the impulse and single-frequency analyses. One method could be termed a truncated single-frequency where the ship attempted to move in a pure oscillatory fashion at a constant frequency and amplitude. The ship motion was practically limited to just a few oscillation cycles though, and so transient effects were bound to contaminate the analysis. Another method was a random towpoint oscillation: the ship meandered more or less about a straight line. Finally, the ship attempted an impulse of sorts trying to maneuver as sharply as possible, in a square wave-like fashion. Square waves can serve duty as

146

T.N. DELMERand T. C. STEPHENS

impulse functions because they can be regarded as the integral of several alternating delta functions. In all three of these experiments the drivers and the responses were sufficiently complicated that DSP techniques were applied to extract the frequency components. They found that the results of the truncated single-frequency and the meandering path experiments produced more or less consistent results and were able to calculate transfer functions. But the results of the square wave maneuver defied analysis. They conjectured that this maneuver drove the system nonlinear, and so the analysis was not valid. The second example is due to Chapman (1984b, referred to here simply as Chapman) at the University of Bath. He studied the motion of a towed body caused by transverse, single-frequency sinusoidal motion of the towing ship. He performed a mathematical analysis of the differential equations describing the towed system. By assuming the sinusoidal motion small compared to the towing motion, he was able to employ perturbation analysis about the steady towing solution. He used a computer to calculate the linear response amplitude, but this was not a simulation in the sense that we will use the term. He neglected phase responses as well as nonlinear effects. Chapman's analysis is important to us because he inspired the form of our numerical experiments, and we use his results as an independent check on the TAAS results. Two papers by Dowling (1988a,b) appeared after the work described here was completed. He described the response of a towed array to an oscillation at the tow point as a series expansion, i.e. analytically. It may be possible to use his results to verify TAAS for the simulation of a towed array. 1.4. TAAS simulation approach TAAS has been used for several years to calculate the large scale motion of towed array systems. We were curious whether TAAS could also be used to simulate small scale oscillations superimposed upon steady towing motion. Accordingly, we modelled a situation similar to Chapman's, calculating the phase and amplitude of the simulated response of a towed weight to an oscillation of the towpoint. We chose to perform a strict single-frequency analysis for the sake of simplicity. KGS had encountered difficulties with the impulse approach in their experiments, so we avoided it in the simulation. In either case, some analysis of the results is necessary to estimate the amplitude and phase, and the single-frequency procedure is much simpler. The ratio of the amplitude at the towed weight to the amplitude at the towpoint, the magnitude of the transfer function, is given by Chapman. Our objective was to attempt to duplicate Chapman's results for the simple towed weight problem and then to move on to a full towed array case. The full simulation, once validated, has an advantage over perturbation approaches like Chapman's: it does not rely on the solution being as expected. The perturbation approach assumes that the main motion of the towpoint can be separated from the small transverse motions. When the assumption is correct, the perturbation approach gives correct results. However, the towed system is nonlinear and perturbation techniques are not designed to handle nonlinearities. In essence, the perturbation approach assumes a form of the solution that may not be valid. The rest of this report consists of a chapter detailing the methods used to simulate an oscillating towed weight and then analyse the output data; a chapter summarizing and discussing the linear and nonlinear results of the experiment; and the conclusions.

Simulation of an oscillating towed weight

147

2. METHODS This section contains descriptions of the methods used to simulate an oscillating towed weight, and then to estimate the amplitude and phase responses. The discussion of the simulation focuses first on definition of the system configuration and scenario, proceeds to consider issues related to model discretization, and ends with descriptions of the two types of steady state calculated. Under the heading "Data Analysis", we present the form of the simulation outputs to be analysed, discuss linear system analyses in general (and impulse vs single-frequency analysis in particular), and finally formulate the amplitude/phase estimation procedures.

2.1. Simulating an oscillating towed weight The response of the towed weight was determined by numerical simulations, or computer experiments, using the computer code TAAS. Detailed descriptions of TAAS models and equations of motion are presented in DSC and DST and are not repeated here. The results of Chapman and of KGS influenced the design of our experiment. The simulation was performed with input parameters that coincided with points in Chapman's tables for bare cables, i.e. cables without fairings. This enabled direct comparison with his tables for the linear part of the amplitude response. The computer code that Chapman published was extended slightly to calculate nonlinear effects and phase for further comparisons. 2.1.1. System configuration and scenario. Figure 1 is a perspective diagram of the towed weight. The towpoint is fixed at the origin of coordinates, and a steady, uniform current flows in the positive x direction. The cable is strung out in a steady state configuration in the x-z plane. The use of a current instead of a forward velocity for the towpoint simplifies the interpretation: velocities are initially zero, and we avoid x coordinate magnitudes that steadily grow, which might effect the error control over the course of the problem. The towpoint oscillation occurs along the y axis, in the direction perpendicular to the plane of the steady state configuration. The form of the oscillation is a pure sine wave. The frequencies were chosen not only to coincide with Chapman's tables, but were also similar to those used by KGS in their experiments, (with one relatively high

Oscillation

J

,

FIG. 1. Towed weight scenario (perspective).

148

T . N . DELMER and T. C. STEPHENS

frequency used to push the simulation). The amplitude of 10 m was less than the 30m amplitude used by KGS. This method of chosing the parameters ensures that they are, with the possible exception of the highest frequency, realistic in the sense that they are attainable in experiment. The towed weight series was comprised of five cases: four had amplitudes of 10 m with periods that range from slightly over 1 min to over 10 min (66.5, 133, 266, 665 sec). The fifth case doubled the amplitude of the 133 sec case to test for the presence of nonlinear effects. Chapman published sufficient details of his towed weight model to allow us to create our own, similar model. (KGS did not give the parameters necessary to repeat their experiment, so their work is not useful for checking a simulation.) Figure 2 is a schematic diagram of the computer model of the towed weight. This system is one kilometer long with a uniform steel tow cable and a 3400-kg weight hanging on the end. The cable itself is three centimeters in diameter and 1.02 kg per m linear density, immersed. The drag coefficient of the cable is 1.07 with a tangential loading function of 0.001 to approximate the value of zero as used by Chapman. The normal loading function is equal to sin20, where 0 is the angle of the cable to the current, to agree with Chapman. No drag properties were assigned to the weight at the end of the cable. In contrast to our simulation, Chapman's calculation ignored inertial effects. This assumption is appropriate for the cable but may not be for the weight. Appendix B discusses the inertial effects of the weight and shows that they should not adversely affect our comparison. Towpoint

Weight, 3400kg

15 Elements: 66.7m..03m dia., 1.02kg/rn

\ L

I I I I t I I L [ I I L I L J 1000m FiG. 2. Towed weight model schematic.

2.1.2. Discrete cable model. As in most finite difference calculations involving wave propagation, the analyst must decide how finely to sample (or zone) the transmission medium (here, the cable) to ensure that waves with the expected range of frequencies propagate correctly. The choice is crucial; waves of a particular frequency simply cannot form if too few reference points (nodes) exist per wavelength. If the zoning is too course, the simulation will fail in the sense that amplitudes and phases will be wrong. Conversely, for a given set of nodes, the simulation is guaranteed to fail if pushed beyond a certain frequency. As a modeling choice we decided it was sufficient to subdivide the kilometer-long cable into 15 elements each 66.7 m long. The basis for this choice was a desire to have five nodes per wavelength at the shortest period (133 sec) for which the simulation was to be accurate; it was intended that the simulation would start to fail at the shorter period, 66.5 sec. The wave was expected to move at about the towing speed, 2.5 m per sec, so the wavelength expected was approximately 333 m (period times tow speed), i.e. three wavelengths per km. This decision was made before the simulations were

Simulation of an oscillatingtowed weight

149

performed and did not take into account either the fact that the cable does not stretch its full length downstream, or the complexity of the waves that exist on the cable. Towed systems support a variety of transverse and longitudinal waveforms, but some subtle points must be noted: • For transverse cable motion, the cable inertia is dominated by the normal fluid drag. Transverse cable oscillation modes thus resemble diffusion in a moving medium more than they do waves. In fact, for immersed cables transverse "waves" are, strictly speaking, fictitious---they persist only so long as they are actively driven; otherwise they simply diffuse away. The simplistic argument above regarding requirements for waveform resolution is then deceptive: there still must be several reference points per "wave-length", but the relation between frequency and wavelength for driven diffusion is different from that of classical wave propagation. The effective wavelengths are longer, and thus fewer nodes per unit cable length are required. The present simulation is thus valid at much higher frequencies than we had planned due to the diffusion. • For longitudinal cable motion, the cable inertia dominates the low tangential drag. Longitudinal cable modes are bonafide waves, but they propagate at near the speed of sound in the steel cable. The transmission speeds are so fast that the propagation times are practically negligible so far as this simulation is concerned. Their high speed means that the waves will travel large distances in a single cycle so that there will be many elements per wavelength, enough to model the waves well. If any doubt arises about the effects of element size on the oscillation results, the remedy is to rezone the cable and rerun the problem, noting any differences. When T A A S produces consistent results at finer and finer resolutions, the problem is considered solved.

2.1.3. Steady state calculations. TAAS calculates the initial steady state configuration by the "dynamic relaxation" method: the system is initialized rather arbitrarily and allowed to dynamically settle in the steady current, until the velocities are negligible. TAAS program inputs expressing the exact initiation are listed in Appendix A. Figure 3 is a plot of the profile of the towed weight system at steady state, subjected to a 2.5 m per sec uniform current. This plot is drawn to scale. The towpoint oscillation will occur in a direction perpendicular to the plane of this figure. Each of the dots in the figure represent one of the 16 modeled nodes. At time zero, the reference time, the towpoint commences sinusoidal movements in the y direction, in the horizontal plane, perpendi~cular to the flow. Inputs for the TAAS program defining this phase are also listed in Appendix A. Figure 4 shows sample response histories for all three coordinates at the weight. The case shown is the high frequency case, period 66.5 sec. Note the initial delay in the y response of the weight. This is the propagation time required for the disturbance to reach the weight. There appears to be a large, positive initial displacement in y, a low frequency component before oscillation takes hold. This decays gradually over the period shown in the figure, after which the system appears to be in a steady state oscillation. The transient phase lasts about 600 sec. In order to compare x, y and z oscillations on the same scale, we had to "unbias" the large x and z coordinate values: We subtracted from x and z their

150

T.N.

-I00

DELMER a n d T. C . STEPHEMS

~ e %o ~e~e ~e

\e \

[,q

\ -400

\ 2.5 m / s

Current --)

\ •\

\ -700

, •

' 400

0

800

Displacement X (rn) FIG. 3. Steady state cable profile.

0.6 Y

v

0.0

VvVv

t~ t'~

0

600

300

900

Time (sec) FIG, 4. X Y Z

displacement time-series at weight, period=66.5 sec.

respective averages over the last few periods during steady state oscillation. (y oscillates naturally about a zero coordinate.) The inset figure is an enlargement of the x, y and z responses during one period of the driver oscillation, near the end of the time history. 2.2. Data analysis

2.2.1. Form of the cable response. One could directly measure on Fig. 4 the period of the y oscillation at the weight. It turns out to be the same as the period of the towpoint oscillation--a linear response. The amplitude of the y response could also be measured on this plot and related to the amplitude of the towpoint oscillation. One could construct the general transfer function from plots such as this. This is why singlefrequency response experiments are attractive.

Simulation of an oscillating towed weight

151

It is evident that the y towpoint oscillation, confined to a line perpendicular to x and z, nevertheless produced x and z reactions. The frequency of the x - z response appears to be twice that of the driver oscillation--a nonlinear effect. The source of these nonlinear x - z oscillations will be explained below. It is also apparent that the x and z oscillations begin displaced somewhat from their final average of zero. There is a net increase in the average value for z of over 0.6 m and about 0.2 m for x. So the x - z planar position of the weight migrates to a point slightly displaced from its initial steady state position. The weight is just one of 15 free nodes modeled on this system. All the other free nodes show similar dynamic responses to this one, and they also show similar net shifts. The average configuration thus distorts slightly as the result of the presence of oscillations. For comparison purposes, Fig. 5 is a plot similar to Fig. 4 but for a different case (133 see--twice the period of the previous plot). The scale in this plot is precisely the same as in Fig. 4 to dramatize the differences. It is obvious that there are half the number of response cycles present on this graph due to doubling the driver oscillation period, and that the x and z responses once again show doubled frequencies. The x and z amplitudes appear to be about the same as for the prior case, but the y response is so much larger that it is off the scale near the beginning of the time series. The lower driver frequency has had a rather large effect on the linear (y) amplitude response at the weight, but its effect on x and z was not so noticeable. There are five separate cases in the towed weight series and each case involves 15 free nodes. The result is 75 plots that resemble Fig. 4. To reduce this large quantity of data, we standardized and automated the process of estimating amplitude and phase response, which is described next. 2.2.2. Linear s y s t e m analyses. Since the oscillation amplitude is small, the linear system response should be comparable with Chapman's linear analysis. For a linear system, there are two methods of analysis that are commonly used: the impulse response and the single-frequency method. The impulse method has some economies compared 0.6

e,

0.0

z

Period = 133 sec

-0.6

0

i

i

3OO

600

900

Time (sec) FIG. 5. X Y Z displacement time-series at weight, period=133 sec.

152

T.N. DELMERand T. C. STEPHENS

with the single-frequency method when used with a simulation. Since we chose singlefrequency, we give a brief justification h e r e preceded by a review of the methods. The salient property of a linear system is that the output v(t) produced by a weighted sum of inputs i(t) = ~Otkik(t), is the sum of the weighted outputs, v~(t), for each input, ik(t), taken separately. The etk are arbitrary weighting coefficients. So v(t) = ~,etkVk(t) for arbitrary input. Given that this mathematical relation holds for any input, one can insert an impulse function into the system at time ~-. This input is the Dirac delta function, ~(t - "r); the output, called the impulse response, h(t,-r) is measured. Using the summation property, but writing it as an integral rather than a sum, and using a weighting function i('r),, one finds that the input is f ~i(t - "r)i('r)d'r and is equal to i(t). The output corresponding to this input is the weighted sum of the individual outputs, i.e.

: v(t) = J h(t,.r)i('r)d'r . In practice there is an upper limit on the value of the input; above this limit the system distorts or simply fails. In addition, there is noise on both the input and the output so that the expression has a lower bound on its applicability. If the system is time invariant, then h(t,'r) becomes h(t - .r) and the absolute value of the time is not important. If the system is causal, there is no output before there is input so h(t,.r) = 0 for t < "r. (The variable t could be replaced by some variable other than time in which case causality would not be so natural a property.) Combining these properties gives the following relation for a linear time-invariant causal system

v(t) = It h(t - ~)i('r)d'r. The integral is called the convolution of h with i, and the notation is v = h • i. Systems that are nonlinear can behave in a linear fashion in some situations. For instance, if the output of the system is some known function of the input, e.g.

v(t) = f f [t,'r,i('r)] d'r then for i small

i

v(t) -~ f(t,'r,O) dx +

fo:

i(x) Oii (t"r'i)li=" d'r

and the deviation Av(t) = v(t) - f f(t,'t,O) d"r is linear in i(t); the form of the relation is the same as for the linear system. For the linear system, f(t,-r,O) = O. However, it is possible that OflOi[i=o = O, which can happen for a number of reasons. An example is a drag function varying as the square of the

Simulation of an oscillating towed weight

153

velocity. If the velocity is initially zero, the derivative is zero so the first nonzero term in the expansion is 0.5izO2f/Mzli= o. In such cases, the system is nonlinear even for small amplitudes. 2.2.3. Impulse response analysis. If the system is linear or if only Av(t) is of concern, it is tempting to consider an input i(t) = ioS(t) to the system where ~(t) is the Dirac delta function. In this case the output is

v(t) = I'~ h(t - "r)ioS(x)d'r = ioh(t) so that by measuring the response to an impulse or delta function input, which is just h(t) (called the impulse response of the system), one can calculate the response to any input i(t) by reinserting h(O and the input into the integral for v(t) and carrying out the integration. The impulse response contains all the information relating the input to the output for a linear system (even if the system is time varying and not causal). If the system is time invariant, the above shows that the output is the convolution of the input with the impulse function. In this case the Laplace transform of the integral yields V(s) = H(s)l(s). Upper case letters are used to denote transforms, e.g.

G(s)

=

e "st

g(t) dt

)

where the assumption is that the input is zero for t < 0. For s = jto, this becomes the Fourier transform if the integral exists. It is possible to invert the Laplace transform; given H(s), h(t) can be calculated. Since h(t) totally describes the system, so does H(s). In fact, H(jto), called the transfer function, totally describes the system and is of central interest for the towed weight problem, Here i(t) is the towpoint oscillation and v(t) is the measured response at the weight, or some point on the system. Once the weight is in a steady oscillatory motion, the amplitude ratio of the weight oscillation to the towpoint oscillation is the magnitude of H(jto) where to is the circular frequency associated with the frequency of oscillation, f, of the towpoint; to = 2"rrf. The phase difference between the oscillation at the weight and the oscillation at the towpoint is the phase of H(jto). If one uses a delta function for the transverse motion at the towpoint, the response is h(t) and one simulation can totally characterize the system, which makes such an experiment tempting. In this case, it is inefficient to do a separate simulation at every frequency of interest. However, to achieve the efficiency, the input must be a delta function. "8(x) is not a function according to the usual mathematical definition of a function," (Dirac, 1947). It is characterized by the relation

I~ f(x)~(x)dx = f(O) for any continuous function f(x). One idea of a procedure to model such a function is to consider a step of width ¢ and height 1/e centered at the origin, where ~ is considered to be arbitrarily small. The transform of the delta function is unity so, for such an input the transform of the input is l(s) = 1, and the transform of the system output is H(s).

154

T . N . DELMER and T. C. STEPHENS

The delta function is a function that has all frequency components represented equally. In fact, since a general system could have substantial response at any frequency, a function to sample a general system must have substantial content at all frequencies. It is impossible to simulate a delta function either computationally or experimentally. If the system is nonlinear and only the variation Av(t) is linear, it is futile even to attempt to simulate a delta function because the large amplitude required will drive the nonlinearities. A reasonable substitute for a delta function is an input that has l(jto) large for those to where the output, V(jto), is substantial. Then l(jto) can be calculated, v(t) measured, from which V(jto) can be calculated finally giving H(jto) = V(jto)/l(jto). The difficulty then comes in designing l(jto). If V (jto) varies greatly, calculating V(jto) from v(t) in the regions where V(jto) is small is difficult because much of the frequency domain will be dominated by the spillover from inaccuracies in the regions where V(jto) is large. KGS attempted to stimulate a step function at the towpoint in an experiment. A step function is the integral of a delta function and has a transform of 1/s. This requires a rapid change in configuration and, as they noted, the system probably went nonlinear in the experiment--a phenomenon that we have observed if the amplitude of the transverse oscillation is too large. The results of that particular experiment could not be reconciled with the results from other types of maneuvers. 2.2.4. Single-frequency response analysis. The alternate method is to drive the system for a sufficient length of time at a single-frequency and constant amplitude for an estimate of the amplitude and phase response to be made. The driving amplitude and phase are known and the calculation of the transfer function, H(jto), is simple, does not suffer from spillover from some other frequency, and does not require division by some calculated quantity, l(jto), which may be small or contaminaied with noise. The method is pedestrian, but it works. 2.2.5. Amplitude~phase estimation. The x,y,z time histories (h la Fig. 4) are analysed to determine the amplitude and phase of the response at the driving frequency and its second harmonic. The x and z data require some pre-processing: • A reference of bias is determined, as in Fig. 4, by averaging the last two cycles of the data and subtracting this from the original data record. • To minimize the influence of the transient interval on the estimates, we discard an initial portion of the data record, defined as that occurring before the first zero crossing after the initial cycle of data. In Fig. 4, this amounts to discarding about the first 500 sec of the z data. The amplitude and phase are then calculated as follows: An estimator (Helstrom, 1975) is defined as

E,(t) = t , S(T)ej''°~ d'r where S is the signal (the x,y or z data), to = 2-rrf is the angular driving frequency, and n = 1 for the linear (y) data; n = 2 for the nonlinear (x,z) data. E,,(t) tends to oscillate while converging towards the representative estimate that we desire. Its convergence

Simulation of an oscillatingtowed weight

155

rate is governed by the length and the quality of the signal, S. After verifying that it has more or less converged, we calculate an average over the last N cycles

E,.=

--~ d'r -At

where T is the end time, and At = 2~rN/nto. E, is a complex number embodying the amplitude and phase that we seek: the amplitude estimate is its magnitude, A, = [ E,[; and the phase estimate is its phase angle, 6 , = arctan(Im(-E,)/Re(E,)). Since this yields a phase of 90 degrees for the sine wave at the towpoint, 90 degrees is subtracted from the phase estimate to get the phase relative to the towpoint for the fundamental frequency, and 180 degrees for the harmonic. The choice of the positive exponent in the integral for E. has the consequence that d~ is not the phase of the transfer function, but its negative. 3. RESULTS AND ANALYSIS In this section, we begin by motivating the form of the graphs used to illustrate the linear and nonlinear system properties. We summarize the linear system behavior, discuss the source of certain nonlinear effects, summarize the nonlinear system behavior (including systematic distortion), and compare our results with Chapman's. We finally make some comments regarding TAAS's somewhat surprising ability to calculate smallscale effects. Figure 6 superimposes the y time-series of four points (nodes 4, 8, 12, 16) evenly spaced along the cable for the case with period 133 sec. The curves are labeled by approximate fractional cable position. The largest of the four curves is for node 4, about 1/4 of the distance towards the weight. The other three curves occur near the mid-point of the cable, 3/4 of the way down the cable, and the smallest curve is at the weight. The smallest curve (labeled 1) is the same y curve that was off the scale in Fig. 5. Figure 6 illustrates not only the y amplitude variation on different points on a cable, but also clearly shows the phase delay as the disturbance moves from the towpoint towards the weight. t/4

Period= 133 see

0

i

0

200

400

Time (sec) FIG. 6. Y time-series at 1/4 cable intervals.

600

156

T.N.

DELMER and T. C. STEPHENS

It was mentioned previously that our series of computer experiments produced 75 curves similar to the four shown here. The problem becomes one of distilling trends out of all this information. The data on Fig. 6 suggest that a useful summary might be a plot of the y response maximum amplitudes at each point on the cable and for all the cases considered. A similar plot summarizing the phase relationships would also be useful. The remainder of this section presents and discusses such plots for both linear and nonlinear data.

3.1. Summary of system linear properties Figure 7 summarizes the y linear amplitude response of the towed system. The points on each curve indicate the y amplitudes at individual cable nodes, calculated as described in Section 2.2. Each curve represents a separate case, labelled by driver period. Rather than plot actual amplitudes, Ay, we plot the ratio Ay/mowhere Ao is the amplitude of the driver oscillation, 10 m. The trends show that the amplitudes attenuate as one moves from towpoint towards the weight, and that the attenuation is greater for higher driver frequencies. One interesting feature is that the two cases with periods 133 sec, one with twice the driver amplitude of the other, produced curves that fall directly on top of one another. This is evidence of a linear response: Doubling the driver amplitude produced precisely double the response amplitude and the ratio remained constant. 1 l~&

Period (aec)

r \'<4~ 'J.o [- \ \ .< o

2 E .<

P

06~i

\

o - o ss.s 4

&~

m_

4

~A

m

~

\\\4 o

\o

•

I

266

o--n

z3a. Ao=20m

A....

4

~,

"4~

~o ~ ~" ~ m - - ~

O.2 O.0

~"~.

4--4

. . . .

4"4

A

'

° - - . o _ _ _ o~"0I - . . I _ _ m - - I - - I - - I . . . . '~ O ~ ° ~ O - - - ' 9 - - '

6

11

16

Cable Node

Fro. 7. Linear amplitude along cable.

Figure 8 summarizes the phase relationships for the linear oscillation component. Again, the points represent individual cable nodes and each curve is labelled by case. All phases are measured relative to the towpoint. During steady oscillation, each node's instantaneous phase of course changes, but all nodes maintain a constant relative phase. Again, the two cases with periods 133 sec fall directly on top of one another. As would be expected, the cases with higher frequencies produce larger phase differences between points, and the largest phase difference of all is that between the weight and towpoint for the highest frequency case. For that case the phase difference approaches 360 degrees, indicating that the weight is approximately in-phase with the towpoint, and the points in between should more or less describe a single oscillatory cycle. If that is

Simulation of an oscillating towed weight

157

360 --

Period (sec) 0 --

270

• --

/

•

/0

865

o

I-]--O 133 A-=20m '

,~

O/

u

/o /

/m /

/Jr

a~

...li/I ~ / A /

°/o" x - - = ~ j ~ j

-

-

T

.

,

i

.

.

.

.

,i.1.

i

.

.

.

-•-,

.

11

6 Cable

FXG. g.

S ~-

/.o/° o//

90

~11 ~1

/

180

~

O/I

/

0 66,5

@ - - 0 133 - - ~ 266

16

Node

Linear phase along cable.

true, then a snapshot of the system from above featuring the y displacements should reveal that the cable forms a single y cycle. Figure 9 is a snapshot of the system from above at time 800 sec for the high frequency case (period 66.5 see). Actual x coordinates of nodes measure hundreds of meters and have been replaced by the cable node index in order to make the y displacements visible. This snapshot occurred at a time when the towpoint (which is node 1 and currently at y displacement zero) was heading in an upward stroke, i.e. outward along the y axis toward the maximum amplitude of 10 m. Figure 8 suggested that the towpoint and the weight would be approximately in-phase, with the intermediate points registering a single cycle. That appears to be the case here. A motion picture of this plot would be interesting. It would show the towpoint sliding back and forth along the y axis between its limits ---10 m in cyclic fashion. The stroke and return in either direction would generate a single large bulge such as the one seen here which would attenuate and propagate along the cable towards the weight.

01

/

>.

/

• -2

At time=600 sec (Commencing upstroke)

./.

-3

/

~o / -4

.

.

.

.

Period 66.5sec a

.

.

.

.

6

u

11 Cable Node

Fno. 9.

Transversc cablc profile.

....

16

158

T.N. DELMERand T. C. STEPHENS

3.2. Source of nonlinear behavior Figure 10 shows a series of snapshots like the one in Fig. 9. Here are cable shapes shown at five distinct times during one half of a complete y period. These shapes show the positive half-cycle: the towpoint is moving away from the origin towards its positive maximum of 10 m and then back to the origin. The curve labelled 1 is the same curve shown in Fig. 9. The curves are labeled in time sequence. Especially evident here is how the displacement magnitudes vary from one end of the cable to the other. In fact for the second half of the cable the oscillation amplitudes are barely visible on this scale. 10

L 2

E

I/2 Period Series

2 0 -2 -4

~

P e . . . .

r

i 6

i

o d . . . . Cable

66.5 sec i . . . . 11 16

Node

FIG. I0. Transverse cable profile sequence. Now is the best opportunity to explain the x-z frequency doubling phenomenon. One can see in Fig. 10 that for each outward towpoint stroke (away from the origin) the entire system must be pulled slightly forward. The familiar "water-pulley" effect guarantees that the system will feel a net tug in the forward ( - x ) direction, the x coordinates will decrease, and the system will shorten to a small degree, relative to its steady state. Conversely, for each inward stroke (toward the origin), the system will tend to relax in the current, the x coordinates will increase, and the system will lengthen. The towpoint moves outward from the origin twice during each period of its oscillation; the tugs forward in x are thus bound to occur at twice the frequency of the y towpoint oscillation. The key point is that it matters not whether the towpoint y velocity is positive or negative so long as it is moving away from the origin for the rest of the system to feel a forward tug. The reason that the z oscillation also experiences a doubled frequency is due to the shape of the cable itself hanging in the current. The cable is bent through nearly a 90 degree curve throughout steady state oscillation. Therefore, any displacements aligned along x near the top of the cable are transformed smoothly into displacements aligned along z near the bottom of the cable. For the towed weight case the x and z reactions must be viewed as being coupled. A simpler way to grasp this nonlinear effect is this: immersed cables vigorously resist perpendicular movements but slide easily along their axis, encased in a "water sheath"

Simulation of an oscillating towed weight

159

(the source of the "water-pulley" effect). This ensures that transvers.e displacements of any point on a cable will induce axial movements of the rest of the cable; the degree of axial motion depends on the cable tension, bending rigidity, distance from the point and of course freedom from outside constraints. The towed cable's nonlinear reaction to the crosstrack towpoint motion manifests mainly as sliding movements along its length. The separate x - z oscillation records are artifacts of the cable shape. 3.3. Summary o f system nonlinear properties

Figure 11 is a summary of the nonlinear amplitude response along the cable. The format is similar to Fig. 6. Since the x - z oscillations are coupled, ,;ve plot the amplitudo ratio of the vector magnitude of the oscillation in the x - z plane to the amplitude of the towpoint. There is no x - z reaction at the towpoint (node number 1) since constraints restrict movement to the y direction. The most notable feature of this plot is that the nonlinear oscillations do not attenuate as one moves towards the weight. Instead they appear to be generated in a region that covers the first half of the cable and then remain more or less constant over the second half of the cable. The magnitude of the oscillations appears small though. The maximum indicated here represents about 1.5% of the driver magnitude. Finally, note that of the two curves labeled 133 sec, the one with driver amplitude 20 m (twice the nominal amplitude, indicated by open squares) lies on a curve that is consistently twice the magnitude of the curve representing the nominal driver amplitude (with filled circles). On the linear amplitude plot, Fig. 6, those two curves overlaid one another indicating a linear response. Here they are clearly distinct. Doubling the driver amplitude produced four times the magnitude of the response. The amplitude ratio is thus 2 to 1 and the response is nonlinear. The nonlinear cable mode has been described as primarily longitudinal, i.e. oriented along the cable. Longitudinal modes propagate much faster than transverse modes: pulling on a cable gets more action to the far end faster than shaking it up and down. Though the amplitudes appear relatively small, the disturbances propagate rapidly and without much attenuation throughout the system. This mode does not show up in the linear analysis.

0.020 Period (~ec)

,<0

0.015

"~"0-...~

0/ /

0.010

O~ °'~°~o--o~o--o--°~°--°-o / e ~ e ~ e ~ e ~ e ~ e

, O/

0.005

133

EI~+~D~D~D__D~O__ D_

/ / ~ 0 ~ 0 ~ o ¢:

86.5

O--e

z ~ - - A 266 • - - • 685 n--o 133. Ao-20m

O--"~O~O--...O /

0--0

.......A ~

~O~O~o--o--o--o--q t, ~ t, ~ t , ~ ,~~ , ~ _ _

~__~__

t, _ _ ~

~

8 <~ 0.0o0

_

"A"'-, A

,

,

L

.

.

.

.

6

I

.

.

.

II

Cable Node FIG. II. Nonlinear amplitude along cable.

.

16

160

T . N . DELMER and T. C. STEPHENS

We calculated phases for the nonlinear x and z motions similar to that for the linear case shown in Fig. 8. But the two phases did not match, indicating that the points on the cable undergo elliptic motion, and making it impossible to define a single phase precisely. Accordingly, we present no plot. A systematic distortion, or shift of the average x-z planar positions of the nodes away from their initial steady-state, was noted for Figs 4 and 5. Figure 12 plots the magnitude of the shift as a function of cable position, for all cases. This is another nonlinear effect, as can be seen by the positions of the two cases labelled 133 sec. Doubling the oscillation amplitude of the nominal case produced a net shift four times the original amount. Note that though the magnitudes of the shifts are small relative to the whole system, they are significantly larger than the oscillation amplitudes in the x-z plane. (In the nominal 133 sec case, the weight shifts 30 cm, but oscillates only about 6 cm, c.f., Fig. 11.)

i .2

0/0 /

Period (see)

El/

O-- 0 86.5 ,--o

~--

0.9

/~

D--Q 133,Ao=20m

/ ~

0.6

//0

133

tx 266

A--. s85

~0~°~°~0-0-0-°-0-°-0-0-~

/o

6

11

16

Cable Node Fro. 12. Systematic distortion magnitudes.

The directions and relative sizes of the shifts are illustrated in Fig. 13. Here, for each case, the distorted cable configuration is plotted as a dashed curve, the original steady state is the solid curve (the same for all cases), and vectors connect the corresponding nodes on both curves. The shift for each case has been magnified by a constant factor for visibility and to maintain about the same size vectors across cases. The cases are ordered from top to bottom according to decreasing magnitude of the shift. All the nodes move upward and downstream except those in the lowest frequency case, which actually shift slightly upstream towards the towpoint. (The slope of the vectors seems to vary smoothly, though, from high frequency to low; between cases 266 and 665, there probably is a case that shifts straight upward. The slopes also appear unaffected by oscillation amplitude.) 3.4. Comparison with Chapman Table 1 is a comparison of the TAAS results with Chapman's published results. This table summarizes the linear or y reactions, only. Shown here for the four primary cases

Simulation of an oscillating towed weight

•

"

"

",.

Period

,-.

[magnificai°n]

161

(aee)

133, Ao 2 0 m

-,,

[lOOx]

66.5

[lOOx]

133

[~00~]

266

[~OOx]

665

[2ooox] FI~. 13. Systematic vector displacement (magnified).

TABLE1. LINEARRESPONSECOMPARISONWITHCHAPMAN Transfer function---Linear case Period (sec) Frequency (msec) Amplitude ratio Phase (deg.)

66.5 15.0 8.32E°3 (8.78E-3) -343 (-329)

(Chapman's results in parentheses).

133 7.52 4.53E-2 (4.71E-2) -236 (-231)

266 3.76 1.45E-1 (1.52E-1) -162 (- 162)

665 1.50 4.80E-1 (4.30E-1) -96 (-98)

162

T . N . DELMER and T. C. STEPHENS

(labeled by period in seconds) are the amplitude ratio and relative phase figures that constitute the transfer function for the weight (the transfer function phases are the negative of those shown in Fig. 8). Chapman's results are in parentheses. The figures are within 5-10% of one another. It is not too surprising that the two approaches agree since both have been validated for various dynamic cases (DSC; DST; Bettles and Chapman, 1985). Table 2 is a similar comparison for the nonlinear case. (The concept of "transfer function", strictly speaking, only applies to the linear case.) Chapman did not calculate nonlinear reactions in his paper. However, we modified his published computer program (described in our Appendix C) to produce estimations of the nonlinear amplitude ratios. These ratios are for 10 m amplitude, and vary linearly with amplitude. Although the phase is ambiguous at other points, at the weight the motion is primarily in the z direction, so the phase is well defined there. Phases for z are presented in the table. Notice that the values are small, very nearly in-phase, which is due to the high speed of the longitudinal waves. Once again the amplitude figures agree to within 5-10%; the phases to within one tenth of a cycle. TABLE 2. NONLINEAR RESPONSE COMPARISONWITH CHAPMAN Transfer function--Nonlinear case Period (sec) Frequency (msec) Amplitude ratio z Phase (deg.)

66.5 15.0 8.5E-3 (7.8E-3) 21 (46)

133 7.52 5.7E-3 (5.4E-3) 34 (47)

266 3.76 4.1E-3 (3.7E-3) 41 (40)

665 1.50 3.7E-3 (3.5E-3) 62 (52)

(Chapman's results in parcntheses).

The systematic distortion is not predicted by the purely geometric considerations described in Appendix C, so no comparison can be made. But some simple physical arguments can be used to shed light on the origin and nature of the shifts: • A curved cable necessarily projects a shorter profile than the same cable straight. Thus the y oscillations cause a shortening of the projected length of the cable in the x-z plane--a pure geometric shift. • The broadside aspect of the cable elements is increased relative to the current by virtue of rotation out of the x - z plane. Since the elements already possess varying degrees of angle-of-attack (and lift), the result is increased lift. • The average increase in the speed of the cable element due to oscillation also contributes to increased drag (and more lift due to existing angles-of-attack). We have not studied the relative effectiveness of the above three processes. All are nonlinear since the direction of the oscillation can have no effect..All would seem to induce upward displacements, consistent with those observed in Fig. 13. They may also be sufficient to explain the upstream/downstream variation of the shift vectors: the latter two processes in addition to increasing lift would tend to increase drag downstream, opposing the upstream tendency of the geometric shift. Higher frequency cases would naturally produce enhanced drag which would tend to obscure a simple

Simulation of an oscillating towed weight

163

geometric shortening. In sum, the systematic distortion is caused by various terms in the equations of motion of the cable representing the interplay of complex, three dimensional effects. 3.5. Comments on response size and error control The response of the weight---especially the nonlinear part--is small, being on the order of a few centimeters. It is surprising that the simulation of a cable one kilometer long can produce reliable results for such small displacements, particularly when the root-mean-square (RMS) error criteria for a single time step is 5 cm. Furthermore, the shape of the cable is a curve and is represented in the simulation as straight lines between nodes. This introduces an intrinsic, i.e. unavoidable, error which is on the order of a meter for the parameters in the simulation. The fact that the simulation produces reasonable results in the face of the large error allowed by the integrator and the intrinsic modeling deserves a brief discussion. The integrator calculates an error for each of the three components of each node position and velocity. It estimates this error by comparing an extrapolation from past data to the update time with the integrated value. The total error estimate is the r.m.s. of the set of individual errors. If the r.m.s, value is too large, the time step is decreased and values at the new time are estimated and calculated. The procedure is iterated. If the error is small, the fact is recorded as a positive factor for increasing the time step (one of several factors influencing the modification of the time step). The real error in a time step, the difference between the calculated and the true solution, may be larger or smaller than estimated because the estimator, as its name implies, is imperfect. It is possible that an excessively small time step will be produced because of the inaccuracy and that the actual results will be much more accurate than expected. The opposite can also occur: the estimated and calculated values can be in error in the same direction compared to the true solution and the integration will diverge from the true solution. Since the error is r.m.s., some variables can contribute more than others to the sum. It can be expected that those variables having the largest change in a time step would contribute most to the error and, in turn, maintain a small time step. The time step would then be such that the error in the slowly changing variables would be small. In addition, the error criteria is local, i.e. reflects the current step only, so the r.m.s, value can recur at each step. The total error in the integration--the global error----can be the sum of the local errors made at each time step. There is no known method for estimating the global error accurately. In reality, the error criteria used in the integration does not have much to do with the error in the final results, at least not directly. Rather, the criteria is a control on the integration process. The control is used to force the solution to behave reasonably. For too large a value of the error control, both the time step determined by the code and the solution vary radically, and the integration generates rubbish. As the error control is tightened, the time step and solution jump now and again, but a solution is obtained. Finally, when the error control is small enough, the time step and solution are smooth. In addition, the variation of the solution and the computer time required to produce the solution vary smoothly and monotonically with the error control. The solution is then approaching the true solution. For values of the error control where such behavior takes place, the solution converges as a function of the decreasing error

164

T.N. DELMERand T. C. STEPHENS

control. If the difference between the calculated and true solutions is called noise, it is clear that this noise is not similar to that typical in experimental measurements. The latter can usually be approximated by filtered white noise. The noise in the integrator would be more similar to experimental noise if the integration were driven by a stochastic process. In the latter case, the solution would not converge as a function of error control but merely converge in some average sense, i.e. converge in the mean, under the optimistic assumption that the integrator would produce any solution at all. The smooth relation of the convergence of the numerical integration with error control is related to the convergence with time step since reducing the error control reduces the time step. Evidence of smooth convergence with time step, which can be taken as evidence that the noise in the integration is not stochastic, is given by the success of an integration method based explicitly on this convergence, the Bulirsh-Stoer method. In contrast to the Gear backward difference method used in TAAS, the Bulirsh-Stoer method (Press et al., 1986) does several integrations between two points in time using different time steps and then uses an approximation to the solution in the form of a ratio of two polynomials in the time step to extrapolate the solution to zero time step. If the solution had a stochastic component, such an extrapolation would not work, and some sort of estimation taking account of the stochastic property of the noise would be necessary. One such method is Kalman filtering, commonly used to estimate signals in noise when the correlation function of the noise is known approximately. The intrinsic error made in approximating a curve by straight lines is a bias in the simulation. The error varies inversely as the product of the square of the radius of curvature and the square of the number of links between the end nodes. Since neither of these parameters changes significantly during the simulation, the bias does not change. By looking only at the changes in the solution between the oscillating and nonoscillating case, the intrinsic error is subtracted out. It is probable that in certain scenarios the small changes in solution would not be observable without reducing the error criteria. This could make the required computation time prohibitive. A more likely cause of failure is to attempt to simulate high frequency response. High frequencies require more nodes and smaller time steps simultaneously, which increases the required computation time radically. Eventually, for high enough frequencies, one must treat the frequency as a parameter as Chapman did, rather than try to integrate the amplitude variations caused by the oscillation.

4. CONCLUSIONS Three major conclusions can be drawn from this work: • Full simulations like TAAS can be used successfully to study linear, first-order system responses. In particular, the simulation of a towed weight oscillating about a steadystate configuration compared favorably with an independent perturbation calculation (Chapman) for a large range of parameters. TAAS can therefore model small changes in a large system as well as the transmission of oscillations throughout a system. • Full simulations are also effective at calculating nonlinear, second-order system responses. Comparison of the TAAS simulated nonlinear reactions of the oscillating

Simulation of an oscillating towed weight

165

towed weight with a slight generalization of Chapman's perturbation calculation showed that second-order effects are modelled correctly. T A A S can therefore illuminate nonlinear b e h a v i o r whereas linear analysis cannot. Nor does a linear analysis reveal when such effects are important. This last point is justification for the use of full simulations even when it appears that a linear analysis would suffice: linear analysis can fail without giving a hint that it has failed. • The primary influence of crosstrack towpoint motion on a towed object may be nonlinear. Our simulation has affirmed that it is the nature of towed cable systems to convert a portion of any transverse cable movement into tangential motion, i.e. motion along the cable axis. This effect is non-linear: it is practically independent of the radial direction of the perturbation, and the spectrum of the response differs radically from the driver spectrum. Since axial cable motion is only lightly damped, the cable serves as an efficient conduit of nonlinear disturbances throughout the system. While the results presented here are for small scale motion, it is important to bear in mind that full simulations do not require that the motion be small or that the system be simple, e.g. have uniform properties, both of which were essential to the perturbation analysis. TAAS has now been validated for large-scale nonlinear motion (DSC, DST), and for small-scale linear and nonlinear motion. This indicates, although it does not prove, that of the complex system. The most commonly occurring sources of complexity in the effects of small oscillations about such motion. The main source of doubt is in the possibility of some qualitative difference between the response of the simple system and that of the complex systgem. The most commonly occurring sources of complexity in the systems that we study are material and/or geometric nonuniformities, including discontinuities in the cable properties. This is already present in the towed weight problem at the weight itself, and the results presented here show that it is treated correctly. The conjecture then is that T A A S could also simulate the effects of small towpoint perturbations on the motion of a towed array. The discontinuity of angle that occurs at the interface between the cable and array is not contained in the towed weight problem, but the inductive leap from towed weight to towed array is credible because there is no fundamental difference in the equations used to simulate the two problems. The fact that T A A S has been shown to handle one type of discontinuity (the weight) indicates that it should handle other discontinuities. We hope to further substantiate this generalization. Acknowledgements--This work was supported by the Defense Research Establishment, Atlantic under contract W7707-6--7832/01-OSC. We wish to thank the contract technical monitor, Jim Tremills, for helpful comments and. for presenting this material at several meetings.

REFERENCES BET'FEES, R.W. and CHAPMAN,D.A. 1985. Experimental verification of a towed body and cable dynamic

response theory. Ocean Engng 12(5), 453-469. CHAPMAN,D.A. 1984a. Towed cable behaviour during ship turning manoeuvres. Ocean Engng 11(4), 327-361.

166

T . N . DELMERand T. C. STEPHENS

CrtAVMAN,D.A. 1984b. A study of the ship induced roll-yaw motion of a heavy towed fish. Ocean Engng 11(6), 627-654. DELMER,T.N., STEPHENS,T.C. and Cot, J. 1983. Numerical simulation of towed cables. Ocean Engng 10(2), 119-132. DELMER, T.N., STEI'HIZNS,T.C. and TREMZLLS,J.A. 1988. Numerical simulation of cable-towed acoustic arrays. Ocean Engng 15(6), 511-548. D i s c , P.A.M. 1947. Quantum Mechanics. Oxford University Press. London. DOWL]SG, A.P. 1988a. The dynamics of towed flexible cylinders, Part 1. Neutrally buoyant elements. J. Fluid Mech. 187, 507-532. DOWLINO, A.P. 1988b. The dynamics of towed flexible cylinders, Part 2. Negatively buoyant elements. J. Fluid Mech. 187, 533-571. HELSTROr,I, C.W. 1975. Statistical Theory of Signal Detection, pp. 262-274. Pergamon Press, Oxford. KEr~NEDV, R,M. and GALL~-~A, F. 1981. Cable dynamics experiment. Technical Report 6467. Naval Underwater Systems Center, New London, CT. KENNEDY,R.M. and S~AHLEN, E.S. 1981. A linear theory of transverse cable dynamics at low frequencies. Technical Report 6463. Naval Underwater Systems Center, New London, CT. PRESS, W.H., FLANr~ERY,B.P., TEUKOLSKY,S.A. and VETrERLIr~6,W.T. 1986. NumericalRecipes. Cambridge University Press, Cambridge. APPENDIX

A

Sample simulation inputs for stationary steady state phase . .

................... ...................

To~ed Weight problem . . . . . . . . . . . . . . . . . . . . . L a t e r a l O s c i L L a t i o n Test . . . . . . . . . . . . . . . . . . . . .

N

" " " .

This t e s t a p p l i e s a sinusoidaL Y o s c i L L a t i o n t o the system t o t l p o i n t , Mhtch i s o t h e r w i s e f i x e d s t t h e o r i g i n . The system i s f i r s t immersed i n s u n i f o r m +X c u r r e n t zmd reaches s t e a d y s t a t e . Then t h e o s c i | t s t i o n begins, and the system i s ettowed t o s e t t l e i n t o " o s c i L L a t o r y " s t e a d y - s t a t e . At t h a t stage, t h e weight i s sampled f o r o s c i L L a t i o n s .

N

STRUCTURAL DEFINITION LOADING FUNCTIONS CLASSaBARE ORAG"1.07 NORNAL=(SIN(T)**2] TANGENTIAL=.O01 TENSION FUNCTIONS CLASS=STEEL STRAIN= [E+EDOT*. 003] STRESS: [3.4E7"S] TOM CABLE LENGTH=IO00 DIAMETER=.03 INERTIAL NASS:I.02 GRAVITATIONAL NASS:1.02 BUOYANCY=I DRAG:BARE . ELAST| C! TY=STEEL ROOES=14" CALL NODE 1 (ToMpoint) CALL LINK 1 (CabLe) CALL NODE 16 (Weight) COUPLl NG:O TORSION=O MOMENT PER LENGTH:I ST I FFRESS=O ANGULAR DAHPING=O

15-(~S. 7m sections

Simulation of an oscillating towed weight

167

SPACING RATE=I IN-LINE PACKAGE EXTRA HASS ON:iRIGHT =3400 SCENARIO N

" The defautt | n i t | a t i z s t i o n has the system aligned, unstretched, stong the +X a x i s w i t h " the towpo|nt s t the o r i g i n . The we|ght w i l t s i n k the cabte dynamicatty t o the proper " depth. n

CONSTRAIN TOM~OIIlT VX=O" CONSTRAIN TO~OINT VY=O CONSTRAIN TOJPOINT VZ=O X CURRENT SPEED --2.5 PROGRAM CONTROLS STOP T IH~=3000 ERROR CRITERION: (0, .05) END

Fix the t o w p o i n t

Sample simulation inputs for oscillatory phase .

...................

taterat Oseitlntion

Test

.....................

N

" " " "

This t e s t a p p l i e s a s i n u s o i d a t Y o s e i t t a t | o n t o the system tOUl:oint, which i s otherwise fixed st the o r i g i n . The system is assumed to be at steady state i n a uniform +X current. The o s c i t t a t i e n begins, and the system i s attowed to s e t t t e i n t o n s c l t t a t o r y steady s t a t e .

N

SCENARIO N

BEGIN: WTSTEAD.FRZ CONSTRAIN TOM~)INT Y=[lO*SIN(2*PI*T/133)] CONSTRAIN TOWPOINT VX=O CONSTRAIN TOUPOINT VZ=O X CURRENT SPEED = 2.5 PROGRAM CONTROLS STOP TINE=1200 ERROR CRITERION: (0, .05) ENO

APPENDIX B This appendix discusses differences in the Chapman and T A A S approaches to the treatment of the boundary condition at the weight. Chapman's work is particularly useful in the present context because his met'hod is completely different from ours. He considers the steady state towing and solves this by integrating an ordinary differential equation with scope as the independent variable. Our solution (DSC) uses the method of lines to convert the partial differential equation for the cable into ordinary differential equations for the motion of points on the cable with time as the independent variable. In Chapman's work the magnitude of the variation in time is considered small and at a single frequency. He considers the amplitude of this motion as a new variable which is a function of scope; time is replaced by frequency, which is a constant. Chapman's method has advantages over other approaches (Kennedy and Galletta, Kennedy and Strahlen, Dowling) in that it integrates the equations in scope numerically rather than assuming an approximate solution, e.g. a straight line. To make the comparison definitive, with no excuses to fall back on as explanations for discrepancies, we chose to compare our results with Chapman's. However, late in the project, we noticed that there is one implicit assumption in Chapman's work: the inertial mass of the weight is negligible. For the assumption to be true, the physical properties of the cable and weight must satisfy a certain relation which is given here, and is satisfied by the parameters used in our simulations. This relation is not necessarily satisfied

168

T.N. DELMERand T. C. STEPHENS

by heavy weights, but is satisfied in many cases. In addition, a new boundary condition at the weight is suggested that accounts for the inertial mass. To have included this condition would have required an additional parameter, or dimension, in Chapman's tables. Had we realized this implicit assumption was present early in the project, we would have used a downward force as the boundary condition rather than simulating a weight. In our work and Chapman's the weight is assumed to have no drag so the steady state solution would have the end of the cable at the weight vertical. This is evident in the boundary condition in Appendix A of Chapman (1984a) where the angle of the cable with respect to the horizontal is ~r/2 or 90°. At the weight the amplitude of the oscillation is set to unity and the change of the amplitude with respect to the scope (measured from the weight) is zero. However, to get the weight to oscillate at frequency to as cos tot requires a force equal to the mass times the acceleration. For a weight, W, having gravitational force on it, using g as the gravitational acceleration, the magnitude of the required force is Wto2/g. This force must be supplied by the horizontal component of the tension. The magnitude of the tension is W, the gravitational force on the weight. If the cable makes a small angle with the vertical, 0 = d13/ds, where s is the scope, the horizontal force is approximately W0. This gives the boundary condition WdfS/ds = Wto213/g. This equation contains a parameter, g, that did not appear previously and adds a new dimensionless parameter to the analysis. One possible choice of this parameter is Vto/g where V is the towing speed; this parameter does not go to zero for the heavy weights that Chapman considers so the parameter can be important. To estimate the importance of this boundary condition, consider the contributions made in the equation for the second derivative, d213/ds2, from the drag and the first derivative at a position near the weight. The drag is in terms of the drag coefficient, Ca, the diameter of the cable, d, and the density of the water, p. The ratio of the terms, which must be small if the boundary condition of zero first derivative is to hold, is w d13/ds 0.59V-C~d '~ 1, where w is the weight in water of the cable per unit length. The denominator here is Chapman's roto/V, so substituting the boundary condition for the first derivative gives the requirement

The left-hand side of this equation is 0.002 for the parameters used here so the inequality is satisfied. The boundary condition can be translated into_Chapman's terminology using dimensionless scope, ~ = sw/W and his complex displacement, 13, to give d~ dtr

Wto2 gw

at the weight. APPENDIX C This appendix describes the extension of the code Chapman (1984b) listed in his appendix B and lists the extended code. His work is generalized to give the first nonlinear, i.e. quadratic, result. Two caveats about the extended code should be made immediately. It is clumsy and slow and it should not be extended further. The code published here is not the code used to calculate the results, but is a compaction of the original using the feature of BASIC that allows several instructions to be placed on a single line; this reduces the number of pages and the empty space on a page. This, along with the idea of publishing the code, is in the spirit of Chapman; his listing has several instructions per line and has one error which spoils the results (the last letter on line 5040 should be R, not S). Due to various exigencies, published codes are almost never the codes used to obtain the published results. We rely on the differential geometry and use the notation of Chapman (1984a,b). The nonlinear response is due to the transverse displacement of the cable as discussed in conjunction with Fig. 10 of the main text. Let s be the scope as a coordinate on the steady state cable and I3 be the

Simulation of an oscillating towed weight

169

transverse displacement due to the oscillation. The change in actual distance from the towpoint will be

where the integral is from beginning to the end of the cable. This expression is good for small 13. Expanding this gives

In terms of Chapman's scaled complex displacement, ~, and scaled scope, ~, the derivative is

d_~= ds

Re e i'°' d13 dtr

where to is the circular frequency of oscillation, t the time, and Re means take the real part of the following quantity. Using the same notation as in our Appendix B, the previous equation gives

+ o., co,

-

-

.,

sin 2tot [_dry dcrJ"

The doubling of the frequency and the presence of a zero frequency component is due to the squaring of the cosines and sines. In terms of the three integrals

(d~,y 1, = I i, = I

\d~r} d~r Ida,/2 \ d ~ ] der

d13,~-~ d13i d~ , l,-i = f ~-~ the amplitude and phase of A are given by

A2,~ = O.25 ( W ) ~/(I, - I,)Z + 41~ tan(~b2~,) = (-~Iri, (It - 1,)). The subscripts_denote the frequencies of the components. The phase of the fundamental is given by tan ~b~, = (13~,13r) evaluated at the towpoint. The notation for the tan is the same as used in the F O R T R A N ATAN2 function, two values being required to invert the tangent and put the angle in the correct quadrant. The three variables, 1,, I,, and Irg were added to Chapman's variables in his code. To compare the results from the integration to those from the simulation, it is necessary to take into account the fact that Chapman starts integrating at the weight. Our phase comparisons are in the opposite

170

T.N. DELMERand T. C. STEPHENS

direction; the phase of the response is calculated relative to the phase of the towpoint so a sign change is necessary and is reflected in the code. Chapman tabulates 10 loglo ~ The extended code calculates the dimensionless response, R, for a unit driving amplitude in the nondimensional variables and prints AMP2=10 logloR. The unit driver corresponds to an amplitude of W/w (m). To get the response of the weight to a driver of amplitude A, the results must be scaled to A by the factor (Aw/W) 2 to account for the difference in the driver's amplitude. To put the dimension of distance into the result, the factor of W/w must be included. The (dimensional) amplitude of the response is related to the code output by Amplitude = (A-~wW)10AMP2~I0 . The complication comes from using dimensionless quantities and having a nonlinear response. The result for the phase does not depend on the amplitude. Note that the usual concept of a decibel would replace the factor of 10 by 20 since the decibel refers to power which is related to the square of the amplitude. The following is the listing of the modification of Chapman's code; it was run using the IBM BASICA interpreter on an IBM PC. 10 20 30 50 60 70 80 90 100 110 120 130 140 150 160 170 180 199 200 210 220 230 240 250 260 280 290 310 320 330 380 410 420 &30 440 450 460 470 480

REM ***LATERAL ATTENUATION PROGRAM*** DIN Y ( 9 ) , Y S ( 9 ) , Y O ( 9 ) , F O ( 9 ) , F B ( 9 ) , F S ( 9 ) , T L ( Q ) , E ( 9 ) EQ=9 : PI-4*ATN(1) : LOGIO-LOG(IO) REM Y(O) IS THE ANGLE BEN Y(1) IS THE TENSION BEN Y(2) IS THE DEPTH REM Y(3) IS BETA (BAR) SUB R REM Y(4) IS BETA SUB I REM Y(5) IS D BETA SUB R D SIGMA REN Y(6) IS D BETA SUB I D SIGMA BEN Y(7) IS INTEGRAL OF (D BETA R D SIGMA)^2 REN Y(8) IS INTEGRAL OF (D BETA [,D SIGMA)^2 REM Y(9) IS INTEGRAL OF D BETA R D SIGMA TIMES REM D BETA I D SIGMA INPUT "NON-DIHENSIONAL PERIODM;P IF P<=O THEN PRINT "ERRORn:GOTO 160 INPUT "CABLE WEIGHT/DRAG RATIOm;GR IF Mt"F" AND TSo"B" THEN PRINT "ERROR":GOTO 200 INPUT "MAX NON-DIMENSIONAL SCOPE";SM IF SMSM THEN PRINT "ERROR":GOTO 240 MSuSQR(WR) : OM=2*Pl/P REM *PRINT HEADINGS* S$="BARE" : IF T$="F" THEN S$="FAIRED" CLS PRINT " PERIOD W/RO ";S$ S$-"" : W=7 : D=4 : VO=P : GOSUB 950 V-16 : VO,4,1t : GOSUB 950 PRINT S$;" CABLE" PRINT PRINT " SCOPE DEPTH A N G L E ATTN(DB) PHASE AHP2(2ND)ANGLE(2ND) AHPO(ATO)" REM *SOLUTION: INITIALIZE** FOR I : 2 TO EQ Y(1)=O NEXT Y(O)=PI/2 : Y(1)=1 : Y(3)=1 : X=O : HE=DX/4

Simulation of an oscillating towed weight 53O 540 550 56O 57O 58O 59O 6O0 650 68O 7O5 7'07 7'10 740 750 7'70 7'74 775 78O 79O 8O0 820 825 830 840 86O 88O 9OO 910 920 930 940 950 96O 1000 1020 1030 1040 1050 1060 1070 1080 1090 1100 1120 1130 1160 1170 1180 1190 1220 1230 1240 1250 1260 1270 1280 1320 1330 17,;0

FOR I:O TO EQ D8.0001 DO=D*ABS(Y(I ) ) IF DO>O THEN D=DO TL(I)=D NEXT REM *PRINT RESULTS* S$:"" : VO=X : W=7 : 0=4 : GOSLIB 950 VO=T(2) : V:lO : GOSUB 950 VOsT(O)*180/PI : D:3 : GOSUB 950 DRIVER=SWt(T(3)^2+Y(4)^2) DRIVERmlO*LOR(DRI VER)/LOGIO VO=ORIVER : D=4 : GOSUB 950 SI NUSsT(4) : COSlNUSsy(3) : GOSUB 1690 VO=ATN2*180/PI : GOSUB 950 SI NUS=+2*T(9)*(T(3)^2-Y(4)^2)-2*(Y(7)-Y(8))*Y(3)*Y(4) REN THIS SIGN REVERSES THE REFERENCE POINT SINUS=-SINUS COSINUSs+4*T(9)*Y(3)*Y(4)+(Y(7)-Y(8))*(Y(3)^2-Y(4)^2) GOS~ 1(#0 PHASE2"ATN2*180/PI : N4POC=Y(7)+Y(8) IF AHPDC > 0 THEN NSPDC=LOG(AHPOC/4)*IO/LOGIO-2*DRIVER N4P2=.O625*(T(7)-Y(8))^2+.25*Y(9)^2 IF'AMP2 • 0 THEN AHP2=LOG(AMP2)*5/LOGIO-2*gRIVER VOsNIP2 : __r~_JIS950 VO8pHASE2 : GOSUR950 VO-NIPOC : GOSUR 950 PRINT S$ IF X>SN-.O0001 THEN PRINT:END REH *MOVE TO NEXT SCOPE* GOSUR 1120 GOTO 530 REM **OUTPUT FORMATTER** VI,,INT(VO*IO^D+.5) : C80 : KS="" : V8ABS(V1) Q81NT(V/IO) : V:V-Q*IO gS'44105("0125450789", V+1,1 )+KS C8C+1 IF C<,4) OR Q<>O THEN V:Q:GOTO 1000 IF VlW THEN 1100 KS:LEFTS(" ",W'C)+K$ S~S$+K$ : RETURN REN **RUNGE'ICUTTA-MERSON METHOD FOR SOLVING ODEnS* * XOsX : KsO : HsHE FOR 180 TO EQ YO( I ) . y ( l ) NEXT 8180 : HOeH : DU=X4"DX-XO IF(DU-H)/DX<.O00001 THEN BI=-I:H=DU EXsXO FOR I=O TO EQ YS(l)-YO(1) NEXT GOSU8 1800 H3:H/3 : H6zH/6 : HS=H/8 : H:'sH/2 FOR I=O TO EQ FO( I )=FS( I ) YS( I )sYO( ! )+H3*FO( I )

171

172

1350 1360 1380 1390 1600 1610 1620 1630 1640 1650 1660 1680 1690 1500 1510 1520 1.%0 1550 1560 1570 1580 1590 1600 1610 1630 164,0 1650 1660 1670 1690 1700 17'10 1720 1730 1740 1750 1770 1780 1790 1800 1810 1840 1850 1880 1890 19"50 1940 1950 1960 1970 1990

T . N . DELMER and T. C. STEPHENS

NEXT EX=XO+H3 : GOSUg 1800 FOR I=O TO EQ YS(1):YO(I)+H6*(FO(I)+FS([)) NEXT _r~um 1800 FOR l=O TO EQ FB(1):FS(I) YS(I):YO(1)+H8*(FO(I)+3*FB(1)) NEXT EX=XO+H2 : GOSUR 1800 FOR I=0 TO EQ FO(1)=FO(I)÷6*FS(!) YS(I)=YO(I)+H2*(FO(i)-3*FO(1)) NEXT EX=XO+H : _~uJ6 1800 FOR [=0 TO EQ Y(I)=YO(I)+H6*(FO(I)+FS(1)) NEXT FOR I : 0 TO EQ E(I)=AOS(YS(I )-Y( [ ) ) / 5 IF E(1)>TL(i) THEN H=H/2:GOTO 1190 NEXT K=K+I : XO=XO÷H [F H1 THEN HE=HO:X=XO:RETURN FOR [=0 TO EO IF E(1)>TL(1)/64 THEN 1160 NEXT H=H*2 : GOTO 1160 REN *ARCTANGENT (ATN2) OF COS]NUS ANO SINUS* REM NEEOS A VALUE OF PI IF COSINUS <> 0 THEN 1750 IF SINUS • 0 THEN ATN2s.5*PI:RETURN IF SINUS < 0 THEN ATN2=I.5*PI:RETURH ATN2sO:RETURN TANG=SlNUS/COSINUS : ATN2=ATN(TANG) IF COSlNUS < 0 THEN ATN2=ATN2÷PI:RETURN |F ATN2 < 0 THEN ATN2=ATN2÷2*PI RETURN REN **CALCULAT|ON OF OERIVATIVES** YSsYS(O) : SAsSlN(YS) : SO=SA IF S8<.5 AND TSsMF" THEN SBz.273+.82?~YS*YS CA=COS(YS) : FS(O)=(CA'SB*SA/IJR)/YS(1) : 01=SA IF TS=mFm THEN 01=SA÷SO*CA/WR FS(1)=D1 : FS(2)sSA : FS(3)sYS(5) : FS(4)=YS(6) 02"Ol*Sa/Mt IF TS='"O" THEN DI"OI+SA *C,A/~ FS(5)8- (01*YS(5)+D2*YS(6))/YS(1) FS(6)s(O2*YS(3)-OI*YS(6))/YS(1 ) FS(7)=YS(5)^2 : FS(8)=YS(6)^2 FS(9)=YS(5)*YS(6) RETURN