Infrasonic Technology Workshop November 3-7, 2008, Bermuda, U.K.

1. Infrasonic Technology Workshop November 3-7, 2008, Bermuda, U.K. Session 5: Detection, Propagation and Modeling Oral Presentation: Acoustic-Gravity Waves from Meteor Entry as well as from Rockets and Missiles D.O. ReVelle EES-17, Geophysics Group Los Alamos National Laboratory Los Alamos, New Mexico 87545 USA

2. Summary of Presentation • Introduction, Overview and Motivations: • Limitations of the PPK Normal Mode Code • AGW’s from very rapidly moving, impulsive sources: • Atmospheric source model: Dirac Delta function type source • Weak shock wave propagation from bolides and rockets: Direct and indirect arrivals • Lamb-edge wave formation: Dispersion and Airy functions • Acoustic wave dispersion: Dispersion and Bessel functions • Ducted acoustic waves: Sound speed and horizontal wind speed sensitivity using Gaussian beam theory “miss” distances, near integral number of hops, etc. • Total signal construction as a function of range, source size, etc. • Applications: The Carancas meteorite fall/crater: 09/15/2007 • Summary and Conclusions

3. The Pierce-Posey-Kinney (PPK) Normal Mode Code • Motivation- New AGW Code work due to significant PPK Code Limitations: • PPK: Designed for stationary point “bomb” sources: Small source limit (Rs << H) and scaling at relatively low heights • Limitations: Relatively large sources at great ranges/long times • Significant sensitivity to upper boundary conditions at ~150 km (where the atmosphere is most poorly known and with highly non-steady conditions AND with an unmixed diffusive separation of light constituents resident on top of heavier ones) • Leaky modes not allowed (Only discrete modes included) • Exclusively linear, but full wave theory (No changes of behavior are predicted at caustics, etc.) • A code was needed for meteors/bolides (small sources at close range and large sources at great range and intermediate ranges, etc.) as well as for rockets and missile sources.

4. Examples: Bolide entry andSecond stage ignition sequence

5. Specific Source Examples: • On the next slide, we have plotted the predicted line source cylindrical blast wave pressure amplitude time series (at the altitude corresponding to the maximum blast wave relaxation radius conditions using ReVelle’s TPFM entry modeling) for: • The Revelstoke Meteorite Fall of March 31, 1965 (Canada) • The Carancas Meteorite Fall of September 15, 2007 (Peru) • On the next slide, we have also indicated a conceptual view of two possible source mechanisms for the generation of infrasound from very rapidly moving rocket and missile sources • On the next slide we present a conceptual view of the construction of a energetics-based, self-consistent solution for AGW from bolides as well as form rockets and missiles

6. Blast Wave Signatures: Revelstoke and Carancas Meteorite Falls 65 km source altitude assumed: Matching of Peru and Paraguay infrasound signal amplitudes 14.3 km sourcealtitude assumed: Based upon TPFM Modeling

7. Rocket Infrasound:Conceptual View Rocket motion Leading shock front Main rocket body Trailing shock front Supersonic/hypersonic shock front regime:  Mach 10-20 for steady state aerodynamic flight Engine High temperature Mach 3 plume with wake turbulence generated “noise” , i.e., acoustic-gravity waves-AGW Nonlinear atmospheric refraction and heating regime at close range

8. Two-Dimensional Source Modeling Model Atmosphere: Horizontally stratified, range independent, steady state, hydrostatic model atmosphere (including seasonal detailed properties and horizontal atmospheric winds, etc.) Geopotential altitude, z Ray-mode skip distance of the ducted wave paths: Intermediate range: Range ~ O({Ro, R1}) Duct height Internal wave generation region Ducted Strato- and Thermo- spheric arrivals Bolide source Lamb wave formation Observer Lamb edge wave guided arrivals Weak shock/linear wave arrivals +x direction

9. Overview:Wave Solutions Expected • Starting from a very narrow short period blast wave pulse whose attributes depend on the source properties (energy, altitude, etc.): • Internal Gravity Waves are launched; Lamb-edge wave solutions (Lowest order gravity wave mode at low frequency) expressed as an Airy function with longer periods traveling faster and arriving earlier than higher frequencies: Normal dispersion- The final shape is a function of range and source energy (Gill, 1982). • Internal Acoustical Waves are launched: Bessel function solutions with shorter frequencies traveling faster and arriving earlier than lower frequencies: Inverse dispersion (Tolstoy, 1973). • Ducted Acoustical Waves (trapped in a waveguide) whose amplitude and wave period are a function of the source properties and whose behavior is critically dependent upon the mean and perturbed atmospheric sound speed and wind speed structure. • Weak shock waves: Direct/indirect arrivals at sufficiently close range whose properties depend on the properties of the source

