THE UNIVERSIT Y OF MICHIGAN COLLEGE OF ENGINEERING Department of Civil Engineering Final Report PROPELLANT LINE DYNAMICS E. B. Wylie W. Zielke Project Director: V. L. Streeter ORA Project 08040 under contract with: NATIONAL AERONAUTICS AND SPACE ADMINISTRATION GEORGE C. MARSHALL SPACE FLIGHT CENTER CONTRACT NO. NAS 8-20312 HUNTSVILLE, ALABAMA administered through: OFFICE OF RESEARCH ADMINISTRATION ANN ARBOR September 1967

This report was prepared by The University of Michigan College of Engineering Department of Civil Engineering under Contract No. NAS 8-20312 PROPELLANT LINE DYNAMICS for the George C. Marshall Space Flight Center of the National Aeronautics and Space Administration. The work was administered under the technical direction of the Propulsion and Vehicle Engineering Laboratory, George Co Marshall Space Flight Center with John R, Admire acting as project manager.

TABLE OF CONTENTS Page LIST OF ILLUSTRATIONS......................................... v LIST OF TABLES................................................. vii I INTRODUCTION.............................................. II METHODS OF ANALYSIS....................................... 3 2.1 Waterhammer Equations............................... 3 2.2 Closed Form Solutions................................ 2.3 Method of Characteristics............................ 10 2.4 Boundary Conditions.................................. 11 2.5 Fortran IV Program for Analysis of Saturn V Propellant Systems................................... 14 2.6 Accuracy and Limitations of the Analysis............. 14 III PUMP OSCILLATION EXPERIMENT............................... 19 3.1 Purpose of the Experiment............................ 19 3.2 Experimental System.................................. 19 3.3 Steady State Measurements............................ 24 3.4 Frequency Response under Non-Cavitating Conditions... 27 3.5 Frequency Response under Cavitating Conditions....... 32 3.6 Constant Volume Joint................................ 37 3.7 Self-Excited Vibration of the System................. 37 IV PROPELLANT LINE WAVESPEED AND PUMP CAVITATION............. 46 4.1 Calculation of Wavespeed............................. 46 4.2 Wavespeed Measurements............................... 47 4.3 Effect of Cavitation Void on Wavespeed and Attenuation........................................ 49 4.4 Experiment Simulating Upstream Pump Cavitation....... 53 4.5 Model Correlating Suction Line Resonant Frequency to Suction Pressure................................. 57 V SUMMARY AND CONCLUSIONS................................... 63 REFERENCES.................................................... 66 APPENDIX A Liquid Propellant Missile Longitudinal Oscillation, R. H. Fashbaugh iii

