APPLICATIONS OF SECOND-ORDER LINEAR DIFFERENTIAL EQUATIONS TO MODEL A HYDRODYNAMIC RAM CAVITY
By: Adam Nesmith, Andrew Lingenfelter, Joshuah Hess, and David Liu
Characterizing cavity dynamics after a hydrodynamic ram (HRAM) event is a vital part of aircraft safety and survivability. In a combat situation, an aircraft’s fuel tank is susceptible to high-energy impacts from bullets, missile fragments, or other projectiles . These impacts cause pressure fluctuations in the fuel tank, which can cause structural and equipment damage [2, 3]. Additionally, when a projectile hits an aircraft’s fuel tank, a cavity is created. This cavity is composed of fluid, fluid vapor, and ambient gases; and it creates a region of low pressure in the aircraft’s tank [4, 5]. As the cavity collapses, fluid is then forced out of the orifice generated by the HRAM event . Accordingly, development of a method to characterize cavity dynamics to accurately predict cavity collapse is needed. A model predicting the cavity phases during an HRAM event could also provide accurate predictions of fluid ejection.
HRAM events are well documented and researched. Disimile and Toy  found spurt characteristics were influenced by the kinetic energy of the projectile. The transient spray was separated into multiple phases, and it was found that at least two of these phases occur during an HRAM event . Previous research has examined the fluid flow rate through the orifice after an HRAM event . Additionally, recent work has determined the variability of discharge coefficients of an orifice created during an HRAM event . It was also found that cavity geometric features are one of the factors contributing to transient spray . These important properties of the orifice and the cavity have provided greater understanding of the dynamics of the cavity. Furthermore, knowing how cavity phases relate to the spray through the orifice has opened the opportunity to further analyze how cavity dynamics relate to fluid ejecta.
A model is needed to predict the oscillation of the cavity as it expands and collapses after an HRAM event, as shown in Figure 1. Figure 1(a) shows the maximum displacement of the cavity before collapse when the projectile hits the tank. Figure 1(b) displays the displacement of the cavity after this initial collapse, before expanding again in Figure 1(c). When the cavity collapses again, it experiences another local minimum displacement, as shown in Figure 1(d). Each of these maxima and minima is plotted with its respective frame letter in Figure 2. The current research models the repeated expansion and collapse of the cavity as a single degree of freedom spring-mass-damper system.
The damping coefficient and natural frequency of the cavity are estimated from this model. Knowing the damping coefficient, the natural frequency, and the cavity displacement allows for the estimation of the magnitude of the system’s forcing function for the cavity’s first expansion and collapse. Additionally, the energy dissipation is calculated from the damping coefficient, natural frequency, and the cavity displacement. Integrating the energy dissipation numerically using the trapezoidal method gives the energy of the cavity’s radial expansion. The percentage of projectile kinetic energy that causes the radial expansion of the cavity is then obtained using the magnitude of the forcing function and the rate of energy dissipation of the cavity. This result gives greater understanding of where the energy of the projectile dissipates during an HRAM event, thereby allowing researchers to better understand HRAM events. Additionally, combining this spring-mass-damper model with another mathematical model, such as the Rayleigh-Plesset equation, could provide researchers a better understanding of cavity dynamics.
Experimental Setup and Measurements
For the HRAM test, a projectile accelerator fired a spherical ball bearing into a tank of water. The experiment setup, methodology, and equipment are detailed in previous work . This current research used the data collected from the 1,200, 1,495, and 1,800-m/s projectile velocity tests. These differing velocities were selected to measure the effect of the projectile’s kinetic energy on the cavity’s geometry . Displacement data of the cavity were acquired using a Phantom high-speed camera at 28,000 frames/second.
The cavity’s maximum and minimum displacements were measured after the data were collected. Sample frames used for Test 3 calculations are shown in Figure 1. The images selected for measuring displacement were imported into numerical analysis software to analyze the number of pixels the cavity radially expanded. These displacement measurements were made at constant time steps. The displacement was then converted from number of pixels to centimeters using a conversion factor. The time stamps where each of these displacements occurred were calculated using the camera’s frame rate. These measurements were made for each of the four different projectile velocity tests. Using this displacement data, the spring-mass-damper model was applied to the cavity oscillation.
Figure 2 shows an example plot of the cavity displacement for Test 3. Examining the plot, it is apparent the cavity displacement curve does not fully exhibit the oscillation expected for a spring-mass-damper system. The local minima of the plot, represented by points b and d, are sharp peaks. This points to the existence of a restoring force, which is not represented by the spring-mass-damper model. To correct this and still use the spring-mass damper model, a few changes to the model were made. The first was using the half-natural frequency for the model rather than the full natural frequency. This was accomplished by using the period from projectile impact to the first cavity collapse. Additionally, since the spring-mass-damper model cannot capture the sharp peaks caused by the restoring force, each peak is modeled separately. For this current research, only the first peak was modeled using the spring-mass-damper model subject to the initial conditions of the projectile impact.
To model the cavity as a spring-massdamper system, a number of simplifying assumptions were made. The first assumption was the cavity model only accounts for one degree of freedom. The cavity’s expansion and collapse thus do not predict the horizontal elongation of the cavity. The second assumption was the cavity’s radial expansion was symmetrical about the center line where the projectile entered. The third assumption was the mass of displaced water was equal to the cavity’s volume multiplied by the density of water. Another assumption was the mass of water displaced remained constant during cavity oscillation. The final assumption was the tank acted as the spring in the vibrating system and the water in the tank provides the viscous damping. Therefore, this model takes a holistic approach to HRAM, capturing tank stiffness, water mass, and ullage.
The governing second-order differential equation for a single degree of freedom, spring-mass-damper system subject to a transient force is:
where x is the cavity displacement, (∙) is the time derivative, ζ is the damping coefficient, is the natural frequency of the system, F(t) is the transient forcing function, and m is the mass of the system . As previously stated, this current research uses the half-natural frequency for ωn because of the sharp peaks in the displacement plot caused by restoring forces acting on the cavity. Solving equation 1 to find the general displacement function, x(t), when subjected to the forcing function, F(t), results in:
where is a dummy variable of integration, t is time, x0 is the initial displacement, v0 is the initial velocity, and ψ is the phase angle shift . The phase angle shift is given by :
To find the half-natural frequency, ωn , the time from projectile impact to the first local minimum was recorded to find sampling period T. This period was converted to ωn using the equation :
Because the period used was from the projectile impact to first cavity collapse, ωn is the half-natural frequency of the cavity oscillation. The first method used to find the damping coefficient, ζ , utilized the decrement of the displacement peaks given by:
where δ is the logarithmic decrement . The term δ was calculated from:
where x1 is the first peak displacement, xn+1 is the displacement at a subsequent local maximum, and n is an integer representing a displacement peak . The accuracy of ζ typically increases by using displacements separated by a period . Thus, δ was calculated using the first peak displacement and the displacement two periods later.
To validate the accuracy of the logarithmic decrement method in determining ζ for this model, a second method was used. A graphical method estimated ζ by plotting a displacement vs. time curve and performing an exponential curve fit on the peak displacements. This provides an internal check of the logarithmic decrement method to find how closely it represents the decay of the peak displacements. Then, by calculating the energy dissipation using both values of ζ, the degree to which the system damping affects the energy results in the model was observable. This in turn shows whether one method of calculating the cavity’s damping is preferable to the other.
To use the graphical method, a graph of the cavity’s oscillation was produced from the cavity’s displacement, shown in Figure 2. A least-squares curve fit was performed on the maximum displacement data points to find the exponential decay function. This is represented by the line in Figure 2. The R2 values of the least-squares curve fits for Tests 1–4 were 0.4723, 0.8373, 0.9820, and 0.9557, respectively. The curve fits thus had differing degrees of accuracy. However, the results were still used to calculate ζ to compare with the displacement decrement method. This exponential decay function is of the same form as the exponential term given by equation 2. Since the half-natural frequency, ωn, is known, ζ is calculated by equating the calculated exponential term from the curve fit with – ζωn and solving for ζ . The curve fit and displacement decrement methods for calculating were applied to the four data sets used and their results compared.
Using the half-natural frequency, the stiffness of the system, k, was calculated using
where mw is the mass of the water displaced by the cavity . The mass of the water in the tank was estimated using numerical analysis software. First, the cavity’s shape was assumed to approximately resemble a cone. The radius of the cavity was measured, and the height of the cavity was estimated by assuming the cavity’s sides decreased linearly to a point. From this radius and height, the volume of the cavity was estimated using:
This volume was then converted to a mass using the density of water. This estimated mass of water displaced in the cavity was then used to calculate the stiffness of the system. The system stiffness was calculated for each test and compared.
To find the energy to displace the volume of water, it was necessary to estimate the magnitude of the forcing function supplied by the projectile. It was assumed the forcing function supplied by the projectile only affected the first expansion of the cavity. This assumption is warranted since the projectile detaches from the cavity after the initial expansion and before the first cavity collapse. Thus, a curve fit of the first displacement peak using the known responses of a spring-mass-damper system subject to an impulse was used. The impulse function was chosen because the projectile inputs a large force into the system over a relatively short period of time before it disconnects from the cavity. From this curve fit, an estimate for the magnitude of the forcing function was obtained.
Next, the energy dissipation rate was found using the various system parameters, the estimated forcing function, and the calculated velocity of the cavity. The velocity of the cavity was calculated using a numerical analysis software from the cavity displacement and the time stamps. The equation used for finding the energy dissipation rate was:
where x is the cavity displacement, is the calculated cavity velocity, u is the magnitude of the forcing function, and f(t) is the type of forcing function assumed . The energy dissipation of the cavity was then plotted with respect to time.
After this, numerical integration using the trapezoidal method was performed on the energy dissipation points corresponding with the first displacement peak of the cavity. The result was the energy required to expand the cavity radially. Finally, this energy was compared with the kinetic energy of the projectile at impact, resulting in the percentage of projectile kinetic energy used in the cavity’s radial displacement. This process was performed on each of the four tests with their differing velocities using the two different calculated ζ values.
DATA AND CALCULATIONS
The resulting water mass, system stiffness, half-natural frequencies, and damping ratios for the four tests are shown in Table 1. The calculated mass of water is similar for each test, with more water displaced the faster the projectile traveled. This result is expected because the greater projectile velocity causes greater radial cavity displacement. The half-natural frequency of the cavity is between 234.57 and 303.33 rad/s for the four tests. From the tests, the half-natural frequency decreases as the projectile velocity increases. The stiffness of the system was calculated using the mass of the displaced water and the half-natural frequency of the system. The stiffness calculated for each test was between 322.5 and 417.16 kN/m. The stiffness of the system decreases as the projectile velocity increases. This result is because the half-natural frequency is squared in the calculation, as previously shown in equation 7; thus higher half-natural frequencies cause higher stiffness calculations.
The system stiffness results confirm the optical observations of the test data. In the tests with faster projectile velocity, the cavity exhibited greater displacement. This larger cavity means more mass was displaced after the HRAM event. The stiffness of the system, k, is fixed for each test, and the system stiffness is related to the displaced mass and the square of the natural frequency, as shown in equation 7. Thus, if displaced mass increases for a test with faster projectile velocity, the natural frequency of the cavity should decrease. The results summarized in Table 1 demonstrate agreement with this theoretical model.
Examining the results further, the damping of the cavity oscillation for each test is between 0.0289 and 0.0565. This means the decrease in cavity displacement after two cycles is small. Each damping coefficient is significantly less than 1, which is expected for an underdamped spring mass-damper system. However, the damping coefficient significantly less than 1 increases the likelihood of error in ζ, since the difficulty of measuring cavity displacement increases . Lightly damped systems also exhibit increased oscillation magnitudes at ωn . The damping coefficients yielded from the log decrement and curve fit methods both agree closely. Test 1 exhibited the greatest difference in damping coefficients between methods, with a difference of almost 9%.
Using the results from Table 1 and modeling the projectile impact as an impulse, the results in Table 2 were generated. Figures 3(a), 3(b), 4(a), and 4(b) display the cavity oscillation plotted alongside the impulse response estimation for the first peak. For each test, the percentage of the projectile’s kinetic energy causing the cavity’s radial displacement is displayed for the different damping coefficients. The cavity’s radial expansion used between 0.66% and 3.1% of the projectile’s energy among the tests. Considering the shape of the cavity and the speed of the projectile in each test, this is a significant amount of energy used merely to expand the cavity radially at its widest section. The fourth column in Table 2 displays the estimated energy (in joules) that causes the cavity expansion. The accuracy of this number is limited because it was calculated using numerical integration of a small number of data points. To improve the accuracy of this energy calculation, more data points around the time of projectile impact are needed. However, it is significant to note the kinetic energy calculations for each test are almost the same for both damping coefficient calculations.
Examining the results of modeling the first displacement peak as an impulse response show different levels of success. For Tests 1 and 4, the theoretical model closely predicts the behavior of the cavity, as shown in Figures 3(a) and 4(b). On the other hand, Tests 2 and 3 exhibit considerable disagreement between the theoretical first collapse of the cavity and the data measured during the first cavity collapse. This disagreement is shown in Figures 3(b) and 4(a). A closer examination of the impulse response theoretical model is needed to analyze possible causes of these differing results.
The impulse response of a spring-mass-damper system after projectile impact assuming x0 and v0 equal to 0 is :
Matching this theoretical model to the different test data required two major aspects. The first was changing the magnitude, u, to create a maximum displacement equal to the measured data. The second was changing the frequency of the sin (ωdt) term of the equation so the maximum theoretical displacement occurred at the location of the maximum measured displacement. This was needed since the measured cavity oscillation did not complete a full cycle before re-expanding into the next peak. This was due to a restoring force acting on the cavity before the cavity oscillation could complete a full cycle. Thus, to properly curve fit the theoretical impulse response of a spring-mass-damper system to the measured data, these two parameters were estimated to produce a shape that matched the measured data.
Examining the different plots, it appears the variation between the theoretical and the measured displacement curves is primarily due to the rate of cavity collapse. In the case of Tests 2 and 3, the cavity collapsed at a slower rate than predicted by the theoretical model. Conversely, Tests 1 and 4 have a more symmetrical first peak around the y-axis. This means the rate of cavity expansion approximately equals the rate of cavity collapse. Therefore, the theoretical model, which assumes a symmetrical sine wave, predicts cavity behavior better for the tests that have symmetrical measured displacement peak.
In summary, there are two main factors contributing to the discrepancy between the theoretical displacement model and the experimental displacement model. The first factor is the requirement of two parameters, the forcing function magnitude and the frequency shift, to match the theoretical model to the experimental data. Because there was no technique in this current research to find these parameters, they were estimated. If a technique to measure or calculate both parameters was developed, the accuracy of the theoretical model in predicting cavity behavior would increase. The second factor is how closely the rate of cavity expansion matches the rate of cavity collapse. As long as these two rates are similar, the displacement peak is symmetrical, and the theoretical model can predict cavity dynamics well. Further work is needed to develop this model to better account for differing rates of cavity expansion and collapse.
To accurately model the system, the cavity’s true radial displacement and true half-natural frequency are needed. Due to the method used in measuring cavity displacement, there are many sources of error that affect the numbers shown in Tables 1 and 2. The variations in the calculations come from two main sources of error. The first source of error comes in selecting the frame exhibiting maximum displacement from the optical measurements. Choosing a frame introduces error because it is difficult to determine the exact frame where maximum or minimum displacement occurs. This error in determining the maximum or minimum also changes the period, T, thus affecting the half-natural frequency of the cavity. Because the half-natural frequency is used in subsequent calculations, this error is propagated through to the final kinetic energy calculations.
The second significant source of error comes from measuring the radial expansion of the cavity. The cavity becomes increasingly asymmetric as it oscillates. This increasing asymmetry makes measuring the radial displacement of the cavity in numerical analysis software increasingly difficult and less accurate. These difficulties in measuring radial displacement introduce error in the damping coefficient calculations for both the logarithmic decrement and graphical methods. Another source of measurement error is found in assuming the cavity is the shape of a cone in order to estimate its volume. The discrepancy between the actual volume of the cavity and the estimated volume assuming a cone reduces the accuracy of the displaced water mass calculations. Because the water mass is used throughout the kinetic energy calculations, the error in estimating volume propagates through the energy calculations.
Quantifying these sources of error proves extremely difficult. To properly quantify the error, a program in the numerical analysis software would have to measure the displacement in pixels of the cavity images. Implementing this in this current research was not possible due to the poor lighting quality of the cavity images. However, in future work, better instrumentation and data collection techniques could make quantifying the sources of error inherent to this model possible. In the case of this present research, the goal is to show the viability of modeling the cavity formed after an HRAM event as a spring-mass-damper model. Although the error propagated through the calculations potentially causes inaccuracy in the magnitudes of the results, the mathematical model has the potential to yield accurate results provided better data analysis techniques are implemented.
When analyzing the results of this model, it is important to note the limitations of this model. One major limiting factor is the method of measuring the cavity’s displacement. Due to the lighting of the high-speed images, displacement measurements were not taken using a numerical analysis software. Thus, the displacement measurements were made on each cavity image using the human eye. Additionally, the measurements were only taken at the point of maximum radial displacement. This limits the model’s ability to account for the entire shape and complex displacement of the cavity as it oscillates.
The radial displacement was also only measured from the side in one direction from the middle of the cavity. It is because of the limits of the displacement measurements in this model that the projectile’s forcing function was assumed as an impulse response. A more sophisticated technique in measuring the cavity geometry as it oscillates could yield better displacement results that capture the cavity holistically. This might result in assuming a different-shaped forcing function to model the projectile, such as a unit step or sawtooth pulse. All of these factors limit the conclusions drawn from the model’s results.
However, there are a number of general conclusions from the results of this model. First and foremost is the viability of modeling a HRAM cavity as a vibrating spring-mass-damper system. With the implementation of more thorough measuring techniques, a spring-mass-damper model could offer a straightforward mathematical approach to cavity dynamics modeling. Increasing the model from a single degree of freedom to multiple degree of freedom potentially would increase the accuracy of the spring-mass-damper model. A second conclusion is that both the log decrement and the curve fit methods for calculating the damping coefficient yield almost the same energy results. Thus, only one method is sufficient for energy dissipation calculations. Finally, it was found that the impulse response theoretical model can predict cavity displacement assuming a symmetrical displacement peak. Using more accurate measurement techniques and the spring-mass-damper model could yield accurate results for cavity oscillation parameters such as the natural frequency and the damping coefficient. Additionally, this model can estimate the percentage of projectile kinetic energy that causes radial cavity displacement.
CONCLUSION AND FUTURE WORK
This research introduces a spring-mass-damper model for predicting HRAM cavity oscillation. By calculating key parameters and assuming an impulse function, a theoretical model for predicting cavity displacement and projectile energy dissipation for the first cavity expansion was developed. Using this model, an estimate of how much projectile kinetic energy results in radial cavity expansion was obtained. Understanding where the energy of the projectile goes during an HRAM allows for greater understanding of cavity dynamics. The viability of modeling cavity dynamics as a spring-mass-damper system was shown; however, there is much remaining theoretical and experimental research needed to continue developing this model.
There are numerous avenues form further research from this current model. From a theoretical modeling point of view, development of the impulse response to better model asymmetric cavity displacement peaks would improve the accuracy of the model. Research into potential causes of differing cavity expansion and collapse rates would help with this. Furthermore, development of a way to find the parameters used in the forcing function would increase the accuracy of the model. From an experimental point of view, further tests can more accurately measure the oscillatory period and cavity displacement. This would reduce the error in the model, thus improving its accuracy and results. Additionally, synchronized camera footage of the cavity at multiple angle would better account for the asymmetry of the cavity as it expands and collapses.
The model could also further develop to account for multiple degrees of freedom, which would better capture the cavity’s dynamics when it splits into multiple sections. Better lighting during image acquisition would allow numerical analysis software to measure the cavity displacement. This would improve the accuracy of the results measuring displacement and oscillatory period as well as make the error in this model quantifiable. Expanding the model to account for the different fluid densities in the cavity as well as cavity mixing as it collapses is another route of potential research. Furthermore, applying this model to characterize the second and third peaks is another avenue of future research.
One possible implementation is using a piece-wise formulation to model each peak as an impulse response subject to different initial conditions. Once the cavity expansion and collapse are modeled, a possible avenue of future research is relating the cavity radius dynamics to the internal pressure of the cavity. This could possibly use the Rayleigh-Plesset equation to further model the cavity dynamics, particularly how the cavity pressure affects the second and third cavity displacement peaks.
ABOUT THE AUTHORS
Mr. Adam Nesmith is an undergraduate research assistant working at the Air Force Institute of Technology (AFIT). His research has primarily focused on vibrational modeling of dynamic systems.
He holds a B.S. in mechanical engineering from Cedarville University.
Maj. Andrew Lingenfelter is an assistant professor in the Department of Aeronautics and Astronautics at AFIT. His research focuses on additive manufacturing, weapons, and aircraft survivability. He has more than 10 years of weapons and aircraft survivability experience in various program office, test, research, and teaching positions. He earned a B.S. in mechanical engineering from the University of Nebraska-Lincoln, an M.Eng. in industrial and systems engineering from the University of Florida, and a Ph.D. in aeronautical engineering from AFIT.
Maj. Joshuah Hess is an assistant professor in the Department of Aeronautics and Astronautics at AFIT. His research interests include orbital mechanics, optimal control, spacecraft attitude determination, relative satellite motion, and estimation theory. He has also investigated adaptive estimation of nonlinear spacecraft attitude dynamics, as well as the relative navigation between satellites conducting proximity operations. Maj. Hess previously worked as a space systems engineer at the National Air and Space Intelligence Center, and he was deployed in support of Operation Enduring Freedom. He holds a B.S. in aerospace engineering from Virginia Polytechnic and State University, as well as an M.S. in astronautical engineering and a Ph.D. in aerospace
engineering from AFIT.
Maj. David Liu currently serves as the Modeling and Simulation Lead at the Special Programs Office of the Air Force Life Cycle Management Center. Previously, he was a research engineer at the Air Force Research Laboratory, before earning his Ph.D. in aerospace engineering at AFIT. Major Liu then served as an assistant professor of aerospace engineering, where he taught graduate-level courses in aircraft combat survivability, rocket propulsion, weaponing, and additive manufacturing. Additionally, he was deployed with the 82nd Combat Aviation Brigade as the Joint Combat Assessment Team member responsible for Regional Command – East and Regional Command – North.
- Ball, R. E. The Fundamentals of Aircraft Combat Survivability Analysis and Design. American Institute of Aeronautics and Astronautics, 2003.
- O’Connell, P. J. “Dry Bay Fire Test and Evaluation for Representative Large Commercial Aircraft.” Technical Report, 46th Test Group, 2012.
- Disimile, P. J., J. M. Davis, and J. M. Pyles. “Qualitative Assessment of a Transient Spray Caused by a Hydrodynamic Ram Event.” Journal of Flow Visualization and Image Processing, vol. 14, no. 3, pp. 287–303, 2007.
- Lingenfelter, A. J., and D. Liu. “Composition Characterization of Cavity Consisting of Multiple Fluids.” 57th AIAA/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, p. 1245, 2016.
- Lingenfelter, A. J., and D. Liu. “Characterization of Hydrodynamic Ram Cavity Dynamics to Transient Spray.” 57th AIAA/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, American Institute of Aeronautics and Astronautics, 2016.
- Disimile, P. J., and N. Toy. “Liquid Spurt Caused by Hydrodynamic Ram.” International Journal of Impact Engineering, vol. 75, pp. 65–74, 2015.
- Lingenfelter, A. J., D. Liu, and M. F. Reeder. “Time Resolved Flow Field Measurements of Orifice Entrainment During a Hydrodynamic Ram Event.” Journal of Visualization, pp. 1–12, 2016.
- Tatman, L., A. J. Lingenfelter, and D. Liu. “Flow Properties Through a Ballistically Generated Orifice During a Hydrodynamic Ram Event.” 2018 AIAA/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, 2018.
- Lingenfelter, A. J., D. Liu, M. Reeder, and E. Brickson. “Synchronized Imagery Assessment of Hydrodynamic Ram Cavity Features to Transient Spray.” Journal of Flow Visualization and Image Processing, vol. 23, no. 3-4, 2016.
- Meirovitch, L. Fundamentals of Vibrations. Waveland Press, 2010.
- Aström, K. J., and R. M. Murray. Feedback Systems: An Introduction for Scientists and Engineers. Princeton University Press, 2010.
Note: This article was previously presented at the American Institute of Aeronautics and Astronautics (AIAA) SciTech Forum in January 2019. Reprinted with permission.