10. Lamb wave (LW) Analyses A.D. Pierce, J.A.S.A., 1963: Isothermal, hydrostatic atmosphere () • R0 = {/(B2)}h = Lamb Wave (LW) formation distance • R1 = {(2 2/(23))h2exp[2B2h]}1/3 = LW dominance distance • R2 = {2/(B2 )}h; R3 = {(2)1/2 /}h where • h = zs /Hp = Dimensionless height of the source • B2 = (2 - )/(2)  0.2143 ;  = 1.40= Ratio of specific heats for air • 2 = 2 - A2 ; 2 = 2 – ¼ ; A2 = ( - 1)/ 2  0.2041 • 2  ( + i) = Scaled wave frequency squared including losses •   o/(cs/Hp) = Non-dimensional (scaled) wave frequency • cs2 = gHp = Adiabatic thermodynamic sound speed squared • Hp = Atmospheric pressure scale height (): Lamb waves can still exist below non-isothermal inversion conditions as shown previously by Kulichkov and ReVelle (2002).

11. LW Analysis Procedure • Using the following set of inequalities we can completely describe the AGW solutions as a function of the horizontal range R: • A. R < {R0 , R1}: Weak shock/quasi-linear acoustic waves only • B. R ~ O{R0 , R1}: Weak shocks/linear and Lamb waves present • C. R > {R0 , R1}: Lamb wave solutions dominate the response • {R0 , R1}: Explicit horizontal length scales developed by A.D. Pierce (1963) for the prediction of the presence and dominance of the Lamb wave after a small explosion (Ro(z) << H(z)), where Ro(z) = Blast wave relaxation radius and H(z) = Density scale height of the atmosphere ~ 10 km {H(z) = - //z}.

12. Lamb Formation and DominanceDistances: Source altitude = 20 km Using a Rayleigh friction viscous decay formalism

13. Lamb Wave Dominance Distance versusWave period and Geopotential height

14. Airy function & its large argument behavior Tunguska bolide observation conditions in Great Britain (June, 1908)

15. Dirac Delta function source:Acoustic, Bessel function solution

16. Summary: Guided Internal WavePropagation Characteristics • Three AGW propagation schemes have been employed: • Downwind- Stratospheric and Thermospheric returns • Upwind- Thermospheric returns only (Diffracted Stratospheric returns were neglected due to small expected amplitudes) • Crosswind- Thermospheric returns only expected • Signal velocities computed internally for the US Standard Atmosphere (1976) model: • csig-Stratospheric (downwind) = 0.2945 km/s • csig-Thermospheric (up-wind) = 0.2667 km/s • csig-Cross-wind = 0.268-0.293 km/s (depending on the phase angle) • Number of hops, slant range and travel time (ducted internal waves) were computed/compared to the horizontal range and travel time for the Lamb wave. • Time delay (internal waves) used to control internal wave onset.

17. Ducted Wave Characteristics • The procedure used to determine if these waves could be ducted was treated “exactly” for an idealized waveguide assuming: • Angle of incidence = Angle of reflection at both the upper and lower duct interfacial boundaries (ignores the detailed bottom topography, especially at the shorter wavelengths). • The vertical gradient of the wave amplitude was assumed to be zero so that the wave frequency and amplitude were unchanged upon reflection at both the upper and lower waveguide boundaries. • A near-integer number of hops, n, must exist between the source and observer. Also if the number of hops < some small limiting value (~0), these solutions were also rejected. • A “miss” distance was computed for each of these “rays” that satisfied the above conditions. The miss distance value was further assigned on the basis of the computed e-folding widths of Gaussian beams (Porter and Bucker, JASA, 1987, etc.) as a function of horizontal range.

18. Summary of predicted ducted wave arrivals vs the integral mode number (range = 1620 km) Miss distance = 20 km, Integral no. of hops +/- 0.20

19. AGW Signal Construction • Combined AGW Pressure Amplitude Response: Using separation of variables, we have constructed amplitude predictions as a function of the blast radius, source height and source energy, etc. in the form (neglecting, a ground reflection factor): • pL(x, z, r, t) = psrc(Ro, z)L(x,t)ZL(z)XL(r) = Lamb wave amplitude • piw(x, z, r, t) = psrc(Ro, z)iw(x,t)Ziw(z)Xiw(r) = Internal wave amplitude • pws(x, z, r, t) = psrc(Ro, z)ws(x,t)Zws(z)Xws(r) = Weak shock wave amplitude • ptL(x, z, r, t) = pi = pL(x, z, r, t) + piw(x, z, r, t) = Total amplitude; r > {R0,R1} • ptws(x, z, r, t) = pi = pws(x, z, r, t) + piw(x, z, r, t) = Total amplitude; r < {R0,R1} where psrc(Ro, z) = Blast wave source amplitude at x = 10viaresults in ReVelle (2005). i(x,t) = Wave shape function (normalized between -1 and 1) Zi(z) = Kinetic energy density conservation (inviscidfluid approximation) Xi(r) = Geometrical wave spreading function for a two-dimensional waveguide Xi(r) = {r/Ro}-;  = Geometric spreading decay factor (constant); ½  ¾ {R0, R1} = Distance scale over which the Lamb wave signals develop

20. Carancas Peru Meteorite Fall- 09/17/2007: Basic Observations • Visual eyewitness accounts: Region of very high altitude (see below) and very rugged terrain in northeastern Peru • Ancillary Observations: Local sounds heard, broken bull’s horn, etc. • 13.6 m diameter crater produced at 3826 m elevation (!) • No satellite data available for this event • Infrasonic wave data recorded at two IMS (International Monitoring System) arrays: • In Bolivia (I08BO) with a large signal/noise ratio (S/N) • In Paraguay (I41PY) with a very small S/N. • Seismic wave data (at several IMS stations) such as LPAZ (LaPaz, Bolivia), UBINS, etc. including direct crater impact arrivals.

21. Carancas: Key InputDirect Entry Modeling Parameters • 0.79 R Initial bolide radius (m) [.000001 - 1000.0] • 12.6 V Initial velocity (km/s) [11.2 - 73.0] • 30.0 ZR Angle of entry with respect to the vertical (deg) [0.0 - 80.0] • 16.0 NMax Maximum number of pieces of fragmentation [1 - 1000] • 1.209 Sf Shape factor (area/volume2/3 ) 1.209 = sphere [1.209 - 2.0] • 0.6667  Shape change factor 2/3 = no change [- 3 to - 0.6667] • 4.605 D Kinetic energy at end height [2.303 - 4.605] i.e. [10% - 1%] • 1.0 BRKTST Allow breakup 0 = no; 1 = yes [0 or 1] • 1.0 FRAGTST Fragments in wake 0 = remain; 1 = Stay with body [0 or 1] • 0.0 PORTST Allow porous materials 0 = no-porosity; 1 = porous [0 or 1] • 1.0 POR Porosity or Fireball group [0 to 1 or 0.0 - 5.0 (uniform bodies) • Initial entry mass = 7.641103 kg • Initial entry kinetic energy = 0.1449 kt = 6.06571011 Joules • Line source blast radius: Top of the atmosphere (initial) value= 61.47 m; Maximum radius value = 156.9 m; Minimum radius value = 27.76 m • Predicted maximum wave period (at x = 10) = 1.469 s • AFTAC source energy (observed maximum 0.62 s wave period)  0.533 t

22. Carancas: Direct Entry Modeling Results

23. Inverse Entry Modeling Summary

24. Predicting the NominalHypersonic boom Corridor • Line source wave normal ray tracing procedure through a specific atmospheric model as a function of geopotential altitude (steady state, hydrostatic and range-invariant) • Assume instantaneous energy release for uniform phase • Specify/measure the adiabatic, thermodynamic sound speed and the mean horizontal winds as a function of height • Specify the entry angle and the azimuth heading angle of the bolide • Specify a fixed entry velocity and predict the cone half angle at al heights and times • Predict the characteristic velocity (Snell’s law constant of Geometrical Acoustics) and wave normal paths valid if R < 2H2/ (H = Vertical duct height; = wavelength) • Predict the arrival timing and arrival angles of the waveforms

25. Temperature/sound speedand wind speed profiles- LaPaz, Bolivia

26. BEST Solution: View Overhead andfrom the West:  = 262 and  = 62 Hypersonic Boom Corridor WEST OVERHEAD

27. Theory and Observations: Infrasound at I08BO: Bolivia Best entry time solution = 16:40:10 UT +/- 5s P.Brown- personal comm. Predicted AGW solution: Propagation delay time ~ 287 s I08BO infrasound arrived at ~ 16:44:20 UT Inferred time delay ~ 250s

28. Carancas: Paraguay Observations and Simulation of Arrivals Downwind propagation conditions Thermo Strato Pa Observed arrivals Time with respect to 1800 UT

29. Summary and Conclusions-I • We have developed a very general computer code in order to calculate the properties of acoustic-gravity wave (AGW) signals from bolides or rockets and missiles as a function of: • Range (and scaled range) • Source energy • Source height • Atmospheric temperature & horizontal wind-Vertical structure • The model incorporates four types of atmospheric responses: • Lamb edge-waves • Internal acoustic waves • Ducted acoustic waves (whose properties depend on the source) • Weak shock waves: Direct (and indirect) arrivals at sufficiently close range

30. Summary and Conclusions-II • The Carancas meteorite impacted in the high mountains of Peru on September 15, 2007. We have since determined that: • No satellite observations are available for this event • Inverse entry modeling (bottom-up) produced only a limited range of size, velocity and angles of entry for an impact at 3826 m  500 m above sea-level (height error assigned). This was also found consistent with a crater diameter= 13.6 m  2 m. Our best, self-consistent modeling result assumed a meteorite density= 3300 kg/m3 and a target “soil” density= 2000 kg/m3. • Direct entry modeling (top-down) produced results that are self-consistent with the inverse (bottom-up) modeling result. • The deduced nominal hypersonic boom corridor is also consistent with existing near-field observations of both infrasonic as well as seismic waves. • Long range modeling of acoustic-gravity wave signals has also been found consistent with the details of the source as discerned by other complimentary approaches.