TABLE OF CONTENTS (CONT'D) APPENDIX B Fortran IV Program for Analysis of Saturn V Propellant Systems, E. B. Wylie APPENDIX C Analysis of Cavitation Extent, P. S. Larsen iv

LIST OF ILLUSTRATIONS Figure Page 2.1 Boundary Conditions in Single Pipeline for Impedance Calculations.............................. 8 2.2 Lox Suction Piping System............................. 8 2.3 Comparison of Experimental and Calculated Results on a Single Pipe Connected to a Reciprocating Pump....... 16 2.4 Head Time Curves, Measured and Calculated, for Uniform Closure of Valve in 8.1 seconds....................... 17 3.1 Schematic Diagram of the Experimental Setup........... 20 3.2 Details of Experimental Apparatus at the Pump......... 21 3.3 Pump, Drive Mechanism, and Piping..................... 22 3.4 Overall View Showing Inlet, Exit Lines and Instrumentation...................................... 22 3.5 Head-Discharge Curve.................................. 25 3.6 Head Rise versus Suction Head......................... 26 3.7 Experimentally Determined Wavespeed in the Discharge Pipe.................................................. 28 3.8 Experimentally Determined Wavespeed in the Suction Pipe.................................................. 29 3.9 Amplitudes of Head Fluctuation versus Frequency....... 31 3.10 Pump Suction and Discharge Pressure Fluctuations...... 33 3.11 The Influence of the Suction Pressure on the Amplitude of Head Fluctuation in Suction and Discharge Side..... 36 3.12 Amplitudes of Head Fluctuation....................... 38 3.13 Constant Volume Joint................................. 39 3.14 Amplitudes of Head Fluctuation with Constant Volume Joint Installed....................................... 40 v

LIST OF ILLUSTRATIONS (CONT'D) Figure Page 3.15 Vibration of Suction Pipe and Discharge Pipe.......... 42 3.16 Self-Excited Vibration.............................. 42 4.1 Pressure at Point in Suction Line Upstream of the Pump Inducer versus Time............................. 0 4.2 Variation of Wave Propagation Velocity in Fluid in Pipeline.............................................. 50 4.3 Cavitation Apparatus Comprising L-shaped Plexiglas Tube Fitted with a Piston Driven by the Mechanical Vibrator.............................................. 54 4.4 Close-up of Plexiglas Tube with Piston. No Cavitation 55 4.5 Close-up of Plexiglas Tube with Piston. Visible Cavitation 5................................... 4.6 TITAN II, Stage I, LOX Suction Line Resonant Frequency versus Cavitating Length Lc and Average Void Fraction a.......................................... 60 4.7 TITAN II, Stage I, LOX Suction Line Resonant Frequency 61 Plate I Pressure Volume Compensating Duct (Outboard Lox).................................... facing page 13 vi

LIST OF TABLES Table Page 3.1 Data Used in Computer Program...................... 32 4.1 Experimental Wavespeeds.............................. 53 4.2 Experimental & Computational Results with Cavitation............................................ 57 vii

I. INTRODUCTION Longitudinal oscillations have been found to occur in a number of liquid-fueled space vehicles, both manned and unmanned. Although the general sequence of events that leads to the self-excited oscillation of the vehicle (POGO stick) is reasonably well understood and documented, analytical methods that predict the entire behavior of the vehicle are not available. The purpose of this study was to investigate a portion of the total vehicle dynamics, namely the propellant feed lines. The dynamics of the Saturn V suction lines were to be analyzed, a Fortran programmed analysis be provided, and an understanding of'the mechanics of unsteady flow in the system be conveyed to the Contracting Officer or his representatives. During the early stages of the Contract these objectives were enlarged to include the analysis of the entire propellant feed system from the feed tank to the engine, including a study of the action of the pump in a dynamic situation. Two experimental programs were proposed and carried through; the first was aimed at providing a better understanding of the extent of cavitation in a pipeline containing a cavitation source and to obtain a measurement of wave speeds in the cavitating zone; the second was to provide some confirmation of the method of analytically modeling the effect of an oscillating pump in a system. All aspects of the original and enlarged objectives have been met and the results are reported herein. Separate reports have been -1

-2 issued on two important phases of the project, and these are included as appendices in this report. The first of these, Appendix A, attempts to convey a basic understanding of the entire problem by discussing the role that each of the basic missile components plays in the development of an oscillation. The second, Appendix B, presents the Fortran computer program to be used for the analysis of the propellant system, together with a description of the program details and the theory underlying the method of analysis. This report first discusses the various methods of analysis with emphasis placed on the computer solution using the method of characteristics. Examples are given showing comparisons between experimental and analytical results. These provide an indication of the accuracy of the method. At the same time the limitations of the method of analysis are clearly set forth. The experimental work in connection with the oscillating pump proved to be very enlightening. It provided a great deal of information pertaining to the accuracy of the method of representing this boundary condition in the program, demonstrated the influence of the P.V.C. on the system, and provided data on the dynamics of the pump during a cavitating condition. These results are reported in detail herein. The cavitating condition in the suction line caused by the pump inducer provides one of the most uncertain elements in the entire analysis. This phenomena and its effect upon wave speed and attenuation are discussed and the experimental work in connection with cavitation are also reported. In addition to the information presented in Chapter IV, a separate description of the various aspects of this topic is presented in Appendix C.

II. METHODS OF ANALYSIS 2.1 Waterhammer Equations The analysis of transients in a complex fluid system such as the Saturn V propellant lines involves the solution of the partial differential equations that describe the flow in each pipeline in conjunction with the appropriate boundary conditions that connect the component parts of the system. A one-dimensional analysis of flow in a pipeline is used with the various characteristics such as friction, compliance, and inertance treated as distributed parameters. In some situations lumped parameters adequately describe the dynamic action of a component part of the system, and some simplification may be realized. The equations of waterhammer written for compressible liquids in elastic pipelines are developed in Appendix B. Since the pressure pulse wavespeed a is very great compared with the fluid velocity V, the general one-dimensional equations of motion and continuity may be simplified. Vt + - + VIVI - g sin a = 0 (2.1.1) p 2D Pt 2 - + av = 0 (2.1.2) P In these equations t and x represent the independent variables time and distance along the pipe, and when used as subscripts, they denote partial differentiation. The density of the fluid is given by p and the pressure is given by p. The angle a is measured between a plane normal to the flight direction and the positive x axis, g is the -3

-4 acceleration at the particular time in flight, D is the inside pipe diameter, and f is the Darcy-Weisbach friction factoro The pressure pulse wavespeed is denoted by a, and is defined by a _/p (2l.13) 1 + KD Ee The bulk modulus of elasticity of the liquid is given by K, the modulus of elasticity of the pipe wall by E, and the pipe wall thickness by e. In the equation of motion the only energy loss is that resulting from the shear stress at the pipe wall, and it is incorporated into the equation by use of the Darcy-Weisbach friction factoro Equation (2.11) includes the nonlinear term for turbulent flow with a velocity exponent of two. A different exponent could be used equally well ( Generally f is considered a constant during the transient calculations although this is not necessaryo Unsteady flow of compressible fluids in elastic pipelines, described by the above equations, can be conveniently handled in either of two different ways, depending upon the type of excitation and the type of results desired from the analysiso The first treats the steadyoscillatory case in which the same cycle is repeated periodically. A closed form solution, referred to as the impedance solution, is used for these studies. The second method, which is used for the Saturn V system analysis, is most useful when a transient condition develops in a system and continues until either another steady-state or steadyoscillatory condition exists. The method of characteristics solution is found to be most convenient since it is able to trace conditions throughout the entire system during the development of the transiento These two approaches are briefly treated in the next two sections.

-52.2 Closed Form Solutions An impedance analysis provides a convenient means of determining the resonant frequencies of complex distributed-parameter fluid systems. In addition, system characteristics such as transfer functions, hydraulic impedance, frequency response, etc., over a wide frequency range, can be economically evaluated for specific system configurations. Other methods of transient-flow analysis of a particular system yield only the system response to a particular excitation at one frequency. The impedance theory assumes the existence of steady-oscillatingsinusoidal waves in the fluid piping system. In the ideal frictionless single pipeline or series system the waves appear as standing waves superposed upon a steady flow condition (which may be zero). In real viscous fluid systems the standing waves appear in combination with traveling waves. The differential equations are taken in simplified form, Eqs. (2.1.1) and (2.1.2) with the slope term removed, average flow relations are removed, the friction term is linearized, and resulting equations are solved for a sine-wave disturbance.() The algebraic equations which result are then used in combination with the boundary conditions to effect a solution. In order to handle a complex system the algebraic equations of one pipe are written in terms of the parameters of that pipe and in terms of its terminating impedance. At a series connection the terminal impedance of one pipe becomes the input impedance of the next pipe. By beginning at a point of known hydraulic impedance and proceeding systematically through the system the hydraulic impedance can be computed at any location. In general a large amplitude hydraulic impedance at a section

-6 in a system indicates that a small discharge fluctuation will produce a large pressure head fluctuation. Thus a look at the variation of the impedance with frequency at a point in a system will indicate the possibility of a severe pressure oscillation, assuming the existence of an exciter. A computer program can be written which will scan a broad frequency range evaluating the impedance at any point in a system. This is the information needed for the determination of the resonant periods of the system. The resonant frequencies of a single pipeline can be calculated directly from the theoretical period, 4L/a. For a single pipeline leading from a constant pressure source the fundamental period is equal to the theoretical period. The higher harmonics in the system are related to the theoretical period by integers. In complex systems the resonant frequencies do not bear these same relationships with the theoretical period. In a series system the fundamental period is at the apparent period which is not Z 4L/a. Also the higher harmonics are not, in general, related by integers to either the theoretical or the apparent period. A series system may be made up of pipes of different diameters connected in series, or a single pipeline of constant cross sectional area but with different wavespeeds or frictional properties at various sections throughout its length. Partial reflections at the junctions of these sections give rise to the apparent period. Since the impedance theory takes into account all forward and reflected waves in the system, it yields the correct resonating periods of the complex system.

-7 Impedance Equations. The hydraulic impedance at a section in a pipeline is defined as the ratio of the oscillatory pressure head to the oscillatory discharge, Zx = Hx/x. Each of these variables is a complex number since each has an amplitude and phase angle associated with it. The solutions for the oscillatory head and discharge in a single pipeline, Fig. 2.1, in terms of the boundary conditions and the pipeline parameters, are as follows: Hx = HR cosh Yx - QR Zc sinh yx (2.2.1) HR Q = sinh yx + QR cosh yx (2.2.2) Zc HR and QR are the complex numbers representing the oscillatory head and discharge at x = 0. The complex constant 7 is given by Y _-2 ( - - + i R) (2.2.3) a2 gA where A is the pipe cross sectional area, g is the acceleration due to gravity, c is the circular frequency, a is the wave velocity in the fluid in the pipeline, and R is the linearized friction constant given, for turbulent flow, by R = (2.2.4) g D A2 In this equation f is the Darcy-Weisbach friction factor, Q is the mean flow in the pipeline in cubic feet per second, and D is the pipe diameter. In Eqs. (2.2.1) and (2.2.2) Zc is called the characteristic impedance of the pipeline, and has the value

-8 F L HR R IL Figure 2.1 Boundary Conditions in Single Pipeline for Impedance Calculations. I LI + L? = 11.2 5" a Lox Tank Pump Pipe 3 ZRo D=20" I L3 =46.2' Pipe 2 D=17"' " Pipe 1 1.25' L210015"s D= 17' o3 =22800 ft/s -3 O-280 200~o3s2800 a3= 2800 ft/s Figure 2.2 Lox Suction Piping System.

-9 Zc = (2.2.5) i wg A Since the hydraulic impedance at any point in a pipeline is given by the ratio of the oscillatory head to discharge, the hydraulic impedance at any point x is Hx ZR-Zc tanh 7x (2.2.6) Zx= --------- (2.2.6) Qx 1 - R tanh yx Zc where ZR = HR/QR. As an example of the application of the impedance method to determine the resonating frequencies of a complex system, the Saturn V Lox suction piping system is considered, Fig. 2.2. The lox supply tank is assumed to behave as an infinite reservoir with a constant pressure. With HR = 0, the impedance at the tank is zero, ZR = 0. The impedance at any point in pipe 3 can be determined using Eq. (2.2.6) with ZR = 0. Zx = - Zc3 tanh 73x (2.2.7) In particular, when x = L3, the impedance at the junction of pipes 2 and 3 can be determined. At the junction, a common pressure exists at any instant and continuity of the flow rate must be satisfied. Therefore the hydraulic impedance at the right end of pipe 3 becomes ZR for the left end of pipe 2. Equation (2.2.6) can now be applied to pipe 2 and the hydraulic impedance can be computed at the right end of pipe 2 (Zc and 7 now pertain to pipe 2). The same procedure is repeated for pipe 1 yielding the hydraulic impedance at the pump. This analysis is set up in the computer program and an impedance evaluation is made

-10 for a number of different frequencies. An examination of these results gives the resonating periods for the particular data analyzed. A change in the data, that is, L2, L3 and a3 will yield a different resonating frequency. In some applications of the impedance theory it is desirable to use ratios different from the head-discharge ratio at a single point in the system. These ratios are called transfer functions and may relate any of the variables in the system, such as the pressure-head variation at one point in the system to the pressure-head variation at another location. These can be computed(l) and the results of a computer program will show the variation of the transfer function over a broad frequency range. 2.3 Method of Characteristics The theory used to solve the partial differential equations in the program for Saturn V is based on the method of characteristics. It leads to a straightforward numerical solution of the differential equations and provides relatively simple relationships which are combined with boundary conditions to provide the transient solution for a complete system. Nonlinear characteristics such as friction losses, pump characteristics, flow through orifices, etc., are treated in their real form without linearization. The method is developed in Appendix B. The partial differential equations, Eqs. (2.1.1) and (2.1.2), are transformed into particular total differential equations by use of the method of characteristics. These equations are then integrated to produce a form suitable for numerical solution on the digital computer. The transient conditions

-11 in each pipeline of a complex system are treated using this procedure. Boundary conditions are used to transfer the effect of one pipeline to another, and to introduce the effect of the pump, the P.V.C., the engine, etc. Lumped parameters can be used where appropriate. In a typical computer analysis of a specific problem, an initial steady-state operating condition is imposed upon the system; then a numerical solution is obtained at finite points throughout the entire system at finite time increments during the life of the transient. The results from the execution of a characteristics method program for a particular system show the complete system response to a particular excitation. This may be contrasted with the results of an impedance program on the same system, which would yield the response at particular points in the system resulting from a periodic input motion at a number of different frequencies. The former method can produce the same results as the latter and is not subject to the introduction of possible errors due to linearization. The characteristics method is, however, more costly if it is necessary to investigate system excitation at a number of different frequencies. 2.4 Boundary Conditions In either of the foregoing analyses it can be observed that the unsteady motion in individual pipelines is treated uniquely, and that boundary conditions are required to transmit the effect of one pipeline to the other, and also to introduce the influence of other component parts of the entire system. The treatment of boundary conditions is well documented in the literature(l) and the particular conditions pertaining to Saturn V

-12 are included in Appendix Bo Therefore, only a brief discussion of some of the elements peculiar to Saturn V is included hereino Minor energy losses at changes in cross-section, elbows9 valves2 etc.9 have not been included in the program. Although it is possible to include these it was felt that the same overall effect could be obtained by adjusting the friction factor in the adjacent pipelineso This distributes the concentrated loss along the pipelineo The pump and its dynamic properties are probably the most important elements in the entire system2 and also the most uncertain due to the limited specific knowledge of the pumpts actionso It is assumed that the steady state characteristic curves are valid, and also that the shape of the curves remains the same during the transient condition. Thus the shape of the head rise versus suction pressure relationship is the same for any flow rate near the operating condition2 and the shape of the head-discharge curve is valid for any suction pressure and pump speed near the operating zoneo The pump is considered as a pipe with a length equal to the mean flow path length through the pump. In this particular pipe the same fluid characteristic equations can be used except the friction term is replaced by the pump pressure rise term. A pressure pulse M-avespeed is assumed through the pump which gives rise to a small time delayo The pump is considered to be one of the possible exciters for the system since a structural motion could cause a physical displacement of the pump with respect to the suction lines and thereby introduce a volumetric change in the system. This excitation is modeled

Plate I

-13 as a periodic discharge input at the pump inlet. Some questions arise as to the ability of the pump motion to excite the system. The experimental program discussed in Chapter III treats this topic in some detail. The pressure volume compensating duct is designed to provide a constant volume in the suction system during a structural oscillation. Thus, a decrease in volume in the suction line above the pump caused by the physical motion of the pump towards the tank results in an increase in volume in the P.V.C.; the net result is that the pipeline above the P.V.C. is not significantly influenced. For pump motion in the opposite direction, a reverse action takes place. The manner in which this condition was programmed is shown in Appendix B and an example of computer results is given to show the stabilizing influence of the P.V.C. A constant volume joint was used in the experimental program discussed in Chapter III and it, too, showed the stabilizing effect (see Plate I). A structural motion of the feed tank was also considered as a possible exciter to the system. This exciter is modeled the same as the pump excitation, that is, the junction is held stationary and a side branch piston is assumed which produces a volume'tric discharge oscillation to yield an equivalent effect. There is some uncertainty at this boundary condition as to just how much of the total feed tank area contributes additional flow during an oscillation. At the engine the flow from the pump discharge line ends in a chamber or manifold before passing through the burner orifice plates to the combustion chamber. The compliance of this chamber is taken into account by lumping its elasticity effects at this point. The burner plate orifices are treated using the nonlinear orifice equation.

-14 A structural excitation is also assumed to be possible at the engine orifice plate; however, this appears to be the least likely point of the structural inputs considered, due to the compliance available in the lox dome and the fuel manifold, and also due to the relatively high losses in the pump discharge piping. 2.5 Fortran IV Program for Analysis of Saturn V Propellant Systems The characteristics method digital computer program for analysis of the propellant systems of Saturn V is presented and discussed in Appendix B. It is written in Fortran IV-G Level and was executed on an IBM 360/67. Complete Saturn V data were not used in clearing the program as portions of the data were not available, particularly pertaining to the discharge piping system and the pump characteristics. The program was checked out using data from Titan II and favorable comparisons were obtained between computed results and experimental results from the literature. (2) An example is given that shows the transients that develop in the lox system as a result of a forced pump motion. The effect of the addition of a P.V.C. in the system is demonstrated by a second example which has the same system excitation but the action of the P.V.C. is included. 2.6 Accuracy and Limitations of the Analysis The best available evidence to demonstrate the accuracy and reliability of an analysis based on the characteristics method is the presentation of experimental results with the corresponding computed data.

-15 The first example was obtained by use of a reciprocating pump. The Triplex pump caused periodic fluctuations in the system which were experimentally recorded and compared with computed data, Fig. 2.3. The test pipeline was a 3-in.-diameter horizontal pipe, 56 feet in length, connected to the suction flange of the pump. The pipe was drawing water from a constant head reservoir which provided the other boundary condition for the system. Pressures shown are at the pipeline connection to the pump. A second example results from a uniform valve closure in a 4,000 feet 1-in.-diameter pipeline connected to a constant head reservoir. The pressures were recorded at the valve, Fig. 2.4. Other examples are available in the literature.(1,2) Chapter III shows further evidence of a good correlation between experiment and theory. As long as the underlying assumptions or restrictions imposed upon the basic equations are kept in mind, it is felt that the proposed solution is accurate. The most important limitations are listed. The equations are based upon the assumption of one dimensional flow. In most systems this does not cause any problems. However, if the length to diameter ratio becomes very small in pipelines and if, at the same time, there are many boundary changes, like elbows, tees, valves, etc., this restriction may become important. The equations assume that the pipeline remains full at all times. If a cavity develops in a system as a result of vapor pressure being established for some duration, then a discontinuity exists and a separate boundary condition must be written to properly describe the

w 60 / \ 40 __ -:o 0 PUMP PERIODS (.2085 SEC) Figure 2.3 Comparison of Experimental and Calculated Results on a Single Pipe g 20Pum 0 t EXPER MENTAL -20 Pe 75 PSI -40 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 PUMP PERIODS (=.2085 SEC) Figure 2.3 Comparison of Experimental and Calculated Results on a Single Pipe Connected to a Reciprocating Pump.

-17 W 80 //!a /F 0LL~ / L EXPERI ENTAL / L/ 0 0 2 4 6 8 10 12 14 16 TIME, Se Figure 2.4 Head Time Curves, Measured and Calculated, for Uniform Closure of Valve in 8.1 sec.

-18 system response. Frictional losses in unsteady flow are assumed to be the same as the losses in steady flow at the same velocityo It is generally known that this is not a completely accurate assumption but, with a few exceptions, reasonable results are usually obtained. High frequency components dampen out much more rapidly than low frequency so energy losses in a high frequency system are not too well described by the steady flow equivalent losso An adjusted friction factor can often be used to compensate for the additional dampingo Systems that are subjected to vapor pressure and possible cavity formation, and also those that may contain undissolved gases are also somewhat difficult to accurately analyze. The wavespeed is a strong function of gas content- in the liquid and also dissipation increases sharply in a liquid-gaseous mixtureo This aspect is discussed in Chapter IVo

III. PUMP OSCILLATION EXPERIMENT 3.1 Purpose of the Experiment While there exists sufficient proof for the reliability in the mathematical treatment of the dynamics of fluid lines, experimental verification was believed to be necessary for the complicated boundary condition at the pump. Simplifications were proposed in the mathematical modeling of the vibrating pump and these needed to be checked. The replacement of the oscillating pump by a fixed pump with a periodic volume input, and the use of steady state pump characteristics during unsteady flow when the pump is operating in the cavitating zone were specific uncertainties. Moreover, the theoretical appearance of a dynamic gain through the pump under cavitating operation seemed to be of such significance for the project that experimental proof was desirable. These initial objectives were exceeded during the experimental work, by demonstrating the effectiveness of a constant volume slip joint, and by establishing self-excited vibrations in the system. The latter were caused by an interaction of structural and hydraulic systems and exhibited significant similarity to the so-called Pogo instability of rockets. 3.2 Experimental System The experimental apparatus was located in the G. G. Brown Fluids Engineering Laboratory of The University of Michigan. Figures 3.1 to 3.4 show the basic components of the experimental setup, which is a closed-loop system consisting of tank, suction line, oscillating -19

-20 Second floor Suction line 71.6' I /2 I.D. Discharge line 74.7' I 1/4 I.D.'I I Air 24.5' Pressure Regulator Basement Figure 3.1 Schematic Diagram of the Experimental Setup.

Swing joints F I 1/4 I. D. Plexiglos L-1.97'"_ _Slip joint H Centrifugal \; -""=^) _~ )pump o Displacement transducer I, I /2 I.D. Plexiglos- apparOs ati L- 4.45 Figure 3.2 Details of Experimental Apparatus at the Pump.

-22 Figure 3.3 Pump, Drive Mechanism, and Piping. i:"aeaaBers s*ag Figure 3.4 Overall View Showing Inlet, Exit Lines and Instrumentation.

-23 pump, and discharge line return to the same tank. The tank is placed in the basement of the laboratory and the tank pressure can be adjusted by a pressurized air supply system to pressures between 0 and 80 psi. The suction line consists of a 1-1/2" steel pipe, 67 feet in length, and a piece of 1-1/2" plexiglas pipe, 4.45 feet in length, which is mounted in front of the pump to permit observance of cavitation bubbles. The discharge line is 1-1/4" steel pipe, 71.50 feet in length and a 2 foot piece of 1-1/4" plexiglas pipe. The centrifugal pump is placed on the second floor of the building and is mounted in such a way as to oscillate in a sinusoidal motion along the axis of the pipes, driven by a pneumatic motor. The double amplitude of the motion is adjustable between 0 and 1 inch and the frequency between 0 and 15 cps. The connection between the moving pump and the fixed pipes is a slip joint on the suction side, thus causing a volume displacement, and swing joints on the discharge side, which have a negligible dynamic input. Provisions were made to take dynamic pressure measurements with strain gauge type transducers at several locations in the system, and static pressures with bourdon gauges and mercury manometers. An orifice flow meter in the discharge line served to determine the steady state discharge, the pump motion was recorded with a displacement transducer and the frequency was measured with an electronic counter. List of Instrumentation: Storage Oscilloscope, Tektronix 564 Optical Recorder, Sanborn 860-4000D Strain Gauge Transducers, Statham, 0-100 psi

-24 Crystal Transducers, Kistler 601 Electronic counter with magnetic pick up, Hewlett Packard Displacement Transducer, Sanborn, 24-DCDT-1000 Centrifugal Pump, 9Sta-rite model JFH 3-51 3o3 Steady State Measurements ao The Steady State Characteristics of the Pump The steady state pump characteristics were assumed to hold for the analysis of the dynamical behavior of the system. They were measured with mercury manometers on both sides of the pump and an orifice flow meter in the discharge lineo The head rise versus discharge curve and the head rise versus suction head curve are presented in Figso 3 5 and 3 60 They are applied in the computer program in a simplified mannero It is assumed that all head rise versus suction head curves can be obtained from the curve for Q = 27.4 gpm by simply shifting it up or down. The amount of shift is obtained from the head rise versus discharge curve given for a suction head of Hs = 50 feet, in Figo 35o. To prove that the suction pressure does not significantly affect the shape of the curve in the region of high discharge, a second curve is shown in Fig. 3~5 for a cavitating pump operation. All curves are related to the suction head instead of the NPSH, since after a sufficiently long running time the water temperature stabilizes on a reproducible level. b -Wavespeed Measurements The propagation velocity of the pressure waves are of fundamental importance for the prediction of the resonating frequencies of

150 1 1 Suction hed =50 feet H20 00 - I I 100 —- | — - ~ t - 0".T —-- - 100 10 20 30 405 I Suction heod I 0 * I I -25.2 feet H201 O *0 I <u OR - 1 0 10 20 30 40 50 Q (gol/min) Figure 3.5 Head-Discharge Curve. I ro Vi

-26 100 0 Iq F 70 LLu z Q:: I 0o 40 00 00 0 0 0 roro rd~ f0 rd ro I ii I I I I I | I I II 30 1 II I I -30 -20 -10 0 10 SUCTION HEAD IN FEET H20 20 30 Figure 3.6 Head Rise versus Suction Head.

-27 piping systems. Equations to calculate the wavespeed are given in Chapter IV, however, since the exact amount of air in the water was not predictable, measurements were carried through in the suction and the discharge pipe. Different experimental approaches are possible, which are discussed in more detail in Chapter IV. The method applied here consisted of generating a pressure pulse in the piping system and measuring the time to traverse the distance between two pressure transducer stations. The test section in the discharge side was 47 feet long and included most of the horizontal pipe section as well as part of the veritical section in the first floor of the laboratory. Fig. 3.7 displays the results as a function of the steady state pressure, which is measured on the downstream side of the flow meter orifice. The test section in the suction side consisted of 20 feet of the horizontal part in the first floor. The pressure to which the wavespeed is related in Fig. 3.8 is measured on the suction side of the pump. When the pressure on the suction side is less than -20 feet of H20 air bubbles appear in the water and the wavespeed decreases rapidly. The same effect can be observed in the discharge side under a correspondingly higher pressure, since the pressure in the discharge pipe is always approximately 40 feet of H20 higher than in the suction side. However, Figs. 3.7 and 3.8 are not quite compatible since they were taken a few weeks apart. In the meantime precautions were made in the system to make the suction side airtight, thus reducing the amount of air in the water. 3.4 Frequency Response under Non-Cavitating Conditions To obtain a set of data the pressure in the tank is set at a

5000 * 0 0* 0 4000 (:3 z 3000 0 LLI a 2000 LJ 3 1000 0 0 F — g,. I ro 0o i __ _ _ _ _ _ _ _ _ _ _ _ I _ _ _ _ _ _ _ _ _ _ _ __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 20 40 PRESSURE 60 IN FEET 80 100 H20 Figure 3.7 Experimentally Determined Wavespeed in the Discharge Pipe.

5000 4z C n LLJ QLLJ <( 4000 300C 2000 o oo o o0 0 0' 00 -0o I 11 I ra) \o I 1000 0 30 - 20 -10 0 10 20 HEAD IN FEET H20 Figure 3.8 Experimentally Determined Wavespeed in the Suction Pipe.

-30 particular level, the pump is moved at a specific frequency' and the resulting pressure amplitudes are recorded. The frequency is then changed and new pressures are recordedo As the frequency is varied over the experimental range the amplitude of the excitation produced by the pump motion varies slightly. In order to place the results on a comparable basis and to account for this variation the experimental results were multiplied by a normalization factor (actual double amplitude/o070). The value.070" was the double amplitude of the exciter at a frequency of 1 cycle per second. The adjusted experimental pressure amplitudes are given in Fig. 3~9 as a function of the exciting frequency. Computed results are also shown in Figo 309 and the agreement between predicted and experimental values is seen to be reasonably good. The theoretical results are based on the impedance method outlined in Chapter IIo Since the mean pressure varies throughout the system the wavespeed for each pipe section had to be chosen dependent upon the mean pressure in that reacho The friction factor of.1 which was used is about five times higher than the steady state valueo The fact that the energy dissipation in a high frequency flow is larger than in steady state flow is experimentally well established but no literature giving quantitative results is available. The value of.1 gives reasonable amplitudes for the resonating case in this system. The resonant frequency is not affected by the friction factor. Experimental and theoretical results show that the oscillating pump can excite the fluid system and establish very high pressure fluctuations under resonating conditions. This occurred at 50-60

-31 0 N CJ II U0 QD LLJ I UiLi 0CL 0 10 20o 30_40 50_ 60 0 9I 1 0 10 20 30 40 50 60 70 80 90 100 CIRCULAR FREQUENCY (rod i ons per second) Figure 3.9 Amplitudes of Head Fluctuation versus Frequency.

-32 radians per second in this example caseo It should be noted from Fig. 309 that the pressure amplitudes on the suction side are equal to pressure amplitudes on the discharge side. Table 3.1 Data Used in the Computer Program Pipe No. Length i 61.7 2 4.1 3 4o45 4 1.3 5 Pump 6 8.3 7.87 8 3.39 9 60.8 Mean Suction Head Mean Discharge 39 Diameter 1.61 1.61 1.61 1.61 (I.D.) Waves peed 2500 2500 1500 2500 Friction Factor.1.1.1.. 1. 1 ol i.1 1.38 1038 1.38 1.38 4000 2000 4ooo 4ooo 4000 4000 46' H20 gpm 3o5 Frequency Response of the System under Cavitating Conditions An intensive study was carried through to gain insight in the dynamics of the experimental system when the pump is cavitatingo In a first series of experiments the pump was shaken with a constant frequency of 3 cycles per second and a double amplitude of motion of.237 inches9 while the tank pressure and thus the mean pressure throughout the system was variedo In Fig. 3,10, a-h the resulting head fluctuations

7Fig. 3.! a Hs =- 28.4 H = 6.i, Q = 342..n. 1~19 3.10,11 s 27 2i!~D = 3205, q 35 U) Oo F 6i3!oc Hs= - 23-9, = 50.2 3 Fg. 3. s - 2!.:: 57.0 40.g...e 3.- P Sucton man Dscharge essure u>ct uat.:ons

Fig. 3 i e. 6 Q -4!z 7..,l9, 3. " Of' H = 7.1 ~, H = 5. 61 q = 4i I+' t -IIi A""g. 3.-Og rig. 3.i"g = 25.7, cD =,-. 4. 3.10h s =5.9, - 88.53, Q = 4A:O 33:o a.igge 3. 0 ('2'= ).r r Zi 25 I *S, D= I-C0-Oq

-35 on the discharge side of the pump, upper curve, and on the suction side of the pump, lower curve, are recorded. The head scale is 1 unit = 11.9 feet H20 and the time scale is 1 unit =.1 sec. The time average values of suction head, Hs, discharge head, HD, and discharge, Q, are given for each picture in feet H20 and gpm, and, to show the points of operation on the head rise versus suction pressure curve, the mean suction pressures are marked in Fig. 3.6. In Fig. 3.10a the head fluctuations on both sides of the pump are very small since the pump is so strongly cavitating that no significant head rise variation is obtained in spite of the considerable discharge fluctuation caused by the moving pump. In Figs. 3.10, b, c and d the point of operation is on the descending portion of the head rise versus suction pressure characteristic of the pump, and a dynamic gain is obtained. The suction pressure is almost constant vapor pressure in Fig. 3.10b and hits vapor pressure in portions of the cycle in Figs. 3.10c and 3.10d. With increasing suction head, wave form and amplitude become almost identical on suction and discharge side, Figs. 3.10e-h. The influence of the suction pressure on the dynamic gain through the pump can be observed also in Figs. 3.11a and b. The mean suction pressure is continuously increased starting from vapor pressure and the amplitude of the oscillations of suction head and discharge head are recorded. The frequency of the exciting pump motion is 3 cps in Fig. 3.11a and 6 cps in Fig. 3.1lb. In another series of measurements the suction pressure at the pump was kept constant on Hs = -26.2 feet H20 and the frequency of the pump motion was varied between 2 and 8 cps. The amplitude of the head

MrMean Head on Discharge Side \ < Head Fluctuation on Discharge Side / Heaa FLuctuation on Suction Siae VMean Hea& on Sucotion Side Fig. 3. a Freuency: 3 cps Fig. 3- lc ixplanat-ir f F ig. 3.2la ana 3.ib.I ~.~..~.~ ~.v~ o ~ 3" a n 3il OI' Figure 3.11 The Influence of the Suction Presszre on he Am_.itzue of Hepa Fuctuationn S uction o an Discharge Side. Head Scae: 1 unit = 11.8 ft of H20 Fig. 3. 1b Frequency: 6 cps

-37 fluctuation on the suction and the discharge side were measured and are plotted in Fig. 3.12e A dynamic gain occurred throughout the whole frequency range varying between 1.5 and 5o The fundamental. frequency of the system is approximately 2 cps in this case, due to the low wavespeed, and a higher harmonic causes high pressure fluctuations at about 8 cps. It should be noted that the system behaves in a highly nonlinear manner, partly due to the dependency of the wavespeed upon the pressure and due to the pump characteristics of the cavitating pump. Distortions of the wave forms are further caused by lateral vibrations of the piping system, which could not be totally excluded during the experiments. 3.6 Constant Volume Joint To demonstrate the effectiveness of a constant volume joint the device shown in Fig. 3.13 was installed into the system replacing the slip joint, which was used normally. The magnitude of the pressure oscillation resulting from the pump motion was observed under conditions equivalent to the experiment described in Chapter III, Section 4~ A comparison of the results in Fig. 3o14 with Fig. 3.9 stresses the substantial reduction of pressure fluctuations observed by use of the device. 3.7 Self-excited Vibration of the System All vibrations described hitherto were forced by moving the pumpo However during the experimental work, self-excited vibrations were observed in the system which continued without external excitation, The oscillation bears a strong resemblance to the POGO instability in that it is the result of an instability encompassing the entire structural

-38 5 4 z 3 2 I 0 N ILL bJ LL z 0 I LL I C) LLJ L0 LJ 0 D 0LJ OD 0 0 30 25 20 15 10 5 0 X-x -- ^ -------—. X |_X'X x x X X x x X 0 0 0 00 * 0 ~ ~ 0 *b o I~ I t - * ~ ___. ____ __L___I - pX ~~~..1- d 2.0 3.0 4.0 5.0 FREQUENCY 6.0 IN cps 7.0 8.0 Figure 3.12 Amplitudes of Head Fluctuation.

SLIDING BEARING CON CONSTANT VOLUME JOINT full size Figure 3.13 Constant Volume Joint

og~ 20 —------- I I!- I I I- I I o 20 0 I ~-J0 o Ie a. ~ L 0 10 20 30 40 50 60 70 80 90 100 CIRCULAR FREQUENCY (RADIANS PER SECOND) Figure 3.14 Amplitudes of Head Fluctuation with Constant Volume.

-41 and fluid system. The instability occurred when the mean suction head was between -27 and -29 feet of H120 and the pump was operating on the curved portion of the head rise versus suction pressure curve. The explanation for this phenomenon can be found in the interaction between piping structure, fluid system, and cavitating pump. As sketched in Fig. 3.15 the horizontal section of the discharge pipe and the suction pipe can move in an axial direction as a mass-spring system with relatively small damping. Both pipes are tied firmly together and execute the same movement. The driving force for the structural vibration is the dynamic head difference between the left end and the right end of the horizontal section. The following sequence of events outlines the growth of the instability. Some initial disturbance creates a pressure pulse in the suction side; this is amplified as it passes through the cavitating pump; it is transmitted into the discharge pipeline where a pressure condition is established that results in a movement of the piping system; the structural movement initiates a new pressure fluctuation in the suction side, thus completing the loop. Although the perturbations are nonsinusoidal in form, a clearer insight into the instability is possible if sinusoidal displacements are assumed, as well as sinusoidal pressure and flow oscillations. Figures 3.10a to 3.10d and Fig. 3.11 show clearly the strong pressure amplification at the pump under certain suction pressure conditions. Under these operating conditions a minor periodic fluctuation in the suction side of the pump will produce a relatively high pressure fluctuation at the pump discharge. If this pressure oscillation is periodic, a standing

DISCHARGE SUCTION PIPE Fi.gure 3.15 Vibration of Suction Pipe and Discharge Pipe. Fig. 3L.6a Self-excited Vibration Upper Curve: Itead on Discharge Side Lower Curve: Head on Suction Side Mean Suction Head - 28.3 ft Ii20 Mean Discharge Head: 1.2.5 ft H120 Time Scale: I unit.2 sec Head Scale:a I unit 5.9 ft HIpO Fig. 3. 16b Self-excited Vibration Upper Curve: Htead on Discharge Side Lower Curve: Head on Suction Side Mean Suction Head: - 27.8 ft H20 Mean Discharge Head: 16.0 ft H20 Time Scale: I unit.^, sec Head Scale: I unit 11. 9 ft Hp2

wave may be created in the discharge piping system which will exhibit significant pressure variations along the discharge pipeline. A further requisite of this standing wave is that the fluid system must respond in such a way as to not nulify the gain produced at the pump. That is, the phase of the reflected wave returning to the pump should add to the incident wave at that point. This will happen, and significant standing waves can be created, at or near the fundamental period or one of the odd harmonics of the discharge piping system. If a standing wave of the fundamental period is assumed in the discharge pipeline, with nodal point at the supply tank and pressure loop near the pump, a periodically-varying pressure field can be visualized at the left and right end of the horizontal section of discharge pipe. These pressures produce the driving force to move this portion of the structure in a periodic manner. The horizontal portion of the discharge and suction piping which are fastened together, as well as the pipe hangers and adjacent connected vertical pipes can be thought of as a spring-mass system. The requisite requirement for this spring-mass system to respond is that one of its natural periods must be near the period of the driving force. When the suction pipeline is moved a pressure perturbation is created at the period of the structural motion. Again, in the suction pipeline, a standing wave can be created which produces a periodic pulse at the pump suction flange. With this complete feed-back loop, the oscillations will continue to develop in magnitude until the dissipative mechanisms within the system just balance the energy being added, or until the magnitudes become sufficiently large that

non-linearities become important and the oscillatory pattern becomes irregular. Figures 3.16a and 3.16b show recordings of the pressure on the suction and discharge side of the pump during two different cases of this self-excited oscillation. The recordings illustrate the high amplification at the cavitating pump. These results show a period of about 0.27 sec which corresponds to the observed natural period of the suspended pipeline system. If the fluid in the discharge pipeline were responding at the fndamental, this period would require an effective wavespeed of about 1100 ft/sec. The third harmonic in the suction system would require an effective wavespeed of 350 ft/sec. Although both of these velocities appear low they are reasonable in the low pressure region in which the system is operating. Furthermore it is not necessary that these harmonic periods be matched exactly. These numbers probably indicate a lower limit of the effective wavespeed in each pipeline. Two important aspects in the auto-oscillation should be recognized. The system can be unstable only when the pump provides sufficient dynamic gain to overcome the fluid friction and the structural damping, and the frequency of vibration is dependent upon the natural frequency of the structure and the response of the fluid system; it is not related to the pump speed. If any one of three major elements of the system (the piping structure, the fluid system, or the cavitating pump) is altered sufficiently, the instability will be destroyed. Definite changes were made in two of these elements to prove this point. (a) The natural period of the structure was altered by

-45 changing the manner in which the pipeline was suspended. The self-excited instability could not be created under these conditions, even by creating fairly large initial pulses. The initial waves quickly dissipated without further response in the system. (b) A strongly oscillating system could be brought to rest merely by either raising or lowering the suction head, thereby decreasing the dynamic gain across the pump. (c) The third element of the system was not so easily varied, however, some indication was obtained that it is equally important to the total instability. It is known that the wavespeed varies greatly with air content in the liquid. In general the self-excited instability could not be created on initial start-up of the equipment. That is, it was only possible to create the auto-oscillation after the closed system had been in operation for some time, and the fluid had the opportunity to reach a certain semi-stable operating condition. On initial start-up the fluid properties were not quite correct to yield a standing wave pattern that would force the structure at or near its natural period. The parallel between this auto-oscillation and the POGO instability is fairly direct. The latter resulted from a closed loop interaction between the engine thrust force, the structure, and the propellant feed system including the pumps. In the present case the discharge piping system is able to provide the driving force to close the loop and cause the instability.

IV. PROPELLANT LINE WAVESPEED AND PUMP CAVITATION 4.1 Calculation of Wavespeed The speed of a pressure pulse wave in an elastic pipeline filled with a liquid is computed using the equation( a.= K 7 (- (4.1.1) J1 + (K/E)(D/e) The numerator is the speed of a pulse wave in an infinite fluid; this speed is reduced by the elasticity of the pipe wall. The propagation velocity of a pressure wave in a pipeline containing a liquid may be reduced significantly if undissolved gas bubbles are dispersed throughout the liquid. A direct calculation shows that even a small quantity of gas greatly reduces the acoustic velocity. Experimental measurements have been published which substantiate the theory.(3 4) The total volume of a section of the pipeline, V, can be expressed as the sum of the liquid volume, Vliqand the gas volume, Vg. A pressure change brings about a volume change which can be expressed AV = AVliq + AV (4.1..1) The gross bulk modulus of elasticity of the fluid is defined K = - nE (4.1.2) AV V and the bulk moduli of the individual components can be expressed

-47 K = _ AP liq Aliq IVliq and K = - g AVg Vg By combining these expressions the gross modulus can be expressed K Kli 1 + ( li 1) V Kg (4.1.3) By expressing the mixture density in terms of the liquid and gas densities P pg pi Vliq P = Pg V+ Pliq v (4.1.4) and by substituting into the wavespeed relation a = K/p, the following expression is obtained a= P g Kliq/PliqK a = [1- (1- g ) 0 girl+ ( -1) g'Kg Pliq v (4.1.5) As an example one percent air content in water under standard atmospheric conditions reduces the wavespeed to approximately 400 ft/sec. If the pipe wall elasticity is included Eq. (4.1.5) becomes a =_...... liq/Pli (4.1.6) [1-(1- Pg) g[l+ (Kl -1) ig + DKli vPliq V g v Ee 4.2 Wavespeed Measurements Two basically different methods can be employed to experimentally determine the wavespeed in fluid lines.

a.) If the combined bulk modulus of the fluid line can be determined -- this includes the compressibility of the fluid and the extensibility of the pipe walls -- then the wavespeed can be calculated employing the relation a2 = K/p, in which K is the combined modulus, p the fluid density and a the wavespeed. A measurement of the change of static pressure with changing volume of a pipe section should provide data for the bulk modulus. However in the experimental work outlined in Chapter III, some difficulties were experienced in completely filling the system with a representative sample, which includes the undissolved gas in the water, and also with preventing leakage from the system. Since the volume changes are small, leakage has a disasterous effect on the accurate determination of the wavespeed using this approach. b.) In most cases it seems more appropriate to measure the speed of a pressure pulse superposed on the actual steady state operating conditions of the system. If such a pulse is generated, a measurement of the time to traverse the distance between two pressure transducer stations delivers the results. The difficulty lies in the generation of a suitable pulse. In the experimental work for this project several devices were examined, such as the opening or closing of a fast acting solenoid valve on a short parallel bypass pipeline, or an automobile spark plug attached inside the pipe. However, the best results were obtained with a simple piston sealed with O-rings and fitted to a tee in the pipe. By striking the piston a pressure wave was produced which gave acceptable signals from the two pressure transducers. The lower the steady state

-49 pressure in the pipe the more pronounced is the influence of the undissolved air in the water, not only on the speed of wave propagation but also on the damping of the pressure pulse. Sharp pulses are strongly distorted as they travel through the fluid and in many cases it is difficult to determine the travel time with sufficient accuracy from the recording. This for example explains the scattering effect of the wavespeed data displayed in Fig. 3.9. 403 Effect of Cavitation Void on Wavespeed and Attenuation When the pump is cavitating,the turbulence created by the pump inducer produces severe random pressure oscillations which may be transmitted upstream into the suction lineo These pressure fluctuations, which are propagated upstream, may be visualized as pressure disturbances oscillating about the mean pressure in the system as shown in Figure 4.1l On the low pressure side of the fluctuation, vapor pressure of the liquid may be reached. When this happens a vapor pocket forms, which subsequently collapses as the pressure increases. Since the turbulence created by the inducer is quite severe, vapor pockets are being created and are collapsing in a very random manner throughout the pipe section. The extent to which these pressure oscillations penetrate upstream is dependent upon the degree of attenuation in the fluid and upon the magnitude of the turbulent condition at the pump. The lox suction pipeline can be visualized as a pipeline carrying the liquid from the tank to some point in the line above the inducer entrance, and a liquid-gas mixture (lox and vapor) in the region next to the inducer. The vapor content is highest at the pump and reduces

-50 w 1r U) LL or CL_ Figure 4.1 Pressure A:n Pressure at Point in Suction Line Upstream of the Pump Inducer versus Time. U) 0. 0 LLI n LL UJ Id 1600 1200 14.7 psia Assumed overage pressure 4.0 in vapor cavity 8001 400 O0L 0.1.2.3.4.5.6.7 VAPOR CONTENT, bg/V- (%) Figure 4.2 Variation of in Pipeline. Wave Propagation velocity in Fluid

-51 to zero at some distance from the pump as the pressure oscillations attenuate. Wavespeed The propagation velocity of a pressure wave through a liquid is reduced considerably if permanent gas bubbles are dispersed throughout the liquid. The theory for the resulting wavespeed of dispersed gas-liquid systems, Eq. (4.1.5), has been verified for the case of permanent gas bubbles for which the gas void fraction can readily be measured and the gas pressure, and hence its bulk modulus, can be inferred (34 5) Figure 4.2 displays the wavespeed variation in the suction line using some assumed values for the many variables. The wavespeed in the liquid in the pipe is assumed to be 2800 ft/sec. The effect of the presence of small quantities of gas in the range less than one percent is apparent in the figure. The results are presented for two different assumed average pressures in the gas bubbles. The presence of cavitation vapor bubbles in a liquid, even though the.bubbles are transient in nature, is expected to have a similar effect of reducing the wavespeed. No data to confirm this was found in the literature. If valid, however, it should be possible to calculate the resulting wavespeed using Eq. (4.1.5) provided the effective cavitation void fraction and gas bubble modulus were known. Both of these factors depend on the time history of the individual bubbles which in turn depends on the imposed fluctuating pressure field causing the cavitation, and these factors are not readily obtained. In the

-52 present study it was assumed that the gas bulk modulus could be based on an effective cavitation bubble pressure equal to the saturation vapor pressure corresponding to the prevailing liquid temperature. A more detailed discussion of the cavitation process is given in Appendix C. Pressure Wave Attenuation For a pure viscous liquid in systems of linear dimensions of the present concern, the bulk attenuation of a pressure wave is negligible except for rapid pressure fluctuations in the ultra-sonic frequency range. However, when the liquid contains dispersed gas bubbles - even in small void fractions - a significant increase in bulk attenuation is observed, particularly at a frequency coinciding with the resonant frequency of the bubbles. The resonant bubble frequency depends inversely on bubble size and for typical cavitation bubbles it is expected to be much greater than any frequency of importance for the overall propellant line dynamics. In this regard bulk attenuation in cavitation zones can therefore by ignored. In evaluating upstream pump cavitation extent, however, where as postulated high frequency pressure fluctuations at the pump inducer inlet propagate upstream to cause the cavitation, the attenuation is believed to be important. For significant attenuation the amplitude of high frequency pressure fluctuations will decrease exponentially with distance from the source near it, that is in the upstream cavitation zone. The theory of attenuation of periodic pressure fluctuations in a bubbly mixture has been verified(5) for the case of permanent gas bubbles. Although cavitation void represents a considerably

-53 more complex system it should be expected that the results for permanent gas bubbles in liquids also apply to cavitation bubbles in liquids provided appropriate average void fraction, gas pressure and bubble size were used. Such information, however, is not available and is difficult to measure. 4.4 Experiment Simulating Upstream Pump Cavitation A simple experiment was constructed to qualitatively demonstrate the reduction in wavespeed in the presence of cavitation. The experiment consists of a piston fitted in a plexiglas tube filled with distilled degased water. Oscillation of the piston produces cavitation in a region near the piston of an extent and void fraction determined by frequency and amplitude of the mechanical vibrator driving the piston. Figure 4.3 shows the experiment and Figs. 4.4 and 4.5 show a close-up of the region near the piston without and with cavitation. Measured resonant frequencies without cavitation were in agreement with theory, Table 4.1. Table 4.1 Experimental Wavespeeds Fundamental Har- Measured Wavespeed Run monic Frequency Frequency fr(hZ) a (ft/sec) Comment 1 (3/4) (a/L) 333 1686 no cavitation 2 (1/4) (a/L) 91 1380 no cavitation 3 (5/4) (a/L) 556 1682 no cavitation 4 (1/4) (a/L) 114 1730 no cavitation 5 (1/4) (a/L) 29.4 445 cavitation The equivalent wavespeed a was calculated from the measured resonance frequency fr by

Figure 4.3 Cavitation Apparatus Corprising L-hshaped Plexiglas Tube Fitted with a Piston Driven by the Mechanical Vibrator.

Figure 4.4 Clse.o.-up of Plexiglas Tube with Piston Cavita-tiAon No Figure 4.) 5 Close-iup of Plexiglas Tube with Piston. Visible Cavitation.

-56 a 4 Lfr n = 0, 1 2 - (2n+l) where L = 45.5 inches in the experiment. The geometry and properties of the test section are D = 1.5 inch inside pipe diameter t = 0.125 inch wall thickness E = 4,5 x 105 psi modulus of elasticity of pipe K = 3.2 x 105 psi bulk modulus of fluid (water) p = 1.94 slug/ft3 density of fluid The predicted wavespeed in the absence of cavitation is thus 1582 ft/sec. by Eq. (4.1,1). Measurements of resonant frequencies and upstream cavitation extent with cavitation were analyzed by the impedance method to calculate the corresponding average void fractions in the cavitation regiono Given the measured quarter wave frequency, the lengths of the cavitating and non-cavitating section, and the wavespeed of the fluid without vapor content, the wavespeed for the cavitation zone was found by an iteration process. Corresponding to these wavespeeds, void ratios were calculated based on the wavespeed equation derived for a liquid-vapor system, including the pipe elasticity term, Eq. (4.1.6), Due to the relatively high elasticity of the plastic pipe used in the experiment, this term has a significant effect on the wavespeed, Table 4~2, Although the void fractions appeared to be reasonable these could not be verified experimentally. It was concluded from this study that the presence of cavitation bubbles causes a reduction in wavespeed and that its order of

-57 magnitude is in reasonable agreement with theory verified for gas bubbles. Table 4.2 Quarter-wave frequency cps 31 31 42 52 63 Experimental Period sec..0323.0323.0238.0193.0159 & Computational Length of Cavitation Section feet.125.167.125.083.083 Results with Cavitation Wavespeed Void Rat in Cavitation % Section fps 138.163 155.128 190.085 200.077 260.045 tio 4.5 Model Correlating Suction Line Resonant Frequency to Suction Pressure Background Pump cavitation in high speed propellant pumps of the inducer type has been observed to occur not only in the inducer and impeller but also upstream of the inducer inlet as much as 3-4 pipe diameters into the suction line.(2) Under such conditions cavitation affects not only the steady pump performance but also the dynamics of the propellant system comprising tank, suction line, pump and discharge line. Cavitation in the suction line increases its compliance hence lowers its resonant frequency which indirectly also affects the pump gain, (see Appendix A). The cavitation extent for a given pump speed and flow rate depends on the pump suction pressure. The effect of suction pressure on the resonant frequency of the Titan II, Stage I, propellant systems with and without cavitation is shown in Figures 4 and 5 of

-58 Appendix A. Such information is necessary to predict overall system dynamic performance and must be obtained experimentally since the cavitation extent cannot be predicted. As described in Appendix A it has been the practice to measure the propellant system dynamics on full scale actual systems using a pulse technique to obtain frequency response data. Detailed information regarding the extent of the cavitation zone upstream from the pump inducer inlet for a given operating condition can be obtained by pressure measurements of standing wave patterns in the suction line or by pulsing the discharge line and measuring the pulse travel time to various positions in the suction line. The latter method gives measured values for local wavespeed from which the local cavitation void fraction can be deduced if the properties of fluid and pipe are known. To facilitate a detailed prediction of a propellant system dynamics by analysis using the method of characteristics or the impedance method, it is necessary to know pump performance, and, ideally, the suction line cavitation extent in terms of void fraction versus position. The latter information is not explicitly available from tests, however, it is implicit in the experimental data expressed in terms of system resonant frequency versus pump suction pressure for a given pressure (see Figure 4 and 5 of Appendix A). Model To analyze the available data on pump suction pressure versus suction line resonant frequency by the impedance method the following simple model is postulated(6, Appendix A) for the cavitation void simple model is postulated for the cavitation void

distribution in and up-stream of a propellant pump of the inducer typeo - Cavitation void is confined to the inducer whenever the suction pressure Ps is greater than a critical value pt for incipience of upstream cavitation - The inducer can be treated as a short section of the suction line of length Lio - For suction pressures ps less than the critical value p' additional upstream cavitation occurs to a length Lc determined by the exponential attenuation of the pressure fluctuations generated at the inducer inlet and causing the cavitationo - The average void fraction in the cavitation region is inverse proportional to the cavitation indexo This model was used together with the impedance method to correlate the data from the lox pump suction line of Titan II, Stage Io First, in the impedance method the inducer and suction line were treated as a pipe having a noncavitating length (L - L - Lc) and a cavitating length equal to Li + L, and of uniform void fraction ao The results of a parametric set of calculations are shown in Figure 4~6. Next, the model was fitted to the data yielding the dashed line in Figure 4.7~ There is conceptual agreement between the data and the results of this semi-empirical analysiso However, more experimental data are required to facilitate prediction of performance, rather than simply correlation of a single set of data. If the model should prove to be generally valid two analytical steps remain. First, to predict from pump geometry, speed and operating conditions the critical suction pressure p' for incipience of upstream cavitation, and secondly, to s

10.0 8.o- \ z0 04.o - 4.0 TA a L, t \ I 2.0 *00 0 1.0.8 --- Or) D 4.3.4.5.6.7.8.9 1.0 f/fo RATIO OF RESONANT FREQUENCIES WITH AND WITHOUT CAVITATION Figure 4.6 TITAN II, Stage I, LOX Suction Line Resonant Frequency versus Cavitating Length Lc and Average Void Fraction a.

1.0 () z Zo LL - z D ~0 Z I 0 F LLI B c 0 LL Z 0 < O I cr 0 4-.8.6 I B - - - - -.- I I I I I I I.4 H / Experiment - - - Correlation with model for upstreom and inducer cavitation.2 m / / 0 1 0 20 40 60 80 100 120 140 160 SUCTION PRESSURE, ps (psio) Figure 4.7 TITAN II, Stage I, LOX Suction Line Resonant Frequency.

-62establish the proper dependence of local cavitation void fraction on suction pressure level and other significant parameters.

V. SUMMARY AND CONCLUSIONS This study is concerned with an investigation of the dynamic behavior of the propellant feed lines of Saturn V. The original Contract objectives have been accomplished and several other relevant topics have been investigated. A Fortran programmed analysis of the propellant feed system, which was one of the primary specific objectives, is complete and has been supplied to the Contractor. A copy of the program and the accompanying explanatory write-up is included in Appendix B. The program was written with the particular configuration of Stage I, Saturn V, in mind, but is adaptable to the other stages of the vehicle by appropriate data changes. A check-out procedure for the programmed analysis involved the use of geometric data from Titan II and experimental data taken on the Titan II test stand. Favorable comparisons were made between experiment and analysis. A second separate report has been issued, and is included as Appendix A herein, in an attempt to provide a basic understanding of the entire missile instability problem. The various component parts of the system were discussed and the influence of each of the components on the stability of the vehicle is explained. A successful experimental program has been completed, in which the uncertainties regarding the analytical modeling of the system excitation due to pump oscillation were removed. Additional data were also obtained. These include system response to pump operation in the -63

cavitating region, dynamic gain across the pump, system instabilities when coupled with structural motions, and the influence of a volume compensating device in a forced vibration. A self-excited oscillation was produced in the pump oscillation apparatus and the characteristics of the instability were very similar to the well-known POGO vibration. Under certain operating conditions the cavitating pump, fluid system, and piping structure provided a complete feed-back loop with a positive gain so that an auto-oscillation occurred. An external forcing function was not required to either initiate or maintain the motion. The actual high speed propellant pumps experience cavitation in the inducer even at relatively high suction pressures and, furthermore, the inducer causes upstream cavitation to extend into the suction line at lower suction pressures. The effect of cavitation on suction line wavespeed for a given system can be inferred from available data on suction line resonant frequency and suction pressure. A simplified experiment was used to gain additional insight into vapor bubbles forming near a forced vibration and their influence on wavespeed and attenuation. Chapter IV and Appendix C discuss cavitation and related problems as well as the experimental program. The inclusion of the pressure volume compensating duct in the suction pipeline greatly reduces the susceptibility of the Saturn V vehicle to possible longitudinal oscillations. In the pump oscillation experiment the addition of the constant volume device virtually eliminated the transients in the system caused by a forced pump motion. Similarly the programmed analysis shows that the P.V.C. is very effective in

-65controlling the amplitude of pressure transients caused by physical displacement of the pump.

REFERENCES 1. Streeter, V. L. and E. B. Wylie, Hydraulic Transients, McGraw-Hill Book Co., Inc., 1967. 2. Streeter, V. L. and R. H. Fashbaugh, "Resonance in Liquid Rocket Engine Systems," Trans. ASME, Jour. of Basic Engineering, Vol 87, Dec., 1965, pp 1011-1017. 3. Kobori, T., Yokoyama, S., and Miyashiro, H., "Propagation Velocity of Pressure Wave in Pipe Line," Hitachi Hyoron, Vol 37, No. 10, Oct., 1955. 4. Pearsall, I. S., "The Velocity of Water Hammer Waves," Symposium on Surges in Pipelines, Inst. of Mech. Engineers, Proc. 1965-66, Vol 180, Part 3E. 5. Silberman, E., "Sound Velocity and Attenuation in Bubbly Mixtures Measured in Standing Wave Tubes," J. Acoust. Soc. of Amer., 29, p 925, 1957. 6. Smith, W., Atkinson, G. L., and Hammitt, F. G., "Void Fraction Measurements in a Cavitating Venturi," The University of Michigan, Industry Program Report IP-581, (Sept. 1962). -66

APPENDIX A LIQUID PROPELLANT MISSILE LONGITUDINAL OSCILLATION R. H. Fashbaugh

TABLE OF CONTENTS Page LIST OF FIGURES.................;.......................... iii I INTRCDUCTION............................................... 1 II DESCRIPTION OF THE MISSILE OSCILLATION PROBLEM............. 2 A. Description of the Missile System........... 2 B. The Relation Between the System Resonant Frequencies.. 5 C. Factors Affecting the' Oscillation Other Than Frequency. 11 D. The Dynamic Gain of the Propellant Pumps.............. 16 E. Suppressing the Oscillation.......................... 19 III TESTING TO DETERMINE THE PROPELLANT SYSTEM RESONANT FREQUENCY AND PUMP'GAIN...................... 24 IV AN ANALYSIS OF THE PROPELLANT SYSTEM USING THE WATERHAMMER CHARACTERISTIC METHOD.................................. 29 V ANALYSIS OF MISSILE IN FLIGHT...................... 31 VI SUMMARY AND CONCLUSIONS.............................. 32 REFERENCES.......... *...... 0 0 6 33 ii

LIST OF FIGURES Figure Page 1. Missile Rocket Engine System............................. 3 2. Diagram of Combined Propellant System, Engine, and Structure................................................. 4 3. Structure Resonant Frequency Versus Flight Time........... 6 4. Oxidizer System Resonant Frequency Versus Pump Suction Pressure................................................. 8 5. Fuel System Resonant Frequency Versus Pump Suction Pressure.................................................. 6. Comparison of the Propellant Systems and Structure Resonant Frequencies............................... 10 7. Characteristic Velocity of the Propellants Versus Mixture Ratio..................................................... 15 8. Typical Steady State Pump Characteristic Curves........... 17 9. Pump Suction System Frequency Response with an Accumulator at the Pump Inlet.......................................... 22 10. Oxidizer Propellant System Test Configuration............. 25 11. Oxidizer Pump Transient Suction Pressure Versus Frequency. 27 12. Oxidizer Pump Gain and Discharge System Pressure Ratio Versus Frequency............................... 28 iii

I. INTRODUCTION This report is concerned with the type of longitudinal oscillation of a liquid propellant missile that is caused by a regenerative coupling between the missile structure and the propulsion system. The purpose is not to give a complete detailed analysis of the problem, but rather to convey a basic understanding of the problem by discussing the role that each of the basic missile components plays in the development of an oscillation. In addition, testing that is necessary for the proper analysis of the problem is discussed, and an analysis approach using waterhammer methods is presented to illustrate how a stability analysis can be developed using test data to substantiate the propulsion system representation. -1

II. DESCRIPTION OF THE OSCILLATION PROBLEM A. Description of the Missile System. The basic missile configuration to be considered is shown in Figure 1 with unnecessary details omitted for simplicity. Two engines are shown with identical propellant systems, and the engines are considered mounted so that their thrust acts at common points on the structure. In an analysis the engines can be considered as acting in parallel and replaced by one engine with double the thrust. This was the case in the Titan II missile oscillation where the flight data showed the pressure oscillations at respective points in each engine system to be of the same magnitude and in phase.(1) If engines are involved that are different or that have non-identical propellant systems, or are mounted differently to the structure, they should be considered separately in an analysis. The turbine gas generator, the turbine, and the pump gear drive assembly are assumed to have no affect on the oscillation. In other words, the pump speed is assumed constant. This also was the case in the Titan II missile oscillation. () To get a clearer understanding of the cause of the oscillation, it is helpful to present the missile system in a block diagram representation as shown in Figure 2. In this figure the components that contribute separately to the development of the oscillation are shown as separate blocks. The oscillation is a structural oscillation so the main block of the diagram is considered the structure block. The remaining blocks then can be considered as acting as a feedback system from the output of the structure block (velocity or acceleration) to the input of the structure block (thrust). The structure vibration causes -2

-3 STRUCTURE THRUST CHAMBER VALVE CHAMBER Figure 1. Missile Rocket Engine System.

STAGE I OXIDIZER SYSTEM _ H _ H _ ~ — ^ —^ ~DISCHARGE V PUMP v LINE | INJECTOR po EXCITATION I CLOSED LOOP /. OPTION MECHANICAL inutTHRUST STRUCTURE -__- __ -____ CHAMBER I STAGE I FUEL SYSTEM Figure 2. Diagram of Combined Propellant System, Engine, and Structure.

-5fluid oscillations in the fuel and oxidizer systems. These fluid oscillations cause an oscillation in thrust chamber pressure which causes an oscillation in thrust. The thrust oscillation is the input to the structure which completes the loop. If this resulting thrust oscillation is of proper magnitude and phase (relative to the structure motion) to reinforce the structure oscillation, a sustaining oscillation will occur. Thus, the problem is of the common stability type. The stability analysis of the complete system is quite complex, however, due to the complicated nature of the structure and to the interaction of the component parts of the propellant systems. In Figure 2, H refers to fluid pressure, V to fluid velocity, T to thrust and x to structure velocity. There are basically three separate systems involved, each with a natural resonant frequency: the structure, the fuel fluid system, and the oxidizer fluid system. The relation between these three resonant frequencies is of prime importance in the problem and will be discussed in the following section. B. The Relation Between the System Resonant Frequencies. The resonant frequencies of the structure, the oxidizer system and the fuel system all vary with the missile flight time. The structure resonant frequency increases with flight time since the propellant mass decreases. A typical variation is shown in Figure 3. Due to the complex nature of the structure, several modes of oscillation, and therefore several resonant frequencies, are possible. In the figure the two lowest frequency modes are shown.

30 20 2N MODE U )I z sL. I0 0 100 T MOmE I I) 110 120 130 I.f% Figure 3. FLIGHT TIME " SEC Structure Resonant Frequency Versus Flight Time. 150

-7 The oxidizer and fuel system resonant frequencies are shown in Figures 4 and 5, respectively. In this case we are concerned with the resonance of the fluid in the tank, pump suction line, pump, and pump discharge line. The resonant frequency is not affected by the level of propellant in the tank since the hydraulic impedance at the tank-suction line junction is practically zero. This has been proven experimentally. ( One would expect the resonant frequency to be constant, however, as Figures 4 and 5 show, this is not the case. The frequency varies with pump suction pressure because of vapor at the pump inlet caused by cavitation of the pump impeller. This cavitation vapor also lowers the frequency a considerable amount from the value one would expect with no vapor in the system. This difference is shown for a typical case in the figures. The pump suction pressure varies with the missile flight time, and therefore, the fuel and oxidizer system resonant frequencies vary with the flight time. These resonant frequency curves have to be obtained experimentally since the extent of the vapor formation by cavitation cannot be predicted. The curves in Figures 3, 4, and 5 are plotted versus the flight time in Figure 6 to show the relationship between the resonant frequencies. This particular plot is typical of the Titan II missile (Stage I)( The missile oscillates at the frequency of the structure, and the oscillation begins when the oxidizer system resonant frequency becomes close to the structure frequency. The oscillation continues over to some point beyond the intersection of the fuel system frequency curve and the structure frequency curve and then stops.

16 - WITH NO PUMP CAVITATION 14 H 12 cn CL z IUJ 0O u.I 10 F Go 8 - WITH PUMP CAVITATION 6 - 4L 30 1 I 50 70 90 OXIDIZER PUMP SUCTION PRESSURE- PSIA Figure 4. Oxidizer System Resonant Frequency Versus Pump Suction Pressure. 120

55 50 20 WITH NO PUMP CAVITATION K c) I ao z w O 15F WITH PUMP CAVITATION 10 5 U'. - I 15 20 25 30 FUEL PUMP SUCTION PRESSURE - PSIA Figure 5. Fuel System Resonant Frequency Versus Pump Suction Pressure. 35

50s H-UEL SYSTEM WITHOUT PUMP CAVITATION 16h SYSTEM WITHOUT PUMP CAVITATION 1) a. U C) I z a L I-. z 4 z 0 U) FUEL SYSTEM 14h 12 STRUCTURE I 0 I IOF OXIDIZER SYSTEM 8 I rOSCILLATION STARTS I I I I I OSCILLATION V STOPS I 6 80 100 120 FLIGHT TIME - SEC. 140 160 Figure 6. Comparison of the Propellant Systems and Resonant Frequencies Structure

-11 The predicted fuel and oxidizer system resonant frequencies without the effect of pump cavitation are shown on Figure 6 to illustrate the importance of accurately knowing the effect of cavitation on the frequency. As seen, the pump cavitation completely changes the problem. Also, any gas that is trapped in the pump suction system (such as in a PVC) can have a similar effect. C. Factors Affecting the Oscillation Other Than Frequency. The possibility of the development of an oscillation in a system depends upon many factors, in addition to the natural frequency of the component parts. The main factors are: the structure gain, the pressure losses in the fuel and oxidizer systems, the dynamic gain of the fuel and oxidizer pumps, the pressure oscillation generated at the bottom of the fuel and oxidizer tanks due to tank motion, the pressure oscillation generated at the fuel and oxidizer pump inlets due to pump motion, the thrust oscillation per unit amplitude of thrust chamber pressure oscillation, and the combustion time lag. The structure gain at a particular point on the structure is the oscillation amplitude at that point per unit amplitude of thrust oscillation. It depends upon the mode shape of the oscillation, and for a given mode of oscillation varies from point to point along the missile. A vibration analysis of the missile structure will yield the structure gain as well as the structure frequencies. In the case of the Titan II missile, the structure gains and frequencies agreed adequately with values measured in flight (l) The fuel and oxidizer system pressure losses are important in the overall gain of the propulsion system. The pump discharge lines

-12 typically have a much higher pressure loss than the pump suction lines. As a consequence, the resonant frequency of the system is essentially determined by the fluid in the pump suction line (including the cavitation vapor), and the pump discharge line (including injector) provides the predominant damping which limits the amplitude of the resonant curve. It should be mentioned that in representing the long suction lines analytically, a distributed representation should be used, the most accurate being a waterhammer type analysis. When the resonant frequencies are known the fuel and oxidizer systems can be represented very well by a waterhammer analysis. The measured steady state pressure losses can be used to determine the pump discharge system friction factors.(2 3) The dynamic gain of the propellant pumps is the oscillation amplitude of the pump discharge pressure per unit oscillation amplitude of the pump suction pressure. This ratio depends upon the steady state characteristic curves of the pump, and can be less than or greater than one. In the case of the Titan II missile, the oxidizer pump gain was about 0.8, whereas the fuel pump gain was in the neighborhood of 4.0.(1) This will be discussed in more detail in section D. In setting up an analysis of the oscillation, it is difficult to properly represent the effect of the motion of the structure on the fluid systems. The main points where the fluid system is excited are at the tank bottom, the pump, and the thrust chamber injector. It is important that the phase angle between the generated pressure oscillations and the structure oscillation be correct as well as the relation between the oscillation amplitudes, since the system stability is sensitive to

-13 this phase angle. One method used is simply to apply the acceleration of the tank directly to the fluid in the tank and the acceleration of the pump directly to the fluid in the suction line, with the fluid system represented by a lumped compliance and lumped inertance. It is felt that the tank excitation is adequately represented in this manner, but that the excitation at the pump is questionable because the change in pump pressure rise due to the change in the fluid velocity relative to the pump caused by the motion of the structure is neglected. Also, because the flow can change through the pump as the pump moves, the fluid in the suction line does not experience the same acceleration as the pump structure. Another approach is to represent the system with the waterhammer analysis, (2,3) and excite the system by varying the pump pressure rise with flow relative to the moving pump, and move the pump boundary with the motion. The waterhammer equations will yield the proper transients developed. A test is currently being conducted at the University of Michigan to establish the latter method of representing the pump motion. The thrust chamber injector motion causes an oscillation of the pressure drop across the injector. This effect is probably small but can easily be introduced into a waterhammer analysis of the system. The thrust oscillation depends directly upon the thrust (4) chamber pressure oscillation, and is given by, T = n At Cf Pc where n is the number of engines considered in parallel, At is the nozzle throat area, Cf is the thrust coefficient, and Pc is the thrust chamber pressure. The thrust coefficient varies slightly with

the thrust chamber pressure and the mixture ratio (flow rate of oxidizer over flow rate of fuel), but very little error is made by treating it as constant. The thrust chamber pressure is given by, () C* (w + f) Pc = At g where, C is the characteristic velocity which is a function of the mixture ratio, wo and wf are the oxidizer and fuel flow rates respectively, and g is the gravitational constant. To illustrate the variation of the characteristic velocity with mixture ratio, a typical curve is given in Figure 7 with a typical operating point shown. Since the slope at the operating point can be quite steep, small variations in the mixture ratio can cause substantial variations in the characteristic velocity. Therefore, this effect should be included in an analysis. The combustion time lag is small compared to a period of oscillation at the frequencies concerned with herein (11 cps in the Titan II case), and therefore causes only a small phase shift. This time lag is due to the time required for the combustion gas to accelerate from zero velocity at the inlet end of the combustion chamber to sonic (4) velocity at the nozzle throat. It can be estimated from, p C L c c PC c where, pc is the average mass density of the gas in the combustion chamber, and Lc is the average combustion chamber length. To include this time lag in a waterhammer digital computer type analysis, it is only necessary to store the Pc values for the time Tc and then compute the thrust. In an analogue analysis the Laplace operator S can be

| ~/1 -^5~~ ok TYPICAL OPERATING POINT II) L 0m I aW E,)or o M- MIXTURE RATIO - Figure 7. Characteristic Velocity of the Propellants Versus Mixture Ratio.

-16 used and the chamber pressure computed from, 1 C (w + w) PC = ------- [ - o f f TC S + 1 At g D. The Dynamic Gain of the Propellant Pumps. The propellant pump pressure rise is a function of the pump suction pressure, the flow rate and the pump speed. A typical set of steady state pump characteristic curves is shown in Figure 8. This set of curves will be used to obtain the dynamic pump pressure rise. Using the steady state curves to obtain the dynamic or transient pump pressures has been varified. (2) The pump speed will be considered constant since the mass of the pump rotor, gear train, and turbine combined is large. (Titan II experience varifies this).(l) The pump discharge pressure, Pd, is given by, Pd- p + PS where P is the pump pressure rise and Ps is the suction pressure. Since the pump pressure rise is a function of the suction pressure and flow rate, so is the discharge pressure, and we can write the change in the discharge pressure as, aPd da d dPd = - dP + where w is the pump flow rate. Also, from above, pd + 1 aPs aPs = m + 1

1200 0. nr e) w Cn CfJ o~: wJ 1000 8: 0. a800 800 LO PE,- R I i-, I 20 60 100 460 500 540 580 (a) PUMP SUCTION PRESSURE, psia (b) WEIGHT FLOW, Ib/sec Figure 8. Typical Steady State Pump Characteristic Curves.

-18 where m is the slope of the pressure rise versus suction pressure curve, Figure 8a. (Note: NPSH is usually the parameter used, but suction pressure is used here since the fluid temperature and thus the vapor pressure is constant). Similarly, d = P = - R where - R is the slope of the pressure rise versus flow rate curve, Figure 8b. Since the slope is negative, the minus sign is put in so that R is a positive number. Combining we obtain, dPd = (m + 1) dPs - Rdw and the pump dynamic gain is therefore given by, ~dP ~d ddw - = (m + 1) - R dPs dPs The first term on the right side represents the effect of the pump pressure rise - suction pressure characteristic, and the second term represents a gain loss due to the pressure - rise - flow rate characteristic. R can be thought of as the pump resistance. The second term also reflects the nature of the fluid impedance of the pump piping system. This shows that the pump gain depends upon the dynamic characteristics of the pump suction and discharge systems as well as the pump characteristics, a fact that should be considered when measuring pump gain in a test. The missile oscillation is sinusoidal so we can represent the pressures and flows as vectors and, Pd = (m + 1) Ps - R Ws

-19 which is a linearlized form where the slopes m and - R are considered constant for a small oscillation around the operating point. The above is useful in gaining an understanding of the pump dynamics and in developing an analogue analysis of the pump. If the waterhammer characteristic method(3) is used the steady state pump curves are used directly in the program. This yields the proper pump gain with the non-linear aspects of the curves present. (2) E. Suppressing the Missile Oscillation. It is of interest to consider some methods of suppressing the missile oscillation. The methods to be briefly discussed are: varying the propellant tank gas pressure, changing the pump suction line, and the addition of an accumulator (surge chamber) to the suction line at the inlet to the pump. Changing the gas pressure in the top of the propellant tanks has two effects on the system. One effect is to change the operating point of the pump on the pressure rise - inlet pressure curve, Figure 8a, and the second is to change the resonant frequency of the propellant system by changing the amount of cavitation and the effective compliance of the cavitation vapor, Figure 4. Changing the operating point changes the pump gain by changing the slope, m. The gain is reduced by increasing the pressure. If the pump is normally operating on the zero slope part of the curve, the pump gain cannot, of course, be reduced by increasing the pressure. In the Titan II case the fuel pump was operating well to the left of the zero slope region, and the oxidizer pump was operating on the zero slope part, so only the fuel tank pressure was changed.

-20 An alteration of the propellant system resonant frequency by changing gas pressure has two effects. This frequency can be moved farther from the structure frequency which reduces the gain of the propellant system, and also changes the phase angle between the structure motion and the thrust oscillation output of the propellant system. The phase angle change might or might not help, depending on whether the phase angle between input thrust and output thrust is increased or decreased. With the frequencies shown in Figure 6, increasing the fuel tank pressure will move the fuel resonance away from the structure resonance and reduce the fuel pump gain. A decrease of the oxidizer tank pressure will move the oxidizer resonance away from the structure resonance, but will not appreciably change the oxidizer pump gain. In the case of the long pump suction lines changing the line diameter, wall thickness or material can change the resonant frequency of the system. The resonant frequency is decreased if the pressure wave velocity in the suction line is decreased, and the pressure wave velocity is decreased if the diameter of the line is increased, the wall thickness decreased, and/or the line material modulus of elasticity decreased. In the example of Figure 6, the oxidizer frequency can be moved away from the structure frequency for instance by changing from a stainless steel line to an aluminum line. The above methods are limited by missile design restrictions and might only succeed in reducing the amplitude of the missile oscillation. A more effective way of suppressing the oscillation is the addition of accumulators at the pump inlets. By properly designing the accumulators, the fundamental resonance of both propellant systems can be

-21 moved considerably below the structure resonance which will stabilize the missile. However, the accumulators introduce two additional resonances to each propellant system, an anti-resonance (zero) which is of no concern (except in phase considerations), and a higher resonance (pole) which needs to be above any prominent structure resonances. Fortunately, at the frequency of the high resonance, the propellant system gain is much lower than at the low resonance. The suction system response curve is shown in Figure 9 with and without an accumulator. The curve is a log - log plot of the ratio of the amplitude of suction pressure oscillation over the amplitude of the fluid acceleration oscillation versus frequency. The accumulators can be either of a piston - spring type design or a compressed gas type (standpipe). The latter is preferred because the mass (inertance) in the accumulator can be low which allows the high resonance to be higher than with the piston type. Also, since a piston has a static breaking force, the piston type accumulator is not effective until an oscillation starts, or unless the noise present (or random vibration) causes large enough pressure fluctuations to overcome the static friction. In designing an accumulator it is advantageous to get the low resonance as far below the structure resonance as possible without having the high resonance too low. The pump cavitation vapor can extend three or four feet into the suction line, and the accumulator will most likely be right in this region since it should be as close to pump as possible to get the low resonance low and the high resonance high. This vapor causes difficulty in analyzing the system with an accumulator. A waterhammer analysis is probably the best since it can approximate the effect of the vapor distribution.

I 0 log (FREQUENCY) Figure 9. Pump Suction System Frequency Response with an Accumulator at the Pump Inlet.

-23As might be expected, putting an accumulator on one propellant system and not the other might cause the missile oscillation to become worse. It is possible for the propellant systems to be out of phase such that one is cancelling the effects of the other, and the removal of one by an accumulator can increase the overall propulsion system gain causing an increased amplitude of oscillation.

III. TESTING TO DETERMINE THE PROPELLANT SYSTEM RESONANT FREQUENCY AND PUMP GAIN. As already pointed out, the oxidizer and fuel system resonant frequencies can only be determined by testing because of the effect of pump cavitation. In fact, each particular pump has to be tested since a correlation between pumps for this factor does not presently exist. It is possible that a correlation might exist through the pump cavitation index parameter (fluid velocity head over equivalent impeller tip speed head), but there is not enough data available on missile pumps to establish any correlation. A typical test configuration is shown in Figure 10 for an oxidizer system. A pulser is used to introduce sinusoidal flow oscillations although other methods could be used such as oscillating a valve or vibrating the pump. It is important that the pump suction and discharge systems be the same from a fluid dynamic standpoint as the actual missile system. The length, the pressure wave velocity, the inertance, and friction losses in each line should be the same as in the missile system. The tank can be other than a missile tank since it does not affect the frequency. The test system should be isolated from the facility return lines by a cavitating venturiL or cavitating orifice placed at the engine nozzle location in the missile system. The thrust chamber injector can be simulated by an orifice. Proper simulation of the missile system is important in order to obtain the correct pump gain as well as the resonant frequency. The resonant frequency can be obtained by pulsing over a frequency range to obtain the frequency at maximum pump suction pressure -24

-25 OXIDIZER TANK PIPE I FEED LINE PIPE 3 VALVE DISCHARGE -- LINE-PIPE 7 VALVE -- OUTFLOW PULSER -PIPE 5'UMP-PIPE 6 CAVITATING (NOZZLE) VENTURI Figure 10. Oxidizer Propellant System Test Configuration.

-26 oscillations. A typical plot is shown in Figure 11 (Titan II stage I (2) oxidizer system). A typical plot of the pump gain is shown in Figure 12. Figure 12 also shows a curve of the pump discharge pressure oscillation amplitude over the orifice outlet pressure oscillation amplitude versus frequency which is helpful in substantiating an analysis of the system. The pulsing should be done at more than one amplitude to seek any dependence of frequency on amplitude that might be present. Other information can be obtained from the test such as the extent of penetration of the cavitation region into the pump suction line. Measuring the standing pressure wave pattern along the suction line and/or pulsing the discharge line with a square pulse and measuring the pulse time delay up the system can accomplish this. It was found that the cavitation effect extended approximately three to four feet into (1,2) the oxidizer pump suction line in the Titan II stage I system.

30 --- 200 c) /,,'^ A'^=- - gw CD W 20 150 w /0 ^~~~~~~~~~~~~ 0 WI I__ ________ _____________l 0 _____ z 0 0 00 23 FREQUENCY CPSz 0 _____ 50,)w -~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - 2 4. 6 8 10 20 30 FREQUENCY CPS I rl I Figure 11. Oxidizer Pump Transient Suction Pressure Versus Frequency.

-28 1.0 I-% w 0 a. 2 Q< w - iL, 0 0 0 <c 0.8 0.6 0.4 0.2' —- ANALYSIS -0 — TEST DATA _. = ^- =s=0 0 1.2 w I10 a. a. 0~'1 crU w z a. a. <U 0.8 0.4 o 0% 8 10 12 14 16 18 20 22 FREQUENCY CPS Figure 12. Oxidizer Pump Gain and Discharge System Pressure Ratio Versus Frequency.

IV. AN ANALYSIS OF THE PROPELLANT SYSTEM USING THE WATERHAMMER - CHARACTERISTIC METHOD. In order to establish a method of analyzing the fluid transients in a missile propellant system a waterhammer analysis of the Titan II stage I test was made, and the results of the analysis were compared with the test data. A detailed description of the analysis will not be given here since it is presented in Reference 2, and it will be assumed here that the reader is familiar with the material in References 2 and 3. The main difficulties in analyzing the propellant system are representing the flow transients through the pump and the pump cavitation. The other parts of the system are represented as explained in Reference 3. In this analysis the pump is considered as a pipe with a length equal to the mean flow path length through the pump. The same fluid characteristic equations are used as are used in a constant area pipe except the friction term is replaced by a pump pressure rise term with a change in sign. This pressure rise term can be derived in the same manner as the friction term is derived in Reference 3 or as was done in Reference 2. For reference, the pump flow characteristic equations are, Vp = 1/2 [(VR + VS) + g (HR - Hs)/ac + 2g (Lp) At] Hp = 1/2 [(HR + HS)+ aC (R - V)/g] where the notation follows that used in Reference 2. The comparison of the test data and computed data are presented in Figures 11 and 12 for -29

-30 the oxidizer system. These data show that for a stationary pump the analysis predicts the pump gain quite well. The AHp used in the above equations is the pump pressure rise predicted from the steady state curves of Figure 8. The pump cavitation is accounted for in the analysis by selecting a short length of the suction line at the pump inlet end in which there is assumed cavitation vapor. The pressure wave velocity is reduced in this section by the vapor. The value of this pressure wave velocity is computed such that the system resonant frequency will agree with the measured value of Figure 4. This is explained in Reference 2. The analysis is not sensitive to the length selected for the cavitation region since the resonant frequency will be computed the same regardless of the length, and the propellant system damping which determines the shape of the resonance curve is primarily in the discharge line. It should be mentioned, however, that if higher order resonances of the suction system are important the cavitation region length could have an effect. In the Titan II analysis a length was used that was determined approximately from the test.

V. ANALYSIS OF MISSILE IN FLIGHT The propellant system analysis presented in Section IV, which deals with the static test stand, forms the basis of a complete missile stability study. The pulser is removed and fluid system excitation caused by the pump motion, tank motion, and injector motion is added. The pump movement affects the relative velocity to the pump which affects the pump head rise, AHp, in the above equation. Flow continuity at the pump boundaries is also affected by the motion. The tank motion is inserted by maintaining flow continuity in a region containing the tank bottom to suction line junction with the tank bottom moving. The effects of the injector motion can be inserted by writing the injector pressure drop in terms of the flow velocity relative to the moving injector. The cavitating venturi in the test analysis is replaced by the engine thrust chamber and nozzle equations given in Section I to convert to the missile system. The structure is represented by a linear 2nd order differential equation with constant coefficients, with the constants determined from the structure vibration analysis, such as, x + ax + bx = cT where x is displacement of a point on the structure such as the pump, a, b, c are constants from the vibration analysis, and T is the thrust inputed to the structure. -31

VI. SUMMARY AND CONCLUSIONS. The problem of longitudinal oscillation in a liquid propellant missile is reviewed in this study. Some of the causes are outlined, remedies are suggested, and a complete analysis is presented. The role played by each of the important component parts of the system on the development of an oscillation is outlined in an effort to convey a broader understanding of the entire problem. A comparison is given between analytical and test stand data. The favorable comparison of these results lead to the conclusion that an accurate stability analysis of the entire missile system can be made. Many of the uncertainties associated with the analysis have been resolved through the experience gained on the Titan II analysis. The use of the steady state pump characteristic curves in the dynamic situation is shown to be satisfactory. A structure vibration analysis can be adequately made (Titan II, Reference 1). The tests underway at the University of Michigan will lead to an improved understanding of the boundary condition at an oscillating pump. With these results in hand, the response of the entire system should be predicted reasonably well. Support for the preparation of this report was obtained from National Aeronautics and Space Administration Contract No NAS 8-20312o -32

REFERENCES 1. Martin Company,'ogo Symposium," Denver, Colorado, February 4, 1964. 2. Fashbaugh, R. H. and Streeter, V. L., "Resonance in Liquid Rocket Engine Systems", Journal of Basic Engineering, (Dec. 1965), pp. 1011. 3. Streeter, V. L. and Lai, Chintu, "Waterhammer Analysis Including Friction", Trans. ASCE, Vol. 128, Part I, (1963), pp. 1491-1552. 4. Ring, Elliot, Editor, Rocket Propellant and Pressurization Systems, Prentice-Hall Space Technology Series, 1964. -33

APPENDIX B FORTRAN IV PROGRAM FOR ANALYSIS OF SATURN V PROPELLANT SYSTEMS E. B. Wylie

This report was prepared by The University of Michigan College of Engineering Department of Civil Engineering under Contract No NAS 8-20312 PROPELLANT LINE DYNAMICS for the George C. Marshall Space Flight Center of the National Aeronautics and Space Administration. The work was administered under the technical direction of the Propulsion and Vehicle Engineering Laboratory, George C. Marshall Space Flight Center with John R. Admire acting as project manager. ii

TABLE OF CONTENTS Page LIST OF FIGURES................................................ iv I INTRODUCTION.............................................1 II DESCRIPTION OF MISSILE SYSTEM............................ 2 III BASIC THEORY.............................................4 IV PROGRAM ORGANIZATION AND DETAILS......................... 11 A. Program Options...................................... 11 B. Program Organization................................. 13 C. Initial Steady State Conditions...................... 14 V BOUNDARY CONDITIONS...................................... 17 A. Series Pipe Connection............................... 17 B. System Excitation at the Feed Tanks.................. 19 C. Turbopump Representation............................. 19 D. System Excitation at the Pump........................ 23 E. Pulser............................................... 24 F. Pressure Volume Compensator.......................... 26 G. Feed Tank Surface.................................... 26 H. Engine for Combined System Analysis.................. 26 I. Engine When Systems are Isolated at the Combustion Chamber..............................................31 J. System Excitation at the Burner Plate Orifice........ 32 VI THE COMPUTER PROGRAM................ o............. 33 VII SAMPLE PROBLEM......................................... 41 VIII SUMMARY................................................. 45 APPENDIX A................................................... 46 REFERENCES.................................................... 52 iii

LIST OF FIGURES Figure Page 1. Schematic Diagram of Saturn'V.....oo.......... 3 2, Control Volumeo........................................ 5 3~ Characteristics in xt Planeo..o...ooooo,........o... 5 4, Characteristics at the Boundary Conditionso, oo....,.. o 5 5. Series Connectiono.oo. ooo......o.o..o...... o.o 18 60 Connection at Feed Tank.........o.................. 18 7o Pump Characteristic Curveo Pressure Rise vso Suction Pressure......................o o o. o o o. 22 80 Pump Characteristic Curve. Pressure Rise vs. Discharge........ OO..o o o................ 22 9o Location of Pulser. o.. o...............o....o.... 26 10o Schematic Diagram of Engine Combustion Chamber, Burner Plate Orifices, and Pump Discharge Lines.......o........ 26 11. Listing of Saturn V Program..o.................. 35 12. Sketch of System Showing Location of Program Output Quantities P,Q........................... o o......... 40 13. Input Data for Example Problem.......................... 42 14. Output from Analysis with Oscillating Pump........... 43 15. Pressure at Pump Inlet of Oxidizer System of Titan II Due to Pump Oscillation, c = 628 rad/se c o.............. 44 iv

I. INTRODUCTION The digital computer program to be used to analyze the propellant system dynamics of Saturn V is presented in this report. The basic theory underlying the method'of analysis is developed and the program details are described. The objective of the program is to provide an analysis of the transient flow condition that may develop in the propellant feed system of a liquid fueled rocket when subjected to longitudinal structural vibrations. The details described herein refer to Stage I, Saturn V, but the program is equally applicable to the other stages of the missile, or to other liquid-propelled missiles of similar configuration. In this report the missile propellant system is first described, the basic theory behind the method of analysis is set forth, then the program is described and details are presented. The report is concluded with the results from an example problem. -1

II. DESCRIPTION OF MISSILE SYSTEM A schematic view of the propellant feed system to one engine of Stage I, Saturn V, is shown in Figure lo Only the component parts considered to be an important influence on the fluid system response are showno The turbine gas generator and turbine are assumed to have no affect on an oscillation in the feed systemo Stage I has five engines, each with separate supply systems leading from the lox and fuel tanks. These are assumed to act independently in parallel without any fluid flow or pressure interaction at the feed tanks. Any differences that exist between the outboard and inboard systems can be treated by using the data appropriate to the system being considered. Figure 1 shows the component parts that are believed to be important in the dynamic response study, namely, the feed tank, the suction lines, the PVC, the pump, the parallel discharge lines, and the engine which includes the injector plate and combustion chamber. The fuel and lox systems are separate units with the combustion chamber beyond the orifice injector plate being the only point in the fluid flow and pressure field that is common to both systemso With this configuration a meaningful analysis can be performed in either system alone, assuming the other steady, or in the two systems together wherein an interaction is permitted at the combustion chamber. If a singlefueled-engine missile were being analyzed only one system would be considered. -2

-3 TA TANK LOX FUEL I TANK ~1 LOX SYSTEM K= SYSTM = 1 FUEL SYSTEM K = 2 SYSTM = 2 COMBINED SYSTEMS SYSTM = 3 I Pulser or PVC Lox Pump Fuel Pump 18 Lox Dome I 1 Fuel Monifold Burner Plate Combustion Chomber Figure 1. Schematic Diagram of Saturn V

III. BASIC THEORY The equation of motion and continuity equation that describe one-dimensional unsteady flow of a compressible fluid in an elastic pipeline are first developed. The equations involve the actual distributed parameters that describe the system and include the nonlinear viscous terms. These partial differential equations are transformed into particular total differential equations by use of the method of characteristics, A finite difference method is then used to place the total differential equations in a form suitable for numerical solution on a digital computer. The transient conditions in each pipeline in the system are treated in accordance with the above description. Boundary conditions are used to transfer the effect of one pipeline to another, and to introduce the effect of the pump, the PVC, the engine, etc. In a typical computer analysis of a specific problem, an initial steady-state operating condition is imposed upon the system; then a numerical solution is obtained at finite points throughout the entire system at finite time increments during the life of the transient. The equation of motion is developed by examining a control volume that is moving with the vehicle, Figure 2o dV Px A Ax + To iD Ax - pg Ax A sin a + pA Ax =0 dt The density of the fluid is given by p and the pressure is given by p. The angle a is measured between a plane normal to the flight direction -4

FLIGHT DIRECTION ArEz-AREA A ro rD Ax-, g pA Ax sin a to to p at 7A SB R S.......... I Figure 2. Control Volume Figure 3. Characteristics on the XT Plane t 4 p 7 S Figure B A A x R 4. Characteristics at the Boundary Conditions

-6 and the positive x axis, and g is the acceleration at the particular time in flight. Other variables are defined in Figure 2. The equation (1) can be simplified to the following form: PX fV2 - + -- + Wx + Vt - g sin a = 0 p 2D The subscripts x and t indicate partial differentiation. The condition of continuity applied to the same control volume at any instant yields: (pAV)x Ax + (pA x)t = 0 which can be simplified to + 1 dA + 1 dp _ 0 Vx+ + = A dt p dt The introduction of the elasticity of the pipe walls and compressibility of the fluid allows the equation to be rewritten in the form: 2 VPx Pt a Vx + + = 0 x P The pressure pulse wave speed in the pipeline is denoted by a, and, for liquid flow in an elastic tube, is defined by a -- KKD 1+ KD Ee The bulk modulus of elasticity of the liquid is given by K, the modulus of elasticity of the pipe wall by E, and the pipe wall thickness by e. The equation of motion and continuity equation can be combined

-7 in a linear manner using an unknown multiplier, X. Thus, by multiplying the continuity equation by X, adding the simplified equation of motion, and rearranging, we have _ I i i t +fV Px (V + ) + Pt + [Vx (V + a 2) + Vt + - g sin a = 0 If dx is restricted to the values V + 1 and V + a2k, then this dt X equation can be written as a total differential equation. X dp dV fV -- + — +- - g sin a = O p dt dt 2D The restriction imposed upon dx requires that X = +. When these dt a values of X are substituted, the following set of equations results 2 1 dp + dV + fV g sin a = (1) ap dt dt 2D C+ =V+ a (2) dt 1 dp dV fV2 - d + dV - g sin a = (3) ap dt dt 2D Cd V-a (4) These are the characteristic equations. The first equation is valid only if the second is satisfied, and the third equation is valid only if the last is satisfied. By use of the first order finite difference approximation to these equations, they can be placed in a form which is suitable for numerical solution on the digital computer. The independent variable (xt) plane is shown in Figure 3. If conditions (p,V) are considered known at points R and S, then they

-8 can be determined at point P at time At later. In Eqs. (2) and (4), the velocity V is generally much less than the wave speed a, and can be neglected. Since a is a constant in any given pipeline, Eqs. (2) and (4) are straight lines on the xt-plane, RP and SP, respectively. It is also convenient to write the equations in terms of mass flow rate (Q = pAV, slugs/sec.). In finite difference form, Eqs. (1) to (4) become QP - QR + a (PP - PR) + f At QR - gpA At sin a = o (5) 20DpA xp - xR = a (tp - tR) (6) QP - Q A - ap PS)a A - gpA At sin a = 0 (7) xp - xS = -a (tp - tS) (8) These four equations enable a solution to be obtained for the four unknowns, p, pp, Xp, and tp, if conditions are known at R and S. An orderly computer solution is obtained if each pipe length is divided into N equal reaches, Ax = L/N. This specifies the time increment that can be used since for stability and convergence of the solution it is necessary that At a < Ax in each pipe. A specified time interval grid is therefore established in each pipeline, Figure 3. If conditions are known along the system at the initial time, to, they can be determined at time At later by a systematic application of Eqs. (5) and (7) at each point P. Simultaneous solution of Eqs. (5) and (7) yields P = QR + QS + A( (P QR-2 + S )] + gpA At sin a (9) 2- L a 2DpA

-9 PP PR + PS + AAt2Dp (QR 2 ) Q 2) (10) PP:7-A R S A2Dp ( (lo Conditions can be found at points R and S by use of a linear interpolation ) between known conditions at points on the specified time interval grid A, C, and B. Equations (9) and (10) permit the calculation of all interior points of the grid for all time, however, the influence of the pipe end conditions must be introduced to yield the complete solution. At the inflow end (left side), Fig. 4, of the pipeline, Eq. (7) can be written in the form A fP t 2A = Q a S - 2A + gpA At sin a + A (11) $=&s~.a PPDAp or p = C1 + C2 pp (12) where C1 is a variable that can be evaluated in terms of known quantities at each new time step, and C2 = A/a. Similarly, at the outflow end (right end) of the pipeline, Eq. (5) can be written ApR f At 2 A QP = QR + fAt QR + gpA At sin a - pp (13) a 2DAp a or p = C3 - C2 pp (14) where C3 can be evaluated in terms of known values. When an external condition that relates Qp and pp, or defines one of the variables, is applied to a pipe end a solution can be obtained in combination with either Eq. (12) or (14) for the inflow or outflow end, respectively.

-10By combining the entire system in an orderly manner with the appropriate boundary conditions, a complete solution is obtained for points throughout the entire system for all time during the life of the transient condition.

IV. PROGRAM ORGANIZATION AND DETAILS Figure 1 displays a schematic diagram of the entire propellant system from feed tanks to combustion chamber. The component parts of the systems have been numbered to enable simple identification in the program. The results from any given execution of the program show the complete system response to a particular excitation. That is, an initial steady state operating condition is defined for the system; then, the transient condition produced by the exciter is followed for as long as may be desired. A summarized discussion of the options in the program are presented below, followed by an itemized account of the major steps in the program. The details involved in the establishment of an initial steady operating condition are then presented. In general, the order of discussion of the program details follows the order of computation in the actual program. A. Program Options The following major options, pertaining mainly to the system to be analyzed and the type of exciter to be used, must be decided upon by the program user and specified in the data. (a) System to be analyzed. Both lox and fuel systems can be analyzed in the same run, in which case the variable SYSTM = 3. Excitation can be in either or both -11

-12 of the systems. If the lox system only is to be analyzed, the variable SYSTM = 1; and if the fuel system only is to be analyzed, the variable SYSTM = 2. (b) Combustion chamber. If one system is excited, it may be desirable to isolate that system from the other in order to study its response alone. This is accomplished by setting the variable SINGLE = 1, which isolates the systems at the combustion chamber and essentially assumes steady-state flow from the second system. The more common situation, wherein one system may influence the other through the combustion chamber, is obtained by the variable SINGLE = O. By using this option, one system can be excited, and the extent to which it influences the other system can be studied. Or, both systems can be excited and the resulting oscillations in the combustion chamber can be observed. (c) Pressure volume compensator. The PVC is included in the fuel and lox systems if the variable PVC = 1. If PVC = 0, the PVC is not in the system. Thus the effectiveness of the PVC for any given excitation can be evaluated. (d) Points of excitation. There are three possible excitation points in each system and one disturbance mechanism which is common to both fuel and lox systems. The desired amplitude and phase angle can be given for each of these elements at the frequency OMEGA. The table summarizes the variables involved. The variables pertaining to a particular location must be set equal to one to use the exciter and an amplitude must be given. If the exciter is wanted in only one of the systems, the variable controlling

-13 the amplitude in the second system must be set to zero. In order to remove any exciter from the system the variable controlling that element is set to zero. Lox Lox Fuel Fuel Location Variable Amplitude Phase Amplitude Phase Junction of supply tank to suction pipeline TANK = 1 AMl(l) ALPHA(l) AM1(2) ALPHA(2) Pulser in suction line PULSER = 1 AM5(l) ALPHA5 AM5(2) ALPHA5 Pump PUMP = 1 AM7 (1) ALPHA7 AM7(2) ALPHA7 Burner Plate Orifice ORIFIC = 1 AM9 ALPHA9 AM9 ALPHA9 B. Program Organization The computing procedure followed in solving a particular problem involving a particular set of data is as follows: a) Read in the data pertaining to the system. b) Calculate steady-state conditions throughout the system and initialize the various variables and constants required for the transient analysis. c) Set up the various column headings for the print out of the results of the calculations and print out the initial steady-state conditions. d) Increment the time variable, T, by At and the integer counter, U, by one. If the lox system is to be analyzed, proceed. If the fuel system only is to be analyzed, the variable K is set equal to 2. e) Calculate the interior points of all the pipelines in the system being analyzed.

f) Caltulate- the boundary cornitiions. atf-each pipe in the system. This includes all the series pipeline connections, the pump, the pulser or PVC, and the feed tank. g) If the lox system was analyzed and the fuel system is going to be analyzed, the program returns to step (e) and repeats the steps for the fuel system. If the fuel system is not to be analyzed the program continues with the next step. h) Calculate the engine boundary condition. i) Substitute all newly computed values of pressure and discharge to a second storage location so they can be reused during the next time increment. j) Compute the pump characteristics and constants to be used for the computations of the next time increment. k) If a print out is desired, transfer to the print statement, otherwise, return to step (d) and repeat the calculations at the new time. A check is made to see if the time exceeds the desired duration of the transient calculations. If so, the program is stopped. C. Initial Steady State Conditions Data are given for the number of reaches in one particular pipeline, as well as its length and pressure pulse wave speed. These data are used to calculate the time increment, DT, to be used for the transient calculation. The number of reaches in each of the remaining pipes in the system are then calculated, based on this DT. A number

-15 of other constants that are needed in later calculations are then established and stored for future use. At any location in the system where identical parallel pipes exist, they are treated as a single pipeline having an area equal to the combined area of the pipes but the actual pipeline diameter is maintained. An example of such a situation would be the two pipes labeled number eight in the pump discharge section of the lox system, Fig. 1. With this representation the velocities are correctly represented in the pipe and the frictional losses are also correctly described. This substitution provides a significant simplification in the program. The given data, in addition to the physical description of the system, that define the initial operating conditions are: the flow in each system, the pressure in each feed tank, and the pressure in the combustion chamber. The factor DPF, which corresponds to the change in hydraulic gradeline in each reach of each pipeline, is computed; then the pressure at the pump outlet is computed by beginning at the combustion chamber and working through the orifice and discharge piping (see Figure 1). All steady state operating conditions (Q,p) are now stored for each point at which transient computations will be performed in the system. The pressures in the suction side are computed by beginning at the feed tank pressure and proceeding towards the pump. When the pump suction pressure is calculated, the pressure rise through the pump (PPUMP) is computed as the difference between the suction and discharge pressures. Since it is unlikely that this value of pressure head rise and given initial discharge will agree exactly with the given pump characteristic curve data, the amount (DPO) which the head

-16 rise vs. suction pressure should be shifted to produce the desired agreement is computed. That is, the pump speed is assumed to be changed slightly to obtain the desired agreement. Also the value of pressure rise, PO, corresponding to the initial discharge value is computed from the given head rise vs. discharge data. This is discussed in more detail under the pump boundary condition. Constants pertaining to the engine operation are also computed. These are discussed in detail under the engine boundary condition. Information dealing with the print out of results is next established in the program. This includes column headings and the print out of some of the newly stored initial conditions. The balance of the program involves the calculations during the life of the transient condition. It,involves repetitive calculations which are performed every At. This is accomplished by performing all calculations to the end of the program, incrementing time and repeating until sufficient results are obtained. The first repetitive calculations involve the interior points along each pipeline. Pressures and flow rates are computed at each time increment by use of the interpolation equations and Eqs. (9) and (10) as described in Section III. At the same time increment as calculations are made at the interior points, computations are also performed at the various boundary conditions. These are discussed in the next section.

V. BOUNDARY CONDITIONS The various boundary conditions are now treated individually in some detail. A. Series Pipe Connection At a junction of two pipes, as pipes 2 and 3, Figure 1, Q2 = P3 and HP2 = Hp3 at any instant. In order to handle all similar series connections the numbers of all adjoining pipes are stored in pairs in the linear array IND2. A simple iteration statement then takes each series connection in sequence and obtains the solution at the new time. Figure 5 shows a typical series connection. Equation (14) is written for pipe X and Eq. (12) is written for pipe J. QPX = C3 - Cx Pp QPJ = C1 + CJ pp Since QPX = QPJ, the solution for pp is C3 - C1 P CJ + CX With pp known, Qp can be found from either of the previous two equations. The pump inlet and discharge connections are treated as special series connections. These are discussed under the pump boundary condition. Also if the excitation mechanism in the system is -17

PIPE J PIPE X FEED I cos (wt+ a,) I H co I Figure 5. Series Connection SUCTION PIPE Figure 6. Connection at Feed Tank

-19 either a structural displacement of the feed tank or pump, an alteration is made to those particular series connections. B. System Excitation at the Feed Tanks A structural motion of the tank with respect to the feed line is assumed if the variable TANK is set equal to one. The amplitude of the sinusoidal motion is AM1, the frequency is OMEGA, and the phase angle is ALPHA1, Fig. 6. The continuity condition at the junction becomes Q1 + P = 2 where QP1 = p (PROP * Al - A2) (AMl) cX sin (wt + a1)' The variable PROP can be used to describe the amount of the tank area that actually contributes additional flow. It would depend upon the shape of the pipe entrance and upon the number of suction lines from each tank. When this continuity condition is used with the earlier boundary condition equations, the pressure at the junction can be determined. C3 - C1 + QP1 pp =- - Cj + CX C. Turbopump Representation In the analysis the pump is considered to be a pipe with a length equal to the mean particle flow path through the pump. The head rise through the pump is considered a linear rise along its length and frictional losses are charged to the pump head rise. A pressure pulse wave speed is assumed through the pump which gives rise to a small time delay. Thus, the pump is treated as any other

-20 pipe in the system with interior point computations performed as described above, and it is connected to the suction and discharge lines using the same procedures as described above for series connections. The equation of motion takes a slightly different form when it includes the head rise produced by the pump. If the total pressure PPUMP head developed by the pump is PPUMP, then an additional force P P A Ax L can be considered to act in the direction of flow on the control volume shown in Figure 2. Then the equation of motion becomes PPUMP dV PA A x Ax - pgA Ax sin a + pA Ax = x L dt which reduces to the following form: Px PPUMP _ --- - g sin a + Wx + Vt = 0 p p L When this equation is combined with the original continuity equation in the same manner as in section III, and placed in finite difference form, the following equations are obtained. - QR + a (Pp - PR) - A At (gp sin a + PPUP)= 0 (15) Xp - xR = a (tp - tR) (16) A PPUIAP)_ O (17) P - S - A (P - PS) - A At (gp sin a + PPP)= (17) a L xp - XS = -a (tp - ts) (18) These equations are used for the pump in place of Eqs. (5) to (8). In order to use this set of equations it is necessary that PPUMP be evaluated. The procedure used for this evaluation is next

-21 described. It is assumed that the steady state characteristic curves for the pump are valid during the transient condition. It is further assumed that the shape of the curves remains the same although the pump speed may vary slightly during the unsteady flow. Tabular data of pressure rise vs. suction pressure (PNP vs. PS) are given, starting at PS = Pi, at equal increments of PS, Figure 7. Since initial steady state conditions in the system are not likely to match the given pump curve exactly, the pump speed is assumed to be altered slightly so that the curve is shifted DP. The initial suction pressure and pressure rise conditions in the system then fall on the newly positioned curve. This shifted curve is used throughout the analysis. At any suction pressure, the constants for a parabola (B1, B2, B3), passing through three adjacent points of the given curve, can be determined. The equation of the adjusted parabola is PRISE = B + DPo + B2 (PS - P ) + B (S P)2 Tabular data are also given for the second pump characteristic curve, pressure-head rise vs. discharge at a constant speed, Figure 8 (PQ vs. Q, data given at equal increments of discharge, DQ). Again for any discharge q, the constants can be determined for a parabola (Z1, Z2, Z3) that passes through three adjacent given data points. Thus for the initial discharge Qo, the parabolic equation is Po =Z + 2 (o Ql) + Z3 (o - ) A slight change in the assumed pump speed is likely to be necessary to shift the curve so the actual initial operating conditions fall on the

-22 PRISE INITIAL OPERATING CONDITION GIVEN DATA DP PNP vs Ps 1,1 — - SUCTION PRESSURE Ps Figure 7. Pump Characteristic Curve. Pressure Rise vs Suction Pressure PPUMP Po I I. PRISE - PO DQ QI QO DISCHARGE Q Figure 8. Pump Characteristic Curve. vs Discharge Pressure Rise

-23 pump curve. The amount the curve is displaced is AP = PRISE - P During the transient calculations, for a given suction pressure, B1, B2, and B3 are found and PRISE can be evaluated. The difference PRISE - P is the amount that the pressure rise-discharge curve must be shifted. With values of Z1, Z2, and Z3 determined with the aid of the most current value of discharge, the relationship between the current head rise and discharge is given. PPMP = Z z - + Z ((Q - l)) + RI ( SE P0 It can be noted that this analysis assumes that the pressure rise vs. suction pressure curve is a constant discharge curve as well as a constant speed curve. For convenience and simplicity in the analysis, values of suction pressure and discharge from the previously computed time increment are used to evaluate PPUMP. D. System Excitation at the Pump. A structural motion of the pump with respect to the suction line is assumed if the variable PUMP is set equal to one. The assumed sinusoidal motion has an amplitude, frequency and phase angle of AM7, OMEGA, and ALPHA7, respectively. The continuity condition at the suction line junction with the pump is QP6 + QP6 = Qp7 where QP6 = pA6 c(AM7) sin (wt + a7) When this continuity condition is used with the boundary equations

-24 presented under A, the pressure at the pump inlet can be determined, C3 - C1 + Qp6 pp = ---- CJ + CX E. Pulser If the variable PULSER is set equal to one in the input data, a sinusoidal pulser is assumed at location 5 of Fig. 1, A piston is assumed in a branch at the junction of pipes 4 and 6 as shown in Fig. 9. The motion of the piston is controlled in amplitude, frequency, and phase by the variable AM5, OMEGA, and ALPHA5, respectively. The continuity condition at the junction is given by QP4 + QP5 = qP6 where QP4 is defined by Eq. (14) 4 = C3 - C4 PP QP6 is defined by Eq. (12) qP6 = C1 + C6 Pp QP5 is defined by QP5 = pA5(AM5) o sin (ct + Q5) The resulting pressure at the junction is given by the relation C3 - C1 + QP5 pp C4 + C6

-25 Figure 9. Location of Pulser LOX DOME MANIFOLD I. * -~~~~~~~~~~~ —— ~~~~~~~~~~~~~~~~~~I Figure 10. Schematic Diagram of Engine Combustion Chamber, Burner Plate Orifices, and Pump Discharge Lines

-26 Fo Pressure Volume Compensator If the variable PVC is set equal to one in the input data, the volumetric response of the PVC is assumed at location 5 of Fig, 1o This volumetric response can be visualized as a piston in a branch, as shown in Fig, 9, that moves in response to the physical displacement of the pipeo The variables that control the PVC motion are the same as those controlling the pump motion, namely AM7, OMEGA, and ALPHA7o The equations at the junction are identical to those for the pulser except QP5 is defined P5 - pA6(AM7) C sin (ot + a7) It is not possible to use the pulser and PVC in the program at one time, G. Feed Tank Su-'-ace The pressure in the fuel and lox supply tanks is given by PTANK and is assumed constant during the transient analysis. With the pressure defined at a boundary the discharge can be computed directly from Eq. (12)o Ho Engine for Combined System Analysis The action of the engine beyond the burner plate orifices has been greatly si.mplified. I:t is assumed that the combined mass flow rate of fuel and lox Q'to the engine nozzle can be related to the combustion chamber pressures P~ by the relation t - "C30 PC (19) The varfiable C30 is dependent upon the nozzle throat area and the characteristic velocity C*-, w.hi3ch in turn is a function of the mixture

-27 (2) ratio Mr (flow rate of lox over flow rate of fuel).(2) Since small variations in the mixture ratio may produce significant changes in the characteristic velocity, the variable C3 is evaluated during the transient as the mixture ratio changes. Tabular data, characteristic velocity vs. mixture ratio, are given for the particular fuel and oxidizer used in the system. Values of C*, CSTAR, are given at equal increments, DMR, of mixture ratio beginning at MR1. For the initial combustion chamber pressure, total mass flow rate, and mixture ratio, the constant ATH is evaluated. ATH = C30 C* = C* Qt/Pc During transient conditions C30 is evaluated from the constant ATH and the C* that corresponds to the current mixture ratio. Figure 10 displays a schematic view of important component parts of the engine that are needed to describe its dynamic response in connection with the fluid pipeline. The notation used to describe the variables in the engine is shown. The fluid in each of the feed lines passes into a volume, then through the burner plate orifices to the combustion chamber. The volume V and effective compliance K' of the space at the orifice entrance are specified by the data, VOL and KP, respectively. The effective compliance, which represents the combined effect of the fluid compressibility and elasticity of the container, can be expressed K' = Ap/(AV/V). For a small time interval the mass flow rate going into storage in the volume is QlO = p. At By combining these two expressions Q10 = Clo Ap

-28 where Cl0 = p V/(K'At). As a first order approximation, Ap = pp - p, where pp is the unknown value of the pressure at the end of the time increment and p is the known value of the pressure at the beginning of the time increment. This expression is used in Eqs. (24) and (28) below. In addition to Eq, (19) and the continuity condition in the combustion chamber QL + QF = tt (20) an equivalent set of equations must be written for the fuel and lox sides of the engine. For the fuel side these equations are as follows: The orifice relationship: F= (AOD)F 2PF(PF - Pc) (21) Continuity at the volume: QF FO10 = F9 (22) Characteristic equation from pipe 9: QF9 = C3F C9F PPF (23) Volume-compliance relationship: FlO lOP ~PF - C10~, p~, (24) F10 = C1F PpF C1OF PF (o4 The corresponding equations for the lox system are %L = (AOCD)L 2PL(PPL - Pc) (25)

-29 L+ L10 = QL9 (26) QL9 = C3L - C9L PPL (27) QL10 = COL PPL - ClOL PL (28) At each new time increment Eqs. (19) to (28) must be solved for the ten unknowns, Qt, Pc, QF, QF10, QF9 QL, QL10, QL9, PPF, PPL Inasmuch as two of the equations involve quadratic terms, a direct solution is not possible. The equations are first combined, then an iteration procedure is used to obtain a solution. Equations (19) and (20) are combined to eliminate Qt C30 P = QF + QL (29) Equations (26) to (28) are combined to eliminate QL9 and QL10 L = L3 - L PPL (30) where L3 = C3L + ClOL PL and L4 = C9L + ClOL Equations (22) to (24) are combined to eliminate QF9 and QF10 F = F3 - F4 PPF (31) where F = 3F+ C O PF and F4 = C9F + C1OF Equations (21), (25), (29), (30) and (31) now involve five unknowns Pc' QF' QL' PPF' and PpL. The method of solution involves the assumption of a small perturbation on each of the three pressures, a linearization of the orifice equations, then iterating and correcting the pressures until the perturbation is arbitrarily small. The final values are true solutions of the

original nonlinear equations, Each pressure is expressed as the sum of a base pressure plus a variation from the base pressure. P P + Pc PPF PF +PF p =P + P PL PL PL The orifice equation for the lox system, Eq. (25) can be written L = ORL - LL - PL + PPL' Pc L -P0 = ORL + PL where ORL = (CDkV) 2pL. If the first two terms of a binominal expansion are used the orifice equation can be written qL = L1 PPL' c ) 2(PL -T) or QL = L1 + L2 PPL L2 P' (25a) where L1= ORL 7PL - Pc and L2 = L1/(2(PL - Pc))' Equations (25a) and (30) can be combined and PpL' eliminated. L1L4 + L2(L3 - L4 PPL) L2L4 QT = p c (32) qL'L2 + L4 L2 + L4 P A similar manipulation of the equations for the fuel system produces a comparable equation.

-31F1F4 + F2(F3 - F4 PPF) F2F4 QF = - - -- PC' (33) F2 + F4 F2 + F4 Equations (29), (32), and (33) are combined to give the solution for Pc' P' = L6/L7 (34) where L1L4 + L2(L3 - L4 PPL) FF4 + F2(F3 - F4 PPF) L6 = C+ - C L2 + L4 +FL- 30 c 0 L2 + L4 F2 + F4 and. L2L4 F2F4 L7 = C30 +........ 7 30 L2 + L4 F2 + F4 Equations (25a) and (30) are combined to find PpL'; and PpF' is found similarly. New base values for each of the pressures are now found by combining the previous base value and the newly computed variation. The process is repeated until the variations are very small. In the program, four iterations are used to establish the new pressures. The four unknown discharges can be evaluated directly when the pressures are known. I. Engine when Systems are Isolated at the Combustion Chamber When the variable SINGLE is set equal to one in the input data, either the lox or fuel system can be analyzed alone without the interference of one upon the other. The fuel system analysis assumes that the lox flow rate to the combustion chamber is constant, QLB. When this assumption is made, Eqs. (19), (20), (21), and (31) can be combined and a direct solution can be obtained for QF. F = (- Fl +/F1 + 4F2)/2

-32 where F1 = ORF + 30 and 2'F3 QLB \ Fa ORF (F QJLB) With QF known, a direct solution is available for the other variables, qt Pc P PpF9 and QF9~ A similar procedure is used when it is desired to isolate the lox system. The fuel flow rate, QFB, is assumed constant for this condition. J. System Excitation at the Burner Plate Orifice A structural motion of the engine with respect to the pump discharge lines is assumed if the variable ORIFIC is set equal to one. The amplitude of the sinusoidal motion is AM9, the frequency is OMEGA, and the phase angle is ALPHA9. In the fuel system, the continuity condition between the end of the pump discharge line and the orifice plate becomes QF9 = QF + qF1O QPF9 (22a) where QpF9 = p A9 o(AM9) sin (ct + c9) Similarly, for the lox system, Eqo (26) becomes QL9 = L + LO -PL9 (26a) where QPL9 = A9 u(AM9) sin (wt + a9).

VI. TI-E COMPUTER PROGRAM A copy of a listing of the program is shown in Figure 11. It is written in FORTRAN IV-G Level and was executed on the IBM 360/67 system at the Computing Center, The University of Michigan. The storage required by the program is approximately 42,000 bytes. The execution time of the program is dependent upon the particular data being processed. The data that is necessary for execution of the program is listed. In all double subscripted data, the first subscript refers to the system (LOX = 1, FUEL = 2) and, except for the pump characteristics, the second subscript refers to the pipe number. D(l1,)...D(1,9), D(2,1)...D(2,9) L(1,l)...L(1,9), L(2,1)...L(2,9) A(1,1)...A(l,9), A(2,1)...A(2,9) F(1,1)...F(1,9), F(2,1)...F(2,9) NO(1,1)...NO(1,9), NO(2,1)...NO(2,9) SLOPE(1,l)... SLOPE((19), SLOPE(2, 1)...SLOPE (2,9) RHO(1), RHO(2), FLOW(1), FLOW(2), CDORIF(l), CDORIF(2), KP(1), KP(2), PTANK(l), PTANK(2), PCO, DQ, Q(l), Q1(2), VOL(l), VOL(2), PQ(l,l)...PQ(1,20), PQ(2,l)...pQ(2,20), DP, PI(1), PI(2), PNP(l,l)...PNP(1,20), PNP(2,1)...PNP(2,20), IND2(1)...IND2(12), J2, DMR, MR1, CSTAR(1).. CSTAR(20), TANK, PULSER, PUMP, ORIFIC, OMEGA, AMl(1), AM1(2), AM5(1), AM5(2), AM7(1), AM7(2), AM9, ALPHA4(1), ALPHAl(2), -33

ALPHA5, ALPHA7, ALPHA9, SYSTM, SINGLE, N(2,7), G, TSTOP, UU, PVC, PROP The program prints out pressure and discharge information at the points indicated by P and Q in Figo 12o The pressures are given in psi although all pressure computations are performed in psf in the program. The discharge unit is mass flow rate, slugs/sec. Definitions of all program variables are given in Appendix A.

-35 FORTRAN IV G.LEVEL O, MOD 0 MAIN CSATURN 5 PROPELLANT SYSTEMS INCLUDING PUMPS AND ENGINE CSUBSCRIPT K OR KK=1 REFERS TO THE LOX SYSTEM CSUBSCRIPT K OR KK=2 REFERS TC THE FUEL SYSTEM CSYSTM=1 LOX SYSTEM ONLY IS ANALYZED CSYSTM=2 FUEL SYSTEM CNLY IS ANALYZED CSYSTM=3 BOTH SYSTEMS ARE ANALYZED CSINGLE=1 SYSTEMS ARE ISOLATED AT COMBUSTION CHAMBER CPVC=1 THE PVC IS INCLUDED IN THE ANALYSIS CALL PRESSURES ARE IN PSFA, ALL FLCWS ARE MASS FLOW RATES, SLUGS/SEC CPTANK(K) IS PRESSURE IN SUPPLY TANKS IN PSF CRHO (K) IS THE MASS DENSITY CTHE VARIABLES FOR THE AMPLITUDES AND PHASE ANGLES AT THE VARIOUS CEXCITATICN POINTS ARE AS FCLLOwS (FREQUENCY OF OMEGA) CLOX TANK TANK=1, AMI(1), ALPHA(l) CFUEL TANK TANK=1, AM1(2), ALPHA(2) CPUMP PUMP=i, AM7(1), AM7(2), ALPHA7 CORIFICE ORIFIC =1, AM9, ALPFA9 CLOX PULSER PULSER=1, AM5(1),ALPHA5 CFUEL PULSER PULSER=1, AM5(2), ALPHA5 0001 PARA3 (DC1,DC2,DC3,DTH )=DC2+. 5*DTH* ( DC3-DC1+DTH* DC DC3 +DC1-2. DC2) ) 0002 QRKRQ (iQl CQ2,DH) =DQQ2-UQQ 1) 0003 PRK(DP1,DP2,DTH)=DP2-CTh*( P2-DP1) 0004 QSS(DQ2,'Q3,DTH)=UD2-UTH* (CQ2-DQ3) 0005 PSS(OP2, P3,DTH)=OP2-TF* (DP2-DP3) 0006 INTEGER UUU, XPULSERZORIFICPUMP,TANK,SYSTMSINGLE,PVC 0007 REAL KP,L,L2,L3, t4,L,L L7,L,MR 1 NO 0008 DIMENSION 0(2,10),F(2,10),L(2,10),A(2l10),AR(2,10),THA(2,10),C(2,1 20 ),N2,10) FF(2,10) QO(2,10),DPF 210)tN(2,10),SLOPE{2,10)1ACCEL( 32, 1O) P4(2,10), P (2,22 ),PNP({220) OC09 DIMENSION RHO(2), FLO~(2), AM5(2), AM1{2), AM7(2), COORIF(2), ORIF 2(2), VOL(2), KP(2), C10(2),Qi(2),PO(2), P{2)tPPUMP(2), DPO(2), PT 3ANK(2), ALPHA1(2), INC2(20), CSTAR(20) 0010 DIMENSION Q(2,10,40), WP(2,10,40), P(2,10,40), PP(2,10,40) 0011 NAMELIST/INSET/U,L,A,F,NCSLOPERHO, FLOWCDORIFKPPTANKPCO,DQ,Q1 2,VOL,PQDPPI,PNP, IN02, J2 DMRMR1,CSTAR, TANKPULSERPUMP,ORIFICOM 3EGAAM1,AM5,AM7,AM9,ALHAALP1,ALPHAALPHA7ALPHA9,SYSTMSINGLE,NGT 3STOPUUPVCt,PROP 0012 10 REAG(5,INSET) 0013 wRITE (6,INSET) 0014 U=0 0015 T=O. 0016 uT=L(2,7)/(N(2,7)*A(2,7)) CINITIAL STEADY STATE CONDITICNS IN SYSTEM 0017 DO 60 K=1,2 0018 PC=PCO 0019 DO 20 J=1,9 0020 N( K,J)=L ( K, J) /( DT *A(K,J) ) 0021 IF(J.NE.5) THA(K,J)=N(K,J)*DT*A(K,J)/L(K,J) 0022 AR( K,J)=.7854*U(K,J)** 2*NO(K,J) 0023 -C(K,tJ)=AR(KJ )/A(K, J) 0024 FF(K,J)=F(K,J)*iT/12.*DIK,J)*AR(K,J)*RHO(K)) 0025 Q0K, J)=FLOW(K) 0026 IF( J.EQ.5QO(K, 5)=0. 0027 IF(J.NE.5)DPFK,J )=F(K,J)*L(K,J)*QO(KJ)**2/(2.*D(K,J )*N(K,J)*AR( 2K,J)**2*RHO(K))-RHO(K)*G*L(KJ)*SLOPE(K,J)/N(K,J) Figure 11. Listing of Saturn V Program

-144 0028 20 ACCEL(K, J)=G*4tHK) *DT*-SLOPE (K, J )*AR (K,J) 0029 DPF (K,5)=0. OC30 OkiF (K) =CUiORI f- K) SiT (2.*RHO(K) ) 0031 iC 0(K)=VUL(K)*RHO(K)/(GT*KP(K)) 0032 P(K,8,1 )=PC+(FLU (K)/CtkIFK) )**2+PF( K,9 )*N(K,9)+PF (K,8)*N(K,8) 0033 P(K,1,1 )=P ANK(K 0034 PP( K, 1, 1) =PTA \K K) 0035 Ri.I TE 1 6, 311 ) DT,(N(K, I ), I=1 9), I AR(K, I), I=1,9) 0035. ill FL;,KAT (iHU,4h CT=,Fb.5,8H N(K,I)=,915,11H AREA(K,I)=,9F6.3) 0037 uL50 J=1,S 0038 NN=N(K,J)+1 0039 bC 30 I=1,NN OC40 P(N,JI )=P(K,J, 1)- I-1)*CPF(K,J) 0041 3C wI K,J,I )=Q(K,J ) 0C42 IF J-6) 5b,40T50 0043 4C PPUMP (K =P(K,d, 1)-P(K,f,N(K,6)+1) 0044 DPF (K,7 )=-PPUMP (K )/NIK,7) 0045 50 P(K,J+1,)=P(K, JN(KJ)+ ) 0046 PP(K,9,N(K,9)+1 )=P(Kc,N(K,9)+1) CEVALUATION OF PUMP CHARACTEt ISTIC CkRVE CONSTANTS 0047 Z= ((K, 6,N(K,6) +i)-i 1(K) )/DQ 0048 P1=Pkw(K,Z) 0049 FP=Pw(K,Z+i) 00C50 P3=P(K,Z+2 ) 0051Z3=(P1+P3-2.*P2 / (2.*UO) 0052 Z2=(P2-P1 )/iD-Z3*D;*( 2.*Z-1. ) 0053 Z I= P2-D -*Z*( Z2+ Z3*D3* LZ 0054. =(P(K,,( K,6) +i )-PI (K) )/DP 0055 Pi=PNP(K,Z) 0056 P2=PNP( K,Z+i ) 0057 P3=PNP(K,Z+2) 0058 b= (P +P3-.*P2)/ (2. L*FPCP) 0059 32=(P2-P )/OP-b3*DP( 2.*Z-1. ) 0060 b FL=P2-DP )Z* ( B2+ E3* P* Z) 0061 PO K)=Z +Z2*(Q(K,6,N( K,.t)+ )-I1 lK ) )+Z3*('IK,6N(KI K6)+1)-01()K) )*2 0062 6C PO(r )=PPUMP(K)-d 1 32 *(P(K,6,N( K,6)+ 1 )-PI (K) )-3*(P(K,6,N(K,6)+1 )2PI(K))**2 0C63 -L=FLC (1) 0064 QF=FLL (2) 0C65 T= L+'~ F 0066 L =(QL/ F- MRiI/D IK+1 0067 TH= (GL/ F-MRt1-Z DMR) /iRK 0068 CSTA=PAtAd(CSTA ( Z-1),CSAR(), CSTAR(Z+1), TH) 0069 ATH=CSTA J]'l/PC 0070 C3u=A1-H/C S-A 0071 IF (SINGLE-1 }8, 70,80 0072 70 6F f=k-1F 0073 Lb=QL OC74 CF30=C30 0075 CL30=C30 OC76 dO F4=CiC(2)+C(2,9) 0077 L4=C10( 1)+C( 1,9) 00C78 JJ2=2J2 0079 fR ITE (6,630 0080 630 FURMAT (iHO,8iH LI'UIC OXYGEN SYSTEM 2 FUEL SYSTEM/111H TIME TANK OUT PULSER P 3UMP IN PUMP OUT ORIFICE TANK OUT PULSER PUMP IN PUMP OUT ORIF 4ICE COMBUSTION) Figure 11. (Continued)

-37 0081 320 DO 321 K=1,2 OC82 DC 321 J=l,9 OC83 321 P4(K,J)=P(K,J,N(K,J)+1)/144. 0084 PC4=PC/144. 0085 NN=N(1,4)+1 OC86 I'=N(2,4)+1 OC87 WRITE (6,640) T,( P4 K, 1),P4( K, 4) P4(K, 6 ),P4(K,7),P4(K,9),K=1,2)tPC4 2,(1,2,1),Q(1,4,NN),Q(1,7,1),tQ1,8,1),QL,Q(2,2, 1),}Q(2,4,1),Q(2,7,1 3),t(2,8,1),QFQT 0088 640 FORMAT (1HO,F6.4,3H P=,11F9.1/7X,3H Q=,llF9.2) 0089 330 T=T+DT 0090 U=U+1 0091 IF(T.GT.TSTOP) GO GT 240 0092 DO 120 K=l,2 0093 IF (SYSTM.EQ.2) K=2 CINl1tRIUR POINT COMPUTATIONS 0094 DO 100 J=l,9 0095 NN=N(K,J) 0096 DO 100 I=2,NN OC97 R=(,KR( (K,J, I-1),Q(K,J, I),THA(K, J 0098 Pt=PKR(P(KJ,1-1),P(K,J,I),THA(KJ)) 0099 iS= SS ( K,J, I),Q(K,J,l+ I),THA(K,J)) 0100 PS=PSS(P(K,J,I),P(K,J,I+1A,THA(K,J)) 0101 WP(K,J, )=.5*(,(Q+S+C(K,J)*(IPR-PS)-FF(K,J)*(QR*ABS(QR)+QS*ABS(QSJl 2)+ACCEL( K,J) 0102 PP(K,J, I )=.5*(PR+PS+(F,-.S-FF(K,J)*(QR*ABS(QR)-QS*ABS(QS)) )/C(K,JI 2) 0103 100 IF (J.EQ.I)QP(K,7,I )=P(K,7,I)+PPUMP(K)*AR(K,7)*DT/L(K,7) CSERIES PIPE CONiNCTILNS 1NCLUU1NG PUMP INLET AND DISCHARGE 0104 00110 M=l t,2,2 0105 X=IND2(M) 0106 I=N(K,X)+1 0107 Wk= RR(d(K,X,{ -1),Q(K,X, I),THA(K,X)) 0108 C3=WR+C(K,X)*RR(P( KX, -1),P(K,X,I),THA(KX))-FF(K,X)*QR*ABS(QR)+ 2ACCEL(K,X) 0109 IF (X.EQ.7)C3=C3+PPUMF(K)"KAR(K,7)*DDT/L(K,7) 0110 J=IN02(M+1) 0111 QS='SS((.(K,J,1),Q(K,J,2),THA(K,J)) 0112 Ci=QS-C(K,J)*PSS(P(K,J, 1)P( K,J,2),HA(K,J))-FF(K,J)*QS*ABS(QS)+AC 2CEL(K,J) 0113 IF(J.E'.7) CI=C1+PPUPF(K)*AR(K,7)*UT/L(K,7) 0114 QP6=0. 0115 1 F PUMP.E. 1. AND. X. EC.6) P6=RHO K ) AR(K, 6)*AM7(K *OMEGA*SIN( OMEGA* 2T+ALPHA7) 0116 (P1=0. 0117 IF( TANK.E. 1.ANUD.X. E.1 )P1=RHO(K)*PROP*AR(K 1)-AR(K,2)) 2*AM1 K)*OMEGA*SIN(UMECA*T+ALPHA1(K)) 0118 PP(KX,I)=(C3-C1+CP6+P1)/(C(K,J)+C(K,X) ) 0119 PP(K,J,1)=PP(K, X,1) 0120 QP(K,X, I)=C3-C(K,X)*PP(K,X,I 0121 110 P(K,J, 1) =P{(K,X, I) CPULSER OR PVC JUNCTI N 0122 I=N(K,4)+1 0123 k=QRKR(,(K,4, I-1),Q(K,4,I ),THA(K,4)) 0124 C3='R+C(K,4)*PRR{P(K,4,i-1),P(K,4,I),THA(K,4))-FF(K,4)*QR*ABS(QR)+ 2ACtL (K,4) Figure 11. (Continued)

-38 0125 QS=(iSS(Q(K,6,1),QK,6,2),THA(K,6)) 0126 C1=QS-C(K,6)*PSS(P(K,6,1),P(K,6,2),THA(K,6))-FF(K,6)*QS*ABSlQS)+AC 2CEL(K,6) 0127 QP K,5,1)=0. 0128 I F ( PULSER.EQ. 1)Q P i(K,5, 1)=RHO(K)*AR(K,5)*AM5(K)*OMEGA*SIN(OMEGA*T+A 2LPHA5) 0129 IF (PULSER.EQ.O.ANO.PVC.EC.1) QP<K,5,1)=-RHO(K)*AR(K,6)*AM7(K)*OMEGA 2*SI N(OMEGA*T+AL PHA7) 0130 P(K,4, I )=(C3-C1+QP(K,5,1) )/(C(K,4)+C(K,6)) 0131 PPK,6, 1)=PP( K,4, I) 0132 PP K,5,1) =PP(K,, 1) 0133 tP (K,4, I)=C3-C( K,4)*PP(K,4,1 ) 0134 P (K,6, 1) =C1+C( K,6)*PP(K,6,1) CBOUNDARY CONDITION AT FUEL ANO LOX TANKS 0135 QS=(SS(Q(K<,1,1),Q(K,1,2),THA(K, 1)) 0136 C1=QS-C(K,1 (PK,1,*SS(P(Kl),P(K,2),THA(K,1))-FF(K,1)*QS*QS+ACCEL(K 2,1) 0137 P (K,i, 1)=C1+C(K, 1 )*PF(K, 1,1) 0138 120 IF (SYSTM.EQ.l) K=2 CBOUNOAKY CONDITIION AT ENGINE 0139 I=N 1,9)+1 0140 QR=RR(Q(1,9,I-1),Q(19,9I ),THA( 1,9)) 0141 C3=-R+C(1,9)*PRR( P(1, t 9,I-1), P (1,9I ),THA( 1,9))-FF( 1,9)*QR*QR+ACCEL 2(1,9) 0142 J=N(2,9)+1 014.3 QR=QRR( Q(2,9, J-1),Q(2,9,J),THA(2,9)) 0144 C2=QR+C(2,9)*PRjR(P(2,9,J-1),P(2,9,J),THA(2,9))-FF(2,9)*QR*QR+ACCEL 2(2,9) 0145 L3=C.3+C10 (1 )*P( 1,9, I 0146 F3=C2+ClO(2) *P(2,9,J) 0147 IF (ORIF IC- ) 140,130,140 0148 130 QPL9=RHU( 1) *AR (1, 9) *Ag9*OMEGA*S N (OMEGA*T+ALPHA9) 0149 iPPF9=RHO (2) AR(2,9)* AM9*OMEGA*SIN(OMEGA*T+ALPHA9) 0150 GO TO 150 0151 140 WPL9=G. 0152 QPF9=0. 0153 150 IFISINGLE-1)190,160,190 0154 160 IF(SYSTM.EQ.2) GO TO 170 0155 L 1=ORIF ( 1)**2*( 1./L4+1./CL30) 0156 L2=CRIF(1) *2*((L3+QPL9)/L4-QF8/CL30) 0157 QL=.5*(-LI+SQRT(LI*Li+4.*L2)) 0158 QT=QL+QFB 0159 PC=QT/CL30 0160 PP(1,9, ) =(-QL+L3)/L4 0161 QP(1,9, I)=C3-C(1,9)*PP(1,9,I) 0162 Z= L CL/QFB-MR1 )/UMR 1 0163 TH=IQL/QF B-MR1- ZL*MR) /MR 0164 CL30=ATH/(PARAB(CSTAR(Z-1),CSTAR(Z),CSTAR(Z+1),TH)) 0165 170 IF(SYSTM.EQ.1) GO TO 210 0166 Fi=F ORIF(2)* 2*( 1./F4+1./CF30 ) 0167 F2=CRIF(2)4**2( (F3+'PL9)/F4-QLBt/CF30) 0168 QF=.5*(-F1+SQRT (Fi*F144.*F2)) 0169 QT= F+QLB 0170 PC=QT/CF30 0171 PP (2,9, J)=(-QF+F3 )/F4 0172 QP (2,9, J)=C2-C(2,9)*PP(2,9,J) 0173 Z=( tL8/CF-MR1 )/OMR+1 Figure 11. (Continued)

0174 TH= ( QLB/QF-MR 1-Z*DMR )/DMR 0175 CF30=ATH/ (PARABICST -)CSTARZ- CSTAR(Z),CSTAR(Z+1, TH)) 0176 GO TO 210 0177 190 DO 200 M=l,4 0178 LI=ORIF (1 *SQRT PP( 1,9, )-PC) 0179 L2=L1/(2.*(PP(1,9,I )-C) ) 0180 F1=ORIF(2)*SQRT(PP(2,9,J)-PC) 0181 F2=F1/( 2*(PP(2,9, J)-PC) 0182 L6=-C30*PC+(L 1L4+L2*( L3-L4*PP(1,9, )+QPL9))/(L2+L4)+(F *F4+F2*(F3 2 -F4*PP(2,9,J)+QPF9 )/(F2+F4) 0183 L7=C30+L2*L4/(L2+L4)+FZ*F4/(F2+F4) 0184 PCP=L6/L7 0185 PLP=(L3-L4*PP(1,9, )-LI+L2*PCP+QPL9)/(12+L4) 0186 PFP=(F3-F4*PP(2,9, J)-F1+F2*PCP+QPF9)/(F2+F4) 0187 PP(1,9,I)=PP( 1,9, +PLP 0188 PP(2,9,J)=PP(2,9,J)+PFP 0189.200 PC=PC+PCP 0190 P(1,9,I )=C3-C( 1,9)*PP(1,9 I) 0191 QP(2,9,J')=C2-C(2,9)*PP( 2,9,J) 0192 QL=CRIF(1)*SQRT (PP( 1,9, )-PC) 0193 QF=CIF(2)*SQRT (PP(2,9, J )-PC) 0194 QT=QF+QL 0195 Z=( L/Qf —MR1 )/DMR+1 0196 TH= ( L/QF-MR-Z*DMR) /CMR 0197 C30=ATH/ (PARAB( CSTARiZ-1),CSTAR { Z},CSTAR (Z+1),TH) ) CSUBSTITUTION STATEMENTS FOR PRESSURE AND DISCHARGE 0198 210 DO 230 K=l,2 0199 IF (SYSTM.EQ.2) K=2 0200 DO 220 J=1,9 0201 NN=N(K,J)+1 0202 DO 220 I=1,NN 0203 Q(K,J,I )=P(K,J,I) 0204 220 P K,J,I)=FP(K,J, ) CPUMP CHARACTERISTIC CONSTANTS AND PUMP PRESSURE RISE 0205 Z=(Q(K,6,N(K,6)+l)-Q1(K) )/DQ 0206 P=1PQ(K,Z) 0207 P2=PQ(K,Z+1) 0208 P3=PQ(K,Z+2) 0209 L3=(P1+P3-2.*P2)/ (2.*DQ*DQ 0210 Z2=( P2-P1 )/DQ-Z3*DQ* ( 2.Z-1. ) 0211 Z =P2-DQ*Z* (Z2+Z3*DQ*Z) 0212 Z=(P(K,6,N(K,6)+1)-PI(K))/DP 0213 P1=PNP(K,.Z) 0214 P2=PNP(K,Z+1) 0215 P3=PNP(K, Z+2) 0-216 B3=(P1+P3-2.*P2)/(2.*DP*DP) 0217 b62=(P 2-P1 )/DP-B3*DP*( 2.*Z-1.) 0218 B1=P2-D P*Z*(B2+B3*DP Z) 0219 PPUMP(K)=Z1+Z2*(Q(K,6,N(K,6)+1)-Q1(K))+Z3*c(Q(K,,NK,6(6)+1)-QLIK) ) 2**2+B1+DPO(K) +B2* P (K,6,N(K,6)+1)-PI(K) )+B3*(P( K 6 NK, 6)+1i-PI (K 3) )**2-PO(K) 0220 230 IF (SYSTM.EO.l) K=2 0221 IF(U/UU*UU-U) 330,320,330 0222 240 STOP 0223 END Figure 11. (Concluded)

-40 TANK LOX + FUEI L_ 2 L TAN K LOX SYSTEM K = 1 SYSTM = I FUEL SYSTEM K=2 SYSTM =2 COMBINED SYSTEMS SYSTM = 3 I PULSER 0 LOX PUMP FUEL PUMP 8 LOX DOME 1 FUEL MANIFOLD BURNER PLATE COMBUSTION CHAMI Figure 12. Sketch of System Showing Location of Program Output Quantities P, Q.

VII. SAMPLE PROBLEM In order to demonstrate the use of the program, a sample problem is included which analyzes a transient condition in Titan II missile. The program analyzes the transients that develop in the oxidizer system as a result of a forced motion of the pump. The motion is assumed to be sinusoidal at a frequency of 62.8 rad/sec and half amplitude equal to 0.007 ft. Figure 13 gives a copy of the listing of the input data and Figure 14 shows a portion of the output from the calculations. The effect of the addition of a PVC in the system can be seen in Fig. 15, where a comparison is made between two runs. Both runs are subjected to the same excitation at the pump, one includes the PVC, the other does not. -41

EXECUTION BEGINS &INSET 0= 10.000000, 10.000000 0.41700000, 0.41700000, 0.41700000, 0.42999995 3.5999994, 0.42999995 2.0000000, 1.5000000 483.00000, 985.00000 1750.0000, 3640.0000 2980.0000, 4230.0000 0. 11000000E-01, 0.11000000E-01, 0 0.0 0.0 NO= 1.0000000, 1.0000000 1.0000000, 1.0000000 1.0000000, 1.0000000 1.c000000, 1.000000, 1.0000000, 0.0 0.0,RHO= 2.8199997 KP= 100000000. 100000000. 5.COOOOOO,VOL= 0.0 189000.00, 224000.00 166000.00, 188000.00 1.0000000 1.0000000 1.0000000, 1.0000000 1.0000000, 1.0000000 176000.00, 190800.00 176000.00, 200000.00 176000.00, 203000.00 176000.00, 204000.00 1.000000, 1.0000000 IND2= 1, 2,, 0.57499999 0.57499999 0.0 1.5000000 2. 5000000 1750.0000 650.00000 0.0 0.72999954E-01, 1.C0000000 ~ 1.0000000 1. 000000 1.0000000 1.0000000 0.0, 1.7599993,PTANK= 10500., 0.0 184000.00 159000.00 1.0000000 1.0000000 1.0000000 176000.00 176000.00 176000.00 176000.00 1.0000000 2,, 0.50000000, 0.57499999, 0.50000000 0.57499999, 0.50000000 0.50000000, 0.57499999, 0.50000000,0.41700000, 0.42999995 0.0,L= 8.0000000, 8.0000000, 23.000000, 0.39999998 1.6699991, 0.0 0.0, 0.50000000, 0.50000000 6.0000000, 3.0000000, 14.000000, 0.0, 0.0,A= 138.00000, 1100.0000, 138.00000, 50.00000 138.00000, 138.00000, 2980.0000, 3000.0000, 2980.0000, 4230.0000, 0.0,F= 0.79999976E-02, 0.79999976E-02, 0.11000000E-01, 0.11000000E-01, 0.12099999, 0.11000000E-01, 0.11000000E-01, 0.72999954E-01, 0.12099999 0.35999998E-01, 1.0000000, 1.4799995, 0.0, 0.0 y 1.0000000, 1.0000000, 1.0000000, 1.0000000, 1.0000000 1.0000000, 1.0000000 1.0000000 1.0000000 1.0000000, 1.0000000,SLOPE= 1.0000000, 1.0000000, 1.0000000 1.0000000, 1.0000000, 0.0, o.0, 1.0000000 0.0, 0.0, 0.0, 0.0, 0.0,FLOW= 17.000000, 8.1699991,CDORIF= 0.36699999E-01, 0.20999998E-01,,000, 2460.0000,PCO= 115000.00,DQ= 1.0000000,QI= 12.000000,PQ= 197000.00, 231000.00, 193000.00, 229000.00, 218000.00, 179000.00, 210000.00, 173000.00, 200000.00 174000.00, 150000.00, 170000.00, 141000.00, 165000.00 1.0000000, 1.0000000, 1.0000000, 1.0000000 1.0000000 1.0000000, 1.0000000, 1.0000000, 1.0000000, 1.0000000 1.0000000 DP= 250.00000,PI= 10500.000, 2000.0000,PNP= 194500.00, 176000.00, 196800.00, 176000.00, 198800.00, 201000.00, 176000.00, 201900.00, 176000.00, 202500.00 203400.00, 176000.00, 203700.00,176000.00, 203900.00 204100.00, 1.0000000, 1.0000000, 1.0000000, 1.0000000 1.0000000, 1.0000000, 1.0000000, 1.0000000, 1.0000000 3, 3, 4, 6, 7, 7, 8, I I 8, 9, O0 O, 6,DMR= 0.999999964E-01,MR1= 1.7999992 5650.0000, 5610.0000 5565.0000 0.0 0.0, 0.0 TANK= O,PULSER= OPUMP= AM5= 0.0, 0.0,AM7= 0.29999 ALPHAS= 0.0,ALPHA7= 0.0,ALi 0, O0, 0 O, 0 0, 1, 1 O, TSTOP= 0.63999951E-01,UU= 8,PVC= &END 0,,CSTAR= 5775.0000 0.0, 0.0 0.0, 0.0,ORIFIC= 998E-02, 0.0 PHA9= 0.0,! 0, O0 O,PROP= 0.19999999 O0 O0, 5756.0000, 0.0, 0.0 OOMEGA= 62.79 0,, 5730.0000, 0.0, 0.0 )9988,AMl= 0.0,ALPHAl= 0.0 1,SINGLE= O, O0 A ^.* Ut Ut J, 5695.0000, 0.0, 0.0, 0.0, 0.0 1,N= O, 0, O, 0, G= 32.199997 2= y I I 9,AM9= 0.0 SYSTM= O0 O0 *0, 0, Figure 13. Input Data for Example Problem

DT= 0.00050 N(K,I)= 33 26 OT= 0.00050 N(KI)= 16 5 6 4 0 6 24 0 I 1 1 2 AREA(K tI)=78.540 7 1 2 6 AREA(K,I)=78.540 FUEL SYSTEM 0.260 0.260 0.260 0.137 0.260 0.260 0.137 0.137 0.196 0.196 0.196 0.137 0.196 0.196 0.145 0.145 TIME 0.0 P= 0.0040 P= Q= 0.0080 P= 00120 P= 0.0120 P= Q= 0.0160 P= Q= C.0280 P= 0.0320 P= 0.0240 P= Q= 0.0280 P= 00520 P= 0.0320 P= 0.0360 P= Q= 0.0400 P= Q= 0.0440 P= Q= 0.0480 P= g= 0.0520 P= Q= Q= 0.0640 P= 0-= LIQUID OXYGEN SYSTEM PULSER PUMP IN PUMP CUT TANK OUT 78.0 17.00 78.0 17.00 78.0 17.00 78.0 1-7.00 78.0 17.00 78.0 16.99 78.0 16.97 78.0 16.94 78.0 16.91 78.0 16,88 78.0 16.86 78.0 16.85 78.0 16.84 78.0 16.84 78.0 16.84 78.0 16.86 78.0 16.90 92.0 17.00 92.3 16.98 92.9 16.96 93.5 16.94 94.1 16.92 94.7 16.91 95. 1 16.91 95.4 16.91 95.5 16.92 95.4 16.93 94.8 16.94 93.6 16.96 92.2 16.97 90.6 16.99 88.9 17.00 87.3 17.02 85.9 17.03 92.0 17. 00 92.4 16.97 92.9 16.95 93.5 16.93 94.1 16.92 94.7 16.91 95.1 16.91 95.4 16.91 95.5 16.92 95.3 16.93 94.8 16.95 93.6 16.97 92.1 16.99 90.5 17.01 88.8 17.02 87.2 17.03 85.9 17.05 1314.4 17.00 1314.2 17.00 1314.0 17.00 1314.0 17.00 1314.1 17.00 1314.2 17.00 1314.4 17.00 1314.6 17.00 1314.8 17.00 1315.0 17.01 1315.1 17.01 1314.7 17.00 1314.3 17.00 1313.7 17.00 1313.1 16.99 1312.5 16.98 1312.0 16.98 ORIFICE 1062.8 17.00 1062.7 17.00 1062.6 17.00 1062.5 17.00 1062.5 16.99 1062.6 17.00 1062.7 17.00 1062.8 17. 00 1C62.9 17.00 1063.1 17.00 1063.2 17. 01 1063,1 17.01 1C62,9 17.00 1062.6 17.00 1062.2 16.99 1061.8 16.98 1061.5 16.98 TANK OUT 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 PULSER PUMP IN PUMP OUT 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 19.8 19.6 1401.2 8.17 8.17 8.17 ORIF ICE 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 1097.2 8.17 COMBUSTION 798.6 25.17 798.6 25.17 798.5 25.17 798.5 25.17 79,1.5 25.16 798.5 25.17 798.6 25.17 798.6 25.17 798.7 25.17 798.7 25.17 798.8 25.18 798.8 25.18 798.7 25.17 798.5 25.17 798.3 25.16 798.2 25.15 798.0 25.15 I I 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 20.2 8.17 19.8 8.17 19.8 8.17 19.8 8.17 19.8 8.17 19.8 8.17 19.8 8.17 19.8 8.17 19.8 8.17 19.6 8.17 19.6 8.17 19.6 8.17 19.6 8.17 19.6 8.17 19.6 8.17 19.6 8.17 19.6 8.17 1401.2 8.17 1401.2 8.17 1401.2 8.17 1401.2 8.17 1401.2 8.17 1401.2 8.17 1401.2 8.17 1401.2 8.17 Figure 14. Output from Analysis with Oscillating Pump

110 (f) 0Q 0r 0: U) Wn U o cu (:) 03 03. 100 90 80 70 I I 0.05 0.1.15.20.25 TIME, SEC.30 Figure 15. Pressure at Pump Inlet of Oxidizer System of Titan II Due to Pump Oscillation,' = 62.8 rad/sec.

VIII. SUMMARY This report includes the digital computer program for the analysis of the propellant systems of Saturn V. Much of the background material necessary to understand the computations performed by the computer during execution of the program is included. Data for Saturn V have not been used in any trial runs of the program as complete data on the geometric properties of the missile are not available. The program was cleared using Titan II data. Computed results have been checked against experimental results obtained (3) on Titan II and reported in the literature, (3) and favorable comparisons have been made. It is our intention that the program can be used directly for any one of the three stages of Saturn V without further modification. However, since the exact configurations are not available, minor changes may be necessary to accommodate the particular data.

APPENDIX A DEFINITION OF PROGRAM VARIABLES ACCEL(k, j) ALPHAl(k) ALPHA5 ALPHA7 ALPHA9 AMl(k) AM5(k) AM7(k) AM9 AR(k,j) A(k,j) ATH B1,B2,B3 C10(k) C1 C2 C30 C3 CDORIF(k) Constant associated with vehicle acceleration in pipe j of system k Phase angle associated with the feed tank motion in system k Phase angle associated with the pulser motion Phase angle associated with the pump motion Phase angle associated with the engine burner plate motion Amplitude associated with the feed tank motion in system k Amplitude associated with the pulser motion in system k Amplitude associated with the pump motion in system k Amplitude associated with the engine burner plate motion Cross-sectional area of pipe j in system k Pressure pulse wave speed in pipe j in system k Constant associated with the engine Variables associated with the pump head rise-suction pressure curve Constant in system k describing the capacitance of fluid and chamber at the engine orifice plate inlet Variable associated with the C- characteristic Variable associated with the C+ characteristic Variable relating combustion chamber pressure and total propellant mass flow rate Variable associated with the C+ characteristic Constant used to describe the discharge coefficient and area of opening of the engine orifice in system k -46

-47 CF30,CL30 C(k,j) CSTAR(i) CSTA DMR DPO(k) DPF(k,j) DP DQ D(k,j) DT F1,F2,F3,F4 FF(k,j) FLOW (k) F(k,j) G IND2(i) I J2 Variables used with combustion chamber boundary condition Constant for system k, pipe j, AR(k,j)/A(k,j) Variable used for storage of tabular data of the engine characteristic velocity vs. the mixture ratio Variable identifying the characteristic velocity at a particular mixture ratio Interval of mixture ratio values at which the characteristic velocities are stored beginning at MR1 Adjustment of the tabulated pump data, head rise vs. suction pressure, to match initial flow conditions Steady-state frictional head loss per section in pipe j in system k Interval of suction pressure values at which the pump head rise values are stored beginning at PI(k) Interval of discharge values at which the pump head rise values are stored, beginning at Ql(k) Pipe diameter, pipe j, system k Time interval at which computations are performed Variables used in connection with the fuel supply to the engine Constant associated with the friction losses in pipe j of system k Given steady-state mass flow rate in system k Darcy-Weisbach friction factor in system k, pipe j Acceleration in the missile axial direction at the time in flight when the transient is being investigated Set of integer constants that identify the pipes, in pairs, by number, that are connected in series Integer variable, generally refers to the pipe section number under consideration Constant identifying the number of series connections in each system

JJ2 J KP(k) K L1,L2,L3,L4 L6,L7 L(k,j) MR1 M NO(k,j) N(k,j) NN OMEGA ORIFIC ORIF(k) PO(k) PlP2,P3 P4(k,j) Constant = 2(J2) Integer variable, generally refers to pipe under consideration Constant defining the effective bulk modulus of volume at the engine orifice plate entrance Integer variable, generally refers to the system, K = 1 refers to lox system, K = 2 refers to fuel system Variables used in connection with the oxidizer supply to the engine Variables used with the boundary condition at the engine Pipe length, pipe j in system k Lowest value of mixture ratio at which characteristic velocities are stored Integer variable, used in series connections and engine boundary condition Integer constant identifying the number of identical pipes connected in parallel that can be treated as a single pipe j in system k Number of reaches into which pipe j of system k is divided Integer variable used in place of N(k,j)+l Angular frequency of excitation Integer used to include or exclude excitation at the engine orifice plate Constant used with the orifice in system k Pump head rise associated with initial steady state discharge, as determined from the given head rise vs. discharge data in system k Variables associated with the pump head rise-discharge curves Pressure at outflow end of pipe j, system k, in psi; used for print out

-49 PCO PCP PC PFP PI(k) PLP PNP(k,i) PP(k,j,i) PPUMP(k) PQ(k, i) PROP PR,PS P(kj,i) PTANK(k) PULSER PUMP PVC QO(k,j) Q (k) Initial steady state pressure in the combustion chamber Variation in combustion chamber pressure at new time increment Pressure in combustion chamber Variation in pressure upstream of orifice plate in fuel system at new time increment Lowest value of suction pressure at which pump head rise values are given in system k Variation in pressure upstream of orifice plate in oxidizer system at new time increment Values of pump head rise stored for equal increments of suction pressure, DP, in system k Newly computed values of pressure at section i, pipe j, system k Pressure head produced by the pump in system k Values of pump head rise stored for equal increments of discharge, DQ, in system k Constant used to define the proportion of the feed tank area to be used when the system excitation is occurring at the tank Interpolated value of pressure head in the pipeline Pressure head at section i, pipe j, system k Pressure level maintained in feed tank of system k Integer used to include or exclude excitation of the system using the pulser Integer used to include or exclude excitation of the system using the pump motion Integer used to include or exclude the pressure volume compensator in the system Initial steady state flow in pipe j of system k Lowest value of discharge at which pump head rise values are given in system k

-50 QFB QF LB QL QP1 Qy6 QPF9 QPL9 QP(k, J~ i) QRQS Q(k,j,i) QT RHO(k) SINGLE SLOPE(k, j) SYSTM TANK THA(k,j) TH T Constant flow rate of fuel to combustion chamber when lox system alone is being analyzed Flow rate of fuel to combustion chamber Constant flow rate of lox to combustion chamber when fuel system alone is being analyzed Flow rate of lox to combustion chamber Variable used when the system excitation is considered to be at the supply tank Variable used.when the system excitation is considered to be at the pump Variable used in the fuel system when the system excitation is considered to be at the engine orifice Variable used in the lox system when the system excitation is considered to be at the engine orifice Newly computed values of flow rate at section i, pipe j, system k Interpolated value of flow rate Flow rate at section i, pipe j, system k Total flow rate to the combustion chamber Fluid mass density in system k Integer constant used to isolate the systems at the combustion chamber Constant used to indicate whether a pipe is in the axial direction of the missile or normal to the axial direction Integer constant used to specify what systems are being studied Integer used to include or exclude excitation of the system using the motion of the supply tank Constant used in pipe j, system k; pertains to the characteristics grid-mesh ratio Variable used with parabolic interpolation Time

-51 TSTOP Constant identifying the time at which the transient calculation should stop U Integer counter, incremented by one when time is incremented by DT UU Integer constant controlling the print out VOL(k) Constant defining the volume at the engine orifice plate entrance X Integer used to identify the inflow pipe number at series connections Z1,Z2,Z3 Variables associated with pump head rise-discharge curve Z Integer used with the pump characteristics Function References PARAB Function used for parabolic interpolation PRR Function used for linear interpolation of pressure at R PSS Function used for linear interpolation of pressure at S QRR Function used for linear interpolation of discharge at R QSS Function used for linear interpolation of discharge at S

REFERENCES 1. V. L. Streeter, and E. B. Wylie, Hydraulic Transients, McGraw-Hill, Inc., 1967. 2. R. H. Fashbaugh, "Liquid Propellant Missile Longitudinal Oscillation," The University of Michigan, NASA Report, Nov., 1966. 3. V. L. Streeter, and R. H. Fashbaugh, -'Resonance in Liquid Rocket Engine Systems," Trans. ASME, Jour of Basic Engineering, Vol 87, (Dec., 1965), pp. 1011-1017. -52

APPENDIX C ANALYSIS OF CAVITATION EXTENT P. S. Larsen

TABLE OF CONTENTS Page LIST OF FIGURES............................................... A. INTRODUCTION......o...1o...ooo...... B. PHYSICS OF CAVITATION.............................. 4 1. Effect of Gas Contents on Pressure Wave Attenuation.. 8 2. The Cavitation Process................................ 10 a. Cavitation Inception............................. 10 b. Cavitation Bubble Size History................... 11 c. Gas Pressure in Cavitation Bubble............... 14 d. Cavitation Bubble Size Distribution.............. 14 e. Cavitation Void Fraction.......................... 15 C. PROPELLANT PUMP CAVITATION AND SUCTION LINE RESONANT FREQUENCY.......o o.................................... 17 1. Backgroundro O......o......................... 17 2. Model Relating Suction Line Pressure and Resonant Frequency......................................... 18 3. Example of Data Correlation.......................... 24 REFERENCES................................................... 27 ii

LIST OF FIGURES Figure Page 1o Schematic for Hypothesized Pressure Amplitude for Fundamental Harmonic with and without Cavitationo Cavitation Experiment.................. 2, Cavitation Apparatus Comprising L-shaped Plexiglas Tube Fitted with a Piston Driven by the Mechanical Vibrator oooo............... oo ooooo ooooo oooo ooo o. 6 3. Close-up of Plexiglas Tube with Pistono No CavitationOo.OOoo... o...............oooooooooooooo 4o Close-up of Plexiglas Tube with Piston. Visible CavitationoOo. o. 000000 0000 00000000 0 0 0 0. 7 5. Cavitation Void Distributions in Suction Line...oo...o 19 60 Titan II, Stage I, Lox Suction Resonant Frequency..... 22 7. Model Including Inducer Cavitation (schematic).,.,.... 23 80 Titan II, Stage I, Lox Suction Line Resonant Frequency versus Cavitating Length Lc and Average Void Fraction o o... o o,, o..... o................. 26 iii

A. INTRODUCTION Pump cavitation in high speed propellant pumps of the inducer type has been observed to occur not only in the low pressure region near the pump inlet and in the inducer, but also several pipe diameters upstream into the suction line measured from the inducer inlet.) Although the mean suction pressure measured at the inducer inlet is considerably higher than the vapor pressure of the fluid (eogo 120-140 psia for the Saturn C-5 LOX line where the saturation pressure is about 4~1407 psia) cavitation may still occur in this region if pressure fluctuations about the mean value of sufficient amplitude exist. Two mechanisms, a random and a periodic, by which such pressure fluctuations can be created are postulated, both originating at the impeller inlet. First, the high turbulence level and secondary fluid motion created at and upstream from the inducer inlet may be visualized as a random source of pressure fluctuations about the prevailing mean pressure. These pressure fluctuations are also propagated upstream and may give rise to cavitation locally and in a region upstream of the inducer of extent determined by the attenuation of the pressure amplitude. Second, the pressure difference existing between the front and back sides of the moving inducer blades near their leading edges as determined from theblade loading may be visualized as radially distributed, rotating pressure disturbances. Relative to a fixed point in the suction line the inducer inlet may therefore be interpreted as a -1

-2 periodic fluctuating pressure source, having a frequency equal to the inducer blade number times its RPM, and an amplitude related to blade loading. While both mechanisms may be important in causing cavitation the latter is presently favored because it is amenable to analysis. Theoretically, given the imposed fluctuating pressure distribution at the inducer inlet, the pressure field upstream into the suction line and the cavitation extent could be predicted provided proper account could be made of nucleation, growth and collapse of cavitation bubbles in terms of an average local void fraction, attenuation and speed of pressure waves in the cavitating zone. In support of the second mechanism described above an estimate is made of possible pressure amplitudes at the inducer blade tip based on the non-conservative assumption that the fluid near the blade is accelerated to the blade speed. With R the radius of inducer, taken as 0.75 ft, and with fluid density estimated at 2.15 slugs/ft3 (liquid oxygen at -290~F) N = 6000 RPM (Tp) - p v2 /2 TIP- TIP VTIP R2xTN 3 26ooo = 470 ft/sec TIP6o 6o (Ap)TIP 2.15 (470)2/(2 x 144) = 1650 psi Although only a fraction of this velocity is attained by the fluid it The presence of cavitation upstream of the inducer inlet may cause inlet flows which deviate significantly from the design in terms of velocity diagrams, hence leading to unfavorable blade loadings.

-3appears entirely possible that a pressure amplitude adequate for initiation of cavitation in a subcooled fluid at 14o psia can exist.

Bo PHYSICS OF CAVITATION The cavitation extent experiment described in Chapter IV simulates a simplified version of the above model. Instead of a radially distributed rotating pressure amplitude imposed on a fluid column the experiment employs a simple, one-dimensional plane pressure source achieved by oscillating a piston in a fluid filled pipe as shown schematically in Figure 1. The experimental set-up is shown in Figure 2. Figures 3 and 4 show a close up of the region near the piston with and without cavitation, respectively. Any analysis of even the simple plane wave case presents a very complex problem in which a coupling exists between the following phenomena, (i) cavitation bubble nucleation, growth and collapse (ii) wave propagation and attenuation in a bubbly mixture of dispersed transient cavitation bubbles in the liquid. The nucleation of cavitation bubbles and their time history depends on the pressure field which in turn depends on the fluid properties in terms of cavitation void fraction and bubble size distribution and the frequency spectrum of the pressure field itself, Although the ultimate objective in this study is only to relate macroscopic observables associated with the cavitation extent (such as pump suction pressure and suction line wave speed) a brief discussion of the various aspects of the cavitation process is given first as a means of introduction. -4

-6. Figure 2* Cavitation Apparatus Comprising,L-shaped P3lexiglas Tube Fitted with a Piston Driven by the Mechanical ViWbrator

7" - Figure 3, Close-up of Plexi.glas Tube with P.iston. Cavi'tatZion. No1 Figure I4 Close-uup of Plexigl.as Cavit at:ion T1 ube with P:ist on. Vis ible

-8 1. The Effect of Gas Content on Pressure Wave Attenuation For a homogeneous single phase viscous fluid the attenuation of pressure waves arises from viscous effects near confining solid surfaces and viscous effect in the bulk of the fluido The latter effect depends on wavespeed, fluid viscosity and the frequency of the pressure oscillations, but is negligible for most practical purposes where frequencies are low compared to the ultrasonic khz and mhz range. In the absence of cavitation, therefore, only wall friction is considered in regard to attenuation, and this effect is essentially negligible for short pipe lengths such as those being considered presently for the cavitation extento The presence of gas or vapor bubbles in a liquid, such as occurs in cavitation, affects both wavespeed and bulk attenuation in a significant way. The strong reduction of the wavespeed for small gas bubble void fractions has been discussed in Chapter IV of this report, and is a direct result of the change in fluid compressibility and complianceo The expression given is valid for low frequency pressure waves Ignoring for the present purpose the effect of wall friction on the attenuation of pressure wave amplitudes, only the effect of friction in the bulk of the fluid is of importanceo The significant effect on attenuation observed in the presence of gas or vapor bubbles is caused by the general compressibility of the gas-liquid mixture --- allowing substantial amounts of compression and expansion work to be converted into internal energyo Since the gas bubbles present in the liquid are brought to vibrate at the imposed frequency of pressure

-9 oscillations a particularly strong attenuation occurs when the imposed frequency of oscillation coincides with the resonant frequency of the individual bubble. When all bubbles present are of the same size a selective attenuation is obtained at a particular frequency. When the bubble size distribution is broad, however, significant attenuation over a range of frequencies results. The bubble resonant frequency for the zeroth order zonal harmonic is, neglecting surface tension, 1 37Fo fr =1 2 (1) r 2irR p where 7 is the polytropic exponent describing the compressionexpansion process of the gas, PO and p denote the liquid pressure and density, respectively. The propagation and attenuation of sound waves in the frequency range of 60 hz - 20 khz through bubbly air-water mixtures in a standing wave tube have been reported by Silberman, showing agreement between theory and experiment. The attenuation may be expressed in terms of a coefficient of attenuation p which describes the exponential decay of the maximum half amplitude of plane waves propagating into a semi-infinite homogeneous fluid Ap(x) = Ap(0) e- Px (2) Both wavespeed and coefficient of attenuation depend on void fraction, bubble size and frequency of imposed oscillations, showing strong resonance characteristics at the resonant frequency of the bubbles when these are of the same size. At moderate to low void fractions of

-10 a = 002 - o01 the magnitude of the void fraction has little effect on the attenuation coefficient, suggesting values about 100 db/ft = 11.5 ft1 in a range of frequencies a decade above the frequency of bubble resonanceo Below the resonant frequency p drops quickly to small values. A typical average value would be about 1 ft o At frequencies low compared to the resonant frequency of the bubbles the bulk attenuation is negligible. Because of the small size of cavitation bubbles the expected bubble resonant frequency in the pump cavitation zone is in the khz range. Bulk attenuation is therefore without import for the low frequency pressure fluctuations characteristic of the overall system resonance, but of great significance for the high frequency (periodic or random) pressure fluctuations responsible for the cavitationo In fact, the bulk attenuation limits the extent of the cavitation zoneo In regard to wavespeed reduction for a given void fraction a similar frequency dependency (2) exists as reported by Silberman(2) This, however, has only a secondary effect on the cavitation extent and need not be consideredo 2, The Cavitation Process ao Cavitation Inception The major obstacle in the analytical and experimental study of cavitation phenomena lies in the problem of cavitation inception. Physically this phenomenon is related to a random nucleation process and depends on a large number of parameters, such as pressure, temperature and thermal properties, size distribution of potential nuclei, gas and impurity content, and turbulence level, No general success has

-11 yet been obtained in the prediction of cavitation inception and extent in flow machinery based on simple dimensionless parameter groups such as the cavitation number or the suction specific speed(3). The explanation may be found in the large number of parameters governing the process. To facilitate any analysis of the extent of cavitation zones resulting from a given pressure disturbance it is necessary to relate in a simple manner the cavitation void fraction or its time-rate of change to the local static fluid pressure relative to the vapor pressure. It is common practice to assume that cavitation will begin at any point in the liquid when the local pressure reaches an arbitrarily assigned rupture or breaking pressure of the fluid. Typical values of -10 psia (4) to -100 psia have been used in bulk cavitation studies for sea water' b. Cavitation Bubble Size History The average cavitation void fraction depends on the time history of the individual bubbles. The problem of growth and collapse of a single bubble in a harmonic pressure field was formulated by Noltingk and Neppiras(5) who also obtained numerical solutions. The analysis assumes the presence of "nuclei" of sufficient size that growth can occur during the period of decreasing pressure. Subsequent collapse may and may not be completed during the following period of increasing pressure depending on the initial size of the "nuclei", the pressure amplitude and level. The transient bubble size distribution in a liquid cavitating owing to a harmonic pressure field is therefore very complex. For the nuclei that complete growth and collapse over

-12 a period of pressure fluctuation their average size over the period is nearly the maximum size since both the growth and the collapse takes place very rapidly. Other bubbles, however, may merely pulsate about a mean size over a period of pressure oscillation. A simple estimate of the effect of suction pressure level on the void fraction may be obtained by considering Rayleigh's solution(6) to the collapse of a single spherical bubble of initial radius Ro in an incompressible liquid resulting from a sudden imposed pressure difference between liquid and saturated vapor APnvo The collapse time to is found to be to p R 915 (3) where p is the liquid density. This analysis ignores surface tension and thermal effects of importance in vaporous cavitationo If the harmonic pressure variation be replaced by a square wave, representing a step change in imposed pressure for the collapse stage, the maximum bubble size Ro is proportional to nApv for a given collapse time which in some manner will depend on the imposed frequency of pressure fluctuationso Therefore, for a constant number of active sites for bubble nucleation it may be expected that the local average void fraction, a c R o3 will depend upon the local liquid pressure as (see Figo 1) 3/2 3/2 a o (Apv) = [Pv - ps + Ap] (4) For a constant local pressure amplitude Ap the local void fraction thus depends on suction pressure ps as

-13 = C1 [C2 - p3/2 (5) where C1 and C2 are constants. The results of the theoretical study of the bubble growth phase by Noltingk and Neppiras(5) suggest a nearly linear increase in maximum bubble radius R with increasing negative liquid pressure amplitude ApQv for constant frequency of the imposed harmonic pressure field, and a practically linear increase in maximum bubble radius Ro with increasing period of pressure oscillations for constant amplitude. The latter is in agreement with Eq. (3) for bubble collapse if to is proportional to the period. If the important parameters for the problem are, a, p, Ro, D, cu, dimensional analysis, however, leads to a form in agreement with Eq. (4), 3 3/2 3/2 3D 3 3/2 /2 p0 = 3 pVIP2 (7) and the cavitation index K is defined by K = (P - PV)/(PVTI/2) (8) Eq. (6) becomes = C4 [C3 P - K]32 (9) In correlating suction line resonant frequency to suction pressure it was found convenient to use a functionally similar, but more simple,

-14 relation for the mean cavitation void fraction a as function of the suction pressure c- ~C5 -- 1 (10) Ps - Pv K Co Gas Pressure in Cavitation Bubbles An important parameter in the description of pressure wave propagation in bubbly mixtures is the gas pressure in the bubbles. As outlined in Chapter IV the wavespeed depends strongly on the value assumed for the gas pressure. A similar effect may be expected for the f2) pressure wave attenuation(2o The experimental data of reference (2) were obtained for air bubbles rising through a water columno In regard to the propagation of pressure disturbances the rising gas bubbles may be considered to be stationary and at the locgal hydrostatic equilibrium pressureo Cavitation bubbles on the other hand are caused by the local pressure fluctuations and furthermore the pressure in a rapidly expanding (7) and collapsing cavitation bubble varies considerably during its lifetime(7)o If the growth and collapse times for transient cavitation bubbles are small compared to the period of pressure oscillation the main contribution to the average attenuation and wavespeed change comes from bubbles near their maximum size and having a gas pressure which deviates little from the saturation vapor pressure corresponding to the liquid temperatureo In the present study it has been assumed that the pressure of the gas phase in the bubbles equals the saturation vapor pressureo do Cavitation Bubble Size Distribution Tn a recent detailed study of ultrasonically induced cavitationt )

-15 an account is given of experimentally observed bubble size distributions in a cavitation field. The study shows that most bubbles appear and disappear over a period of the imposed pressure oscillations at 20 khz. The majority of bubbles have maximum sizes which are about one half or less that of the resonant size corresponding to the imposed pressure oscillationso This is important information in appraising the coefficient of attenuation which as mentioned above is a strong function of bubble size for a given frequency. The validity of superposition of the effect of several bubble sizes on the attenuation coefficient for a range of frequencies is suggested to a first approximation by Figure 10 of Reference 2 involving two bubble sizes. Conversely, for a fixed frequency the effect of a range of bubble sizes may be evaluated by weighing the attenuation coefficient as function of bubble size with the actual bubble size distribution. e. Cavitation Void Fraction Apparently no thorough experimental investigation has been made yielding quantitative data on the magnitude and distribution of cavitation void fractions in turbopumps, although visual observations have appeared. A single study(9) reports on the measurements of cavitation void fraction distributions in a cavitating venturi using a gamma ray densitometer. Typical profiles show a void fraction of 2 - 10% in the core region of the flow accounting approximately for 60 - 80o% of the flow, while void fractions up to 90 - 100% were observed near the wallso The severe cavitation conditions of the venturi, thus typically producing average void fractions in the order of 10 - 15% over a length

of about two throat diameters, would be expected to far exceed the cavitation void fraction at the pump inlet presently considered. As already pointed out both the wavespeed and the coefficient of attenuation are relatively weak functions of the void fraction once this exceeds a value of about 02 - 0o3% (see also Figure 4.2, Chapter IV).

C. PROPELLANT PUMP CAVITATION AND SUCTION LINE RESONANT FREQUENCY 1. Background Pump cavitation in high speed propellant pumps of the inducer type has been observed to occur not only in the inducer and impeller but also upstream of the inducer inlet as much as 3-4 pipe diameters into the suction line 9. Under such conditions cavitation affects not only the steady pump performance but also the dynamics of the propellant system comprising tank, suction line, pump and discharge line. Cavitation in the suction line increases its compliance hence lowers its resonant frequency which indirectly also affects the pump gain. The cavitation extent for a given pump speed and flow rate depends on the pump suction pressure. The effect of suction pressure on the resonant frequency of the Titan II, Stage I, propellant systems with and without cavitation is shown in Figures 4 and 5 of Ref. (11). Such information is necessary to predict overall system dynamic performance and must be obtained experimentally since the cavitation extent cannot yet be predicted. As described in Ref. (11) it has been the practice to measure the propellant system dynamics on full scale actual systems using a pulse technique to obtain frequency response data. Detailed information regarding the extent of the cavitation zone upstream from the pump inducer inlet for a given operating condition can be obtained by pressure measurements of standing wave patterns in the suction line or by pulsing the discharge line and measuring the pulse travel time to various positions in the suction line. The latter method gives measured values for local -17

-18 wavespeed from which thle local cavi-tation void fraction can be deduced if the properties Of fluid and pipe are knowno To facilitate a detailed prediction of a propellant system dynamics by analysis using the method of characteristics or the (12) impedance method(2) it is necessary to know pump performance, pp(w9ps)9 and, ideally, the suction line cavitation extent in terms of void fraction versus position. The latter information is not explicitly available from tests, however, it is implicit in the experimental data expressed in terms of system resonant frequency versus pump suction pressure for a given pressure (see Figure 4 and 5 of Refo (11)). A model is hypothesized in order to attempt a more generalized correlation of such data leading eventually to a reliable method of predicting the cavitation extent from a minimum. of test data. 20 Model Relating Suction Line Pressure and Resonant Frequency To relate suction line resonant frequency to the pump suction pressure or the cavitation index, K = (Ps - v)/(pVTIP/ it is necessary to kinrow both the geometric configuration of the cavitation void and how this is influenced by the suction pressure. Ignoring any non-uniform radial distribution in the local cavitation void the three simplest void configurations that may be hypothesized are (Fig. 5) (i) Lumped cavitation bubble at inducer.inlet (ii) Constant cavitation void fraction, a < 1 of extent L. (iii) Variable, decreasing void fraction5 a(x), of extent L, Based on the simple lumped cavitation bubble configuration

-19 -ILcl \,...._i E K-il' O o o O o o o O Od 0 0 000 o 0 o o 0 oO o 00 0 L ooo D 0o 00 oo o o i I — Lc-I a~(x) I 0 m 0~ ( o ~,'' 0 O 0, 0'. n 0 0... I L 0 0 Inducer Pump Inlet Propellant Tank Figure 5. Cavitation Void Distributions Line. in Suction

-20 (i) Thorson derived an expression for the resonant frequency of the suction lineo Using dimensional analysis, assuming that geometrical similarity applies to the cavitation bubble, the ratio9 C1 - Lc/D, of cavitation length to pipe diameter was derived as the characteristic parameter which might correlate uniquely with the cavitation indexo Measuring then the actual suction line resonant frequency for given cavitating conditions and knowing the suction line geometry and properties the effect of cavitation was isolated in a "cavitation compliance number", Q = 1/(2jC /2) = Dc As seen from the geometry (Fig. 5 (i)) 3/2 / a D Vc (11) where Vc is the cavitation bubble volume. The experimental data analyzed by Thorson showed Q to depend on the cavitation index in a similar manner for several different propellant pumps, although no unique correlation appeared to existo Furthermore, no theoretical prediction of a possible relation between the two parameters was forwardedo It may be expected that a cavitation void configuration more in agreement with reality is required to predict cavitation extent and its effect on suction line resonant frequency. While the cavitation bubble model (i) involves only one parameter) Lc, the uniformly distributed model (ii) involves two parameters, a and L o Also it is important to note that cavitation is not confined to the suction line upstream of the pump inleto In fact appreciable cavitation is known to exist in the inducer region even at high suction pressures without causing a noticeable reduction in pump performance. Cavitation in the inducer is also expected to cause a reduction of the suction line resonant frequency.

-21 Examination of liquid propellant inducer type pump data for suction line resonant frequency versus suction pressure indicates consistently low suction line resonant frequencies for relatively high suction pressures. This tends to confirm the effect of inducer cavitation, known to exist even for high suction pressures. For the Titan II, Stage I Lox pump, for example, treating the inducer as a 1-foot pipe length with a uniform cavitation void fraction of about a = 0.0015, results in a reduction of the resonant frequency of the 28.6 ft. long suction line to 720 of its value in the absence of cavitation. Cavitation void fractions of this order of magnitude could well exist at suction pressures in excess of 100 - 200 psia as suggested by the data reproduced in Figure 6. In addition to the cavitation in the inducer, upstream cavitation caused by pressure fluctuations generated by the inducer inlet blades will further contribute to a reduction of the effective wavespeed of the suction line, hence its resonant frequency. Depending on the magnitude of the effective amplitude, Apo, of the pressure fluctuations at the inducer inlet, this phenomenon will not appear except at suction pressures ps below a certain critical value pI, which must depend on pump design and operating conditions. Once upstream cavitation occurs, its extent may be evaluated as suggested by Eq. (2) for Ap(L,) =Apo e Lc = Ps - Pv (12) The model proposed below attempts to relate cavitation extent (hence suction line resonant frequency) to pump suction pressure. The model shown schematically in Figure 7 is based on the following assumptions:

1.0 U) LU Oz UL> cr< LL U ZD ZI LLI J Of C] LL Z 0< I 4%4i. I.81 - -.6-.4 I row I / / / Experiment Correlation with model for and inducer cavitarion upstream 2 — / / OL 0 20 40 60 80 100 120 140 160 SUCTION PRESSURE ps(psia) Figure 6. Titan II, Stage I, Lox Suction Resonant Frequency.

f/fo 1 I REGION OF INDUCER AND UPSTREAM CAVITATION I REGION OF INDUCER CAVITATION ONLY C L I - a 0~~~~~~~~~~~~~~~~~~~~~~~~~~~I 0 PS APo Lc X X Figure 7o Model Including Inducer Cavitation (schematic).

-24 (1) For suction pressures Ps greater than a critical value, p, where s Ps Pv +A (13) s v 0 cavitation occurs, only in. the inducer, equivalent to a pipe length Li. The average void fraction over this equivalent pipe length depends on suction and vapor pressures according ro Equation (10)o (2) For suction pressures Ps less than the critical value9 p' additional upstream cavitation occurs to a length Lc determined s' by Equation (12), Ps - Pv Lc =-.in (14) Ps - Pv (3) The average void fraction in the upstream cavitation region depends in general on both pressure level, Ps - Pvy and the amplitude, Ap, of the pressure fluctuations causing the cavitation, c = f(ps - pv, Ap) (15) 3. Example of Data Correlation To illustrate the nature of the proposed model, the Titan II, Stage 1, Lox data was examinedo For simplicity Equations (10) and (15) were assumed to be identical. Then the cavitation extent is described by (a, L) where _= c — (16) Ps - Pv Lf (P < Ps ) (17) L1 Ps - Pv Li + - n — p, >p) (18) P - Pv

-25 Using this model together with the results of resonant frequency calculations based on the impedance method, as shown in Figure 8, leads to the dashed curve in Figure 6 showing agreement with the data from the Titan II, Stage I, Lox pump suction line. The cavitation extent was determined from Equations (16), (17) and (18) taking c = 0.20 psi, -1 Pv = Psat(T = 800K) = 4.36 psia, O, Li = 1 ft., P = 0.56 ft, p' = 92 psia. In view of Equation (13) the physical significance of this value of pT s is that the hypothesized pressure fluctuations generated at the inducer inlet have an amplitude in the order of Apo \ 92 psi. This value and the magnitude of P is in agreement with the general parameter ranges discussed above. It should be stressed that the numerical values of the parameters in this example were obtained by correlation with the data and thus do not have a theoretical basis. Conceptually, however, the model presented appears to be in agreement with the physics of the problem.

10.0 8.0 C.) -J I -J z 0 LI I ui C) Q: OD 6.0 4.0 2.0 1.0.8.6.4 I I I3' I, I I II I I I I I I i.3.4.5.6.7.8.9 1.0 f/fo RATIO OF RESONANT FREQUENCIES WITH AND WITHOUT CAVITATION Figure 8. Titan II, Stage I, Lox Suction Line Resonant Frequency versus Cavitating Length Lc and Average Void Fraction a.

REFERENCES 1o Stripling, L, B,, Acosta, A. J., "Cavitation in Turbopumps Part 1 and 2," Transo ASME, Series Do, p, 326 (1962). 2, Silberman, Eo, "Sound Velocity and Attenuation in Bubbly Mixtures Measured in Standing Wave Tubes," J. Acoust, Soc. of Amer., 29, p. 925, (1957). 3. Holl, J. W., Wislicenus, G. F., "Scale Effects on Cavitation," Trans. ASME, Series D., 83, ppo 385-398 (1961). 4. Cushing, V., Bowden, Go. and Reilly, D., "Three Dimensional Analysis of Bulk Cavitation," Interim Report for ONR, Contract NONR-3709(00) EPCO Project No. 1o06 Sept 24, 1962. 5. Noltingk, B. Eo, Neppiras, Eo A., "Cavitation Produced by Ultrasonics," Proc. Phyo Soc., (London), B63, p. 674-685, (1950). 6. Lamb, H., Hydrodynamics, Dover, 6th edition (1945), p. 122. 7o Hickling, R., Plesset, Mo So, "Collapse and Repound of a Spherical Bubble in Water," The Phys. of Fluids, 7, p. 7 (1964). 80 Olson, H. G,, "High Speed Photographic Studies of Ultrasonically Induced Cavitation and Detailed Examination of Damage to Selected Materials,' Ph D Thesis, Univ. of Mich., College of Engineering (Augo 1966), 9o Smith, W., Atkinson, Go L., Hammitt, F. G., "Void Fraction Measurements in a Cavitating Venturi," The University of Michigan, Industry Program Report IP-581 (Sept. 1962). 10o Fashbaugh, R, H., Streeter, Vo Lo, "Resonance in Liquid Rocket Engine Systems," Jo Basic Eng., 87, (Dec. 1965), p. 1011. 11o Fashbaugh, Ro Ho, "Liquid Propellant Missile Longitudinal Oscillations," University of Michigan (Nov. 1966). 12, Streeter, V. L., "Computer Solution of Surge Problems," Proc. Inst. Mecho Engo, 180, Part 3E (1965-66). 13. Thorson, Mo H., "Analysis to Determine Correlation of Propellant Suction Line Resonant Frequencies with Turbopump Cavitation Compliance," Feb. 11, 1964, (Private communication from NASA)o -27

UNIVERSITY OF MICHIGAN 3 9015 03529 7889