1 Introduction

The history of shock wave theory hearkens back, at least, to Siméon-Denis Poisson’s “A Paper on the Theory of Sound,” published in 1808 (Poisson 1808, 1998; Salas 2007). In this work, Poisson, who himself built on the theories of Joseph-Louis Lagrange, used physical laws to formulate a differential equation describing the vibrations of air molecules and the resulting transmission of sound waves through the air. Poisson demonstrated that for a one-dimensional spatial domain, the solution to his model consists of wave-like profiles of the air velocity v(x) satisfying the implicit expression \(v = f[x - (a+v)t]\), where a is the speed of sound and f is an arbitrary function.

In 1848, George Gabriel Stokes’ “On a Difficulty in the Theory of Sound” picked up on the threads of Poisson’s results (Stokes 1998). Stokes demonstrated that if f is taken as a sine wave, v(x) will develop points of infinite slope in finite time. Moreover, after that time, the solution cannot be written as an algebraic expression of the form familiar in Stokes’ era. Stokes dealt with the problem of infinite slope by proposing a “jump” in the solution, which we now know as a shock. Stokes struggled with the idea of this new type of solution, ultimately disclosing to the reader, “It seems to me to be of the utmost importance, in considering the application of partial differential equations to physical, and even to geometric problems, to contemplate functions apart from all ideas of algebraical expression.” Throughout his life, Stokes vacillated on whether or not he believed his own hypothesis regarding shocks, ultimately including a disclaimer in his collected works (Salas 2007). One of Stokes’ principal concerns was that the mathematical model failed to preserve conservation of energy in the system.

Belief in the existence of shocks found firmer scientific footing with developments in the field of thermodynamics. In 1870, William Rankine proposed “jump conditions” at shocks which did, in fact, enforce conservation of energy. An insight that followed was that the thermodynamic process within the shock must be non-adiabatic. In 1885, Pierre-Henri Hugoniot demonstrated that in order for the system to conserve energy, there must be a jump in entropy across the shock (Hugoniot 1889, 1998). The two aforementioned observations led to the so-called Rankine–Hugoniot (RH) equations, which describe the conservation of mass, momentum, and energy across the shock, and which remain fundamental to understanding the properties of the discontinuity.

Our discussion thus far has drawn heavily from the history outlined in Salas (2007). The works of Poisson, Stokes, and Hugoniot that we have cited are available in a compilation of republished historical papers key to the development of shock wave theory (Johnson and Chéret 1998), which additionally contains relevant works of Samuel Earnshaw, Georg Friedrich Bernhard Riemann, Lord Rayleigh, Geoffrey Ingram Taylor, Hans Bethe, and Hermann Weyl. The valuable work of Clifford Truesdell contextualizes Stokes’ oeuvre and gives a broad view of challenges in the development of shock wave theory and thermodynamics Truesdell (1980, (1984).

Towards the end of the nineteenth century, the mathematical research on shock waves found application to explosions. Simple models of explosive shock waves often mention the eponymous Chapman–Jouget (CJ) state, which is the hypothesized state of the chemical products behind the reaction zone of a detonation. In this zone, products move at the minimum possible velocity greater than the speed of sound in the detonation medium in order to propagate the wave forward without interference from the shock’s rarefaction waves (Dremin 2012; Guo et al. 2016; Chéret 1999). Also known as the “zero-reaction zone” model, CJ theory implicitly supposes that the transformation of explosive materials into products via a shock wave is instantaneous (Dremin 2012).

Proposed by Chapman (1899) and built upon by Jouguet (1905, (1906), the hypothesis coalesced into the first theory of detonations at the start of the twentieth century (Chéret 1999). While the theory is useful in predicting detonation states, it is also limited in its application due to its breakdown in the case of non-ideal gases (Keshavarz and Warey 2007). Despite ample experimental evidence that the law was far from universal, mathematicians and physicists seeking to build generalized models of detonations during this time assumed that CJ theory was not just locally true, as was the original context of its proposal, but useful in determining the state of a flow at any point downstream of the shock (Chéret 1999).

One of the most ambitious detonation modelers was Hungarian-American John von Neumann, a titan of scientific discovery who made numerous and influential contributions to mathematics, physics, and computer science during the first half of the twentieth century. Von Neumann’s published works comprise six volumes on topics ranging from functional analysis to game theory to linear programming. During and after his involvement in the Manhattan Project, von Neumann published ten papers that use analytical and numerical methods to investigate the hydrodynamics of shock waves and detonations (von Neumann 1963a, c, d, e, f, g, h, i, j, k). We summarize these in Table 1.

Table 1 John von Neumann’s previously known research on shock waves and detonations, listed in chronological order
Full size table

The next major leap in detonation theory was driven by necessity during World War II. Between 1940 and 1943, Yakov Borisovich Zel’dovich, John von Neumann and Werner Döring independently contributed to the development of a model of detonation that would come to be known, by their initials, as ZND theory (Dremin 2012). Although ZND theory uses the hypothesized CJ state, it assumes a finite (rather than infinitesimally thin) reaction zone that follows the shock wave. In this zone, materials are compressed before they undergo chemical transformation (Dremin 2012).

In 2017, a private collector shared with us an unpublished, handwritten manuscript of John von Neumann’s (von Neumann unpublished). We present a transcription and analysis of this previously unknown and unstudied work and contextualize it within the history of shock waves and detonations as described above. The transcription is a faithful reproduction of von Neumann’s writing, with exception of several commas that have been excluded for readability.

In the manuscript (von Neumann unpublished), von Neumann studies bombs comprising a primary explosive charge along with explosive booster material. His goal is to calculate the minimal amount of booster needed to create a sustainable detonation, presumably because booster material is often more expensive and more volatile. In service of this goal, he formulates and analyzes a partial differential equation-based model describing a moving shock wave at the interface of detonated and undetonated material. Section 2 provides an overview of the manuscript’s origin and contents. In Sect. 3, we present each section of von Neumann (unpublished) interspersed with our own analyses and explanations, including the correction of two small errors in von Neumann’s work. Finally, in Sect. 4, we conclude by situating (von Neumann unpublished) in the context of other work of the era, and hypothesizing the motivation behind the manuscript and the possible reasons it was never published. We suspect that the document was written during World War II, around the same time as (von Neumann 1963j, k), both of which reference a pending publication relating to boostered detonations that never materialized (or was never declassified).

2 Overview of the manuscript

Von Neumann’s manuscript (von Neumann unpublished) was obtained by a private collector through a rare book dealer, who authenticated the manuscript. The document was part of a larger archive of working manuscripts, reports, notebooks and letters produced by von Neumann and his collaborator Raymond J. Seeger, who originally compiled the archive. The manuscript we investigate, (von Neumann unpublished), is not dated. The top right corner of the first page contains von Neumann’s handwritten initials; see Fig. 1. In the pages following von Neumann’s primary manuscript, there are additional pages entitled “Stable Equilibrium of a Plate or Membrane under an External Pressure f” very likely written by Seeger; see Fig. 2. We did not investigate these pages and restrict our attention to the work confidently attributable to von Neumann. Based on the complexity of operations between steps in von Neumann (unpublished), we suspect this document was a write-up of other work, rather than the original derivation.

Fig. 1
figure 1

Top portion of the first page of von Neumann (unpublished) in von Neumann’s handwriting. Here, he begins by stating several assumptions and defining certain quantities in his detonation model. See Sect. 3.1 for a typed transcription

Full size image
Fig. 2
figure 2

Image of the unrelated section of the manuscript in Seeger’s handwriting. We did not study this portion of the document

Full size image

Von Neumann aims to calculate the ratio of booster material to primary charge necessary to create an operational bomb. In this model, the more powerful (but typically more expensive and more volatile) booster material helps create a successful detonation in the primary charge. Von Neumann assumes the detonation process creates an ideal gas with uniform composition as the shock wave progresses through the solid explosive material. He further assumes that the isentropic reaction reaches thermodynamic equilibrium instantaneously with a non-adiabatic jump discontinuity at the line of the shock wave. Rather than implement a discontinuity in liberated energy between the booster and main charge, von Neumann simplifies by modeling the energy output moving away from the center of the explosive as a monotonically decreasing function. The solution von Neumann reaches is meant to describe a boostered detonation along a line, cylinder, or sphere, corresponding respectively to dimensions \(q=1,2,3\).

Table 2 All variables and parameters in von Neumann (unpublished)
Full size table

Table 2 summarizes the variables and parameters used in von Neumann (unpublished). Figure 3 provides a schematic of the boostered detonations under consideration. In the following section, we present a transcription of all 19 sections of von Neumann’s manuscript (von Neumann unpublished), along with our accompanying analyses. The transcription is von Neumann’s wording verbatim, with the exception of erroneous commas that we have excluded to improve readability.

3 Manuscript transcription and analysis

3.1 Sec. 1 transcription

1. Consider a solid explosive of specific volume \(v_0\), which explodes, liberating energy \(e_0\) per unit mass and transforms it into gas that we assume to be ideal, of adiabatic exponent \(\gamma \).

Let p, v be pressure and specific volume of this gas, D the detonation velocity, V the gas velocity behind the detonation front. We treat the explosive reaction as instantaneous.

Then the Rankine-Hugoniot equations are:

$$\begin{aligned} \sqrt{\frac{p}{v_0-v}}= & {} \frac{D}{v_0} = \frac{D-V}{v} \end{aligned}$$
(1)
$$\begin{aligned} \frac{1}{2}p(v_0-v)= & {} \frac{1}{\gamma -1}pv - e_0. \end{aligned}$$
(2)

(2) is solved by

$$\begin{aligned} v = \frac{\gamma -1}{\gamma +1}v_0 + \frac{2(\gamma -1)}{\gamma +1}e_0p^{-1}. \end{aligned}$$
(3)

Now (1) gives

$$\begin{aligned} D= & {} \sqrt{\frac{\gamma +1}{2}v_0} \ \frac{p}{\sqrt{p-\frac{(\gamma -1)e_0}{v_0}}} \end{aligned}$$
(4)
$$\begin{aligned} V= & {} \sqrt{\frac{2}{\gamma +1}v_0} \ \sqrt{p-\frac{(\gamma -1)e_0}{v_0}} \end{aligned}$$
(5)

Thus p is indeterminate, but it determines v, D, V by (3), (4), (5).

Fig. 3
figure 3

Explosive devices that von Neumann studies in von Neumann (unpublished). a Boostered explosive for a one-dimensional device (\(q=1\)), a cylindrical device (\(q=2\)), and a spherical device (\(q=3\)). b Schematic of the detonation process. In this schematic, we move into the reference frame of the detonation front. In panel (a), the front moves to the right for \(q=1,2\) and radially outward for \(q=3\). In (b); however, the front is stationary by design. In Region I, the charges and, eventually, exterior gas, contact the detonation front at supersonic velocity. In Region II, the detonation front, the material compresses and decelerates to subsonic speed. Region III consists of exploded material, now in a gaseous phase and moving at sonic speed

Full size image

3.2 Sec. 1 discussion

Eq. (1) is derived from the RH equation for conservation of mass and the Rayleigh line for solid explosives, and (2) is derived from the RH equation for conservation of energy, also known as the Hugoniot equation (Le Roy 2000).

In this section, von Neumann’s derivation requires that \(p=p_1-p_0\) in (1) and \(p=p_1+p_0\) in (2), where \(p_1\) is the pressure ahead of the shock wave and \(p_0\) is the pressure behind. Thus, von Neumann seems to make an unstated assumption that \(p_0\) is negligible. In von Neumann (1963g), he explicitly states that the difference between \(p_1\) and \(p_0\) is so great that \(p_0\) is negligible, so he assumes \(p_0\) is zero. However, the model in von Neumann (1963g) assumes a gaseous point source, such that “as the original high pressure sphere shrinks to a point, the original pressure will have to rise to infinity.” Von Neumann does not make this assumption in any of his other published works (von Neumann 1963b), and this assumption is not made in detonation models today.

Von Neumann follows the same initial steps in von Neumann (unpublished) as he does in von Neumann (1963g), solving for specific volume, detonation velocity, and gas velocity. Eqs. (3),(4), and (5) in von Neumann (unpublished) are analogous to (2.9\(^\prime \)), (2.10\(^\prime \)) and (2.11\(^\prime \)) in von Neumann (1963g) (with the latter expressed here using the notation of von Neumann (unpublished)):

Expressions from von Neumann (unpublished) Expressions from von Neumann (1963g)

These expressions differ because in von Neumann (unpublished), von Neumann includes the term \(e_0\) in (2). He eliminates this term in von Neumann (1963g) because \(e_0 = p_0v_0/(\gamma - 1)\) and \(p_0\) is negligible. Perhaps \(e_0\) is included in von Neumann (unpublished) because it describes the energy of a solid, rather than a gas as in von Neumann (1963g), so gas laws do not apply.

Following this calculation, von Neumann (unpublished) and von Neumann (1963g) diverge. The latter relies on an energetic calculation that takes into account thermal and kinetic energy, while the former uses a power law to describe energetic output; see Sect. 3.5.

3.3 Sec. 2–3 transcription

2. According to Chapman and Jouguet, the stationary detonation is characterized by that value of p, which renders D a minimum. By (4) this occurs when

$$\begin{aligned} p = \frac{2(\gamma -1)e_0}{v_0}. \end{aligned}$$

Affixing an index s to all quantities referring to this state, we now have:

$$\begin{aligned} p_s= & {} \frac{2(\gamma -1)e_0}{v_0}, \end{aligned}$$
(6)
$$\begin{aligned} v_s= & {} \frac{\gamma }{\gamma +1}v_0, \end{aligned}$$
(7)
$$\begin{aligned} D_s= & {} \sqrt{2(\gamma ^2-1)e_0}, \end{aligned}$$
(8)
$$\begin{aligned} V_s= & {} \sqrt{\frac{2(\gamma -1)e_0}{\gamma +1}}. \end{aligned}$$
(9)

3. In what follows we shall study detonation processes, in which p and D are higher than the Chapman-Jouguet (stationary) value. We put therefore

$$\begin{aligned} p = w \frac{2(\gamma -1)e_0}{v_0} = w p_s, \quad w > 1. \end{aligned}$$
(10)

Then (3), (4), (5) give

$$\begin{aligned} v= & {} \frac{\gamma +w^{-1}-1}{\gamma +1}v_0 = \frac{\gamma +w^{-1}-1}{\gamma }v_s, \end{aligned}$$
(11)
$$\begin{aligned} D= & {} \frac{w}{\sqrt{2w-1}}\sqrt{2(\gamma ^2-1)e_0} = \frac{w}{\sqrt{2w-1}}D_s, \end{aligned}$$
(12)
$$\begin{aligned} V= & {} \sqrt{2w-1}\sqrt{\frac{2(\gamma -1)e_0}{\gamma +1}} = \sqrt{2w-1}V_s. \end{aligned}$$
(13)

3.4 Sec. 2–3 discussion

The RH Eqs. (1) and (2) are underdetermined. Von Neumann uses the CJ hypothesis to find the value of pressure \(p_s\) corresponding to the minimum detonation velocity \(D_s\) necessary to prevent deflagration (collapse) of the detonation. The solutions for \(p_s\) and \(D_s\) appear in (6) and (8). Given these solutions, he further solves for the values of \(v_s\) and \(V_s\), given in (7) and (9). Finally, von Neumann introduces \(w>1\), defined by (10). Here, w is a convenience of notation meant to convey that the wave propagates faster than the CJ value. Eqs. (10–13) restate the values of \(p_s\), \(D_s\), \(v_s\), and \(V_s\) using this new notation.

3.5 Sec. 4–7 transcription

4. Consider now an explosive of variable composition. Using the point at which the detonation begins as origin, and denoting the distance from it by x, we assume that the power of the explosive decreases as x increases.

The purpose of this model is to describe the mechanism of a booster. Indeed, a charge with booster consists of two zones: The inner one, \(x<x_0\), being occupied by a more powerful explosive than the outer one, \(x>x_0\). (These two explosives are the “booster” and the “main charge.”) We replace these two zones by a continuous decrease, mainly for the sake of mathematical convenience: It will appear, that for a certain power-law of decrease the differential equations of the problem are amenable to numerical treatment.

In order to avoid complications which are unconnected with our main problem, we assume that the density of the explosive is the same everywhere. We assume further that the gas resulting from the explosion is always ideal and has the same adiabatic exponent \(\gamma \) everywhere.

Thus the only quantity which can be used to express the power of the explosive is the energy of the explosion \(e_0\). Therefore \(e_0\) must be a monotone decreasing function of x. As mentioned above, we assume a power-law:

$$\begin{aligned} e_0 = e_0(x) = ax^{-2n}, \quad n > 0. \end{aligned}$$
(14)

5. We can treat this problem in \(q = 1, 2, 3\) dimensions, corresponding to plane, cylindrical, spherical detonation waves, respectively. Our considerations are best carried out without specializing this q.

6. When an explosive of this nature detonates, it is to be expected that the detonation will be more powerful at any particular place x, than would correspond to a homogeneous explosive, having everywhere the composition that exists at x. This must be so, since the explosive at x is backed, and hence “boosted”, by the more powerful explosive at \(x^\prime <x\). Hence we assume

$$\begin{aligned} p = wp_s, \quad w > 1. \end{aligned}$$
(15)

Owing to the similitude within the entire arrangement, we may expect that

$$\begin{aligned} w \text { is independent of }x. \end{aligned}$$
(16)

Accordingly we make this assumption.

7. In the subsequent paragraphs, we will determine the stationary detonation process which is compatible with these assumptions. It will appear that for a suitable concept of stationarity, precisely one such process exists, if \(v_0\), \(\gamma \), a, and n are given. This determines w.

More specifically: We shall obtain a functional relation between n and w, which can be written as

$$\begin{aligned} w = w(n) \text { or }n = n(w). \end{aligned}$$
(17)

3.6 Sec. 4–7 discussion

Von Neumann introduces the physical description of a booster, explaining that it composes the inner layer of a two-layer explosive. For “mathematical convenience,” von Neumann describes the energy difference between the booster and main charge of the explosive as a monotonically decreasing function and assumes uniform density of the explosive rather than incorporating a jump discontinuity in the velocity of evolved material as the shock wave progresses between the two materials. Perhaps von Neumann never returned to this subject in his published works because this simplification compromised the usefulness of the model.

He proposes that the model describes an explosive of any dimension given that the system is constant across all angular dimensions. The generalization to several dimensions in this subsection poses some issues for the model based on previous assumptions, particularly the CJ hypothesis. Later expansions on the CJ hypothesis by Taylor and Zel’dovich showed that the hypothesis was not applicable to spheres due to an infinite gradient in the solution behind the detonation front (Bach and Lee 1972).

Here, von Neumann also expands on the difference between a boostered and standard detonation, and he explicitly states his goal, namely, to “determine the stationary detonation process which is compatible with these assumptions.” Restated, von Neumann hopes to use the more easily experimentally determined quantities \(v_o\) (specific volume of solid), \(\gamma \) (adiabatic exponent of gas), a, and n to derive w (see Table 2), which in turn would allow him to compute the ratio of booster to charge necessary for a successful explosive, as he does later in the manuscript. Through the rest of the manuscript, von Neumann solves a Lagrangian partial differential equation for w.

3.7 Sec. 8–9 transcription

8. Before we undertake these computations, let us consider the interpretation of our procedure somewhat more closely.

Consider the detonation at the moment when it has progressed as far as x. The total energy liberated up to that moment is (\(S_q\) is the area of the surface of the unit sphere in q-dimensional space)

$$\begin{aligned}&\int _0^{x}{e(x') \cdot S_qx'^{q-1}}\mathrm{d}x' \\&\quad = S_qa\int _o^{x}{x'^{q-2n-1}}\mathrm{d}x'. \end{aligned}$$

In order that this expression be finite, we must require \(q-2n>0\), \(n<\frac{q}{2}\). In order to have the outward decrease of e(x), mentioned in 4., we must require \(n>0\). Hence

$$\begin{aligned} 0<n<\frac{q}{2}. \end{aligned}$$
(18)

Now the above integral becomes

$$\begin{aligned} \frac{S_qa}{q-2n}x^{q-2n}. \end{aligned}$$

This energy is liberated in the volume

$$\begin{aligned} \int _0^x{S_qx'^{q-1}}\mathrm{d}x' = \frac{S_q}{q}x^q. \end{aligned}$$

Hence the average energy concentration per unit mass is

$$\begin{aligned} e_a(x) = \frac{q}{q-2n}ax^{-2n}, \end{aligned}$$

and this exceeds the concentration e(x) at x by a factor

$$\begin{aligned} b = \frac{e_a(x)}{e(x)} = \frac{q}{q-2n}. \end{aligned}$$
(19)

This factor b represents in our model the ratio total (i.e., booster plus charge) : charge. Hence

$$\begin{aligned} \text { booster : charge }= b - 1 = \frac{2n}{q-2n}. \end{aligned}$$
(20a)

9. The formulae of 7. can now be used in this way:

An explosive may be considered “sluggish,” if its detonation is not self-supporting. This means that the stationary detonation wave of this substance is not sufficient to initiate in it the explosive reaction in an adjacent layer. That is, the pressure \(p_s\) or the mass velocity \(V_s\)—whichever may be the significant quantity—is not sufficient to produce that result. Denote the values which are actually needed to that end by p and V. Express these quantities by (10) or (13), thus introducing the quantity \(w>1\). This \(w>1\) gives then a quantitative measure of the “sluggishness” of the explosive.

Now (17), (20a) show how w determines the relative amount of booster, which is required to detonate a given quantity of this explosive.

These formulae express that a given quantity of booster will only detonate a definite amount of the explosive, and no more. And this amount will also depend on the value of \(q = 1, 2, 3\), i.e., on the character of the detonation wave: Plane, cylindrical, spherical.

3.8 Sec. 8–9 discussion

Using his assumed power-law \(e_0=ax^{-2n}\), von Neumann computes the amount of energy liberated in the volume of the explosion \(e_0\) up to a point x. This exceeds the energy liberated at point x by a factor b, which is the total amount of material (booster and charge) divided by the amount of charge. If the energetic output of the explosive materials and the desired dimension of the explosive, q, are known, then one can find b, the minimum necessary amount of booster material, by rearranging (19) to obtain (20a).

3.9 Sec. 10–13 transcription

10. We now pass to the computations.

Denote the time by t, and the position of the detonation front at that time by

$$\begin{aligned} \bar{x} = \bar{x}(t). \end{aligned}$$

Then by (12), (14)

$$\begin{aligned} \frac{d\bar{x}}{dt}= & {} D = \frac{w}{\sqrt{2w-1}}\sqrt{2(\gamma ^2-1)e_0} \\&= \frac{w}{\sqrt{2w-1}}\sqrt{2(\gamma ^2-1)a} \cdot \bar{x}^{-n}. \end{aligned}$$

Hence

$$\begin{aligned} \bar{x}^n\frac{d\bar{x}}{dt} = \frac{w}{\sqrt{2w-1}}\sqrt{2(\gamma ^2-1)a} \end{aligned}$$

and (assuming \(\bar{x}(0) = 0\), i.e., that the detonation began at \(t = 0\)),

$$\begin{aligned} {\left\{ \begin{array}{ll} t = \alpha \bar{x}^{n+1} \text { or } \bar{x} = \alpha ^{-\frac{1}{n+1}}t^\frac{1}{n+1} \\ \text {with} \\ \alpha = \frac{\sqrt{2w-1}}{w} \frac{1}{(n+1)\sqrt{2(\gamma ^2-1)a}}. \end{array}\right. } \end{aligned}$$
(20b)

Further (14) and (10), (11), (13) give

$$\begin{aligned} p= & {} w\frac{2(\gamma -1)a}{v_0}\bar{x}^{-2n}, \end{aligned}$$
(21)
$$\begin{aligned} v= & {} \frac{\gamma +w^{-1}-1}{\gamma +1}v_0, \end{aligned}$$
(22)
$$\begin{aligned} V= & {} \sqrt{2w-1}\sqrt{\frac{2(\gamma -1)a}{\gamma +1}}\bar{x}^{-n}. \end{aligned}$$
(23)

(21), (22) determine the adiabatic coefficient, with which the gas originates at this point of the detonation wave:

$$\begin{aligned} k = pv^{\gamma } = w \cdot 2(\gamma -1)a \cdot (\frac{\gamma +w^{-1}-1}{\gamma +1})^{\gamma } \cdot v_0^{\gamma -1} \cdot \bar{x}^{-2n}. \end{aligned}$$
(24)

11. The gas behind the detonation wave, i.e., in the interval

$$\begin{aligned} 0<x<\bar{x} = \alpha ^{-\frac{1}{n+1}}t^{\frac{1}{n+1}}, \end{aligned}$$
(25)

is governed by the Lagrangian differential equation

$$\begin{aligned} \text {X}_{tt} \text {X}_x = -vp_x. \end{aligned}$$
(26)

Here

$$\begin{aligned} \text {X} = \text {X}(x, t) \end{aligned}$$
(27)

is the position at the time t of that gas particle which was at x at the moment \(t = 0\). v, p obtain from the formulae

$$\begin{aligned} v= & {} v_0x^{-(q-1)} \text {X}^{q-1} \text {X}_x, \end{aligned}$$
(28)
$$\begin{aligned} p= & {} kv^{-\gamma }, \end{aligned}$$
(29)

where \(k = k(\bar{x})\) obtains from (24), since the motion of the gas is adiabatic throughout the zone (25) behind the detonation wave.

Hence (26) becomes

$$\begin{aligned} \text {X}_{tt}= & {} -w \cdot 2(\gamma -1)a \cdot (\frac{\gamma +w^{-1}-1}{\gamma +1})^{\gamma } \cdot \\&\cdot x^{-(q-1)} \text {X}^{q-1} \cdot (x^{-2n+\gamma (q-1)} \text {X}^{-\gamma (q-1)} \text {X}_{x}^{-\gamma })_x, \end{aligned}$$

i.e., using (20),

$$\begin{aligned} {\left\{ \begin{array}{ll} \text {X}_{tt} = -A \cdot x^{-(q-1)} \text {X}^{q-1}(((n+1)\alpha x^n)^{-2}x^{\gamma (q-1)} \text {X}^{-\gamma (q-1)} \text {X}_{x}^{-\gamma })_x \\ \text {with} \\ A = \frac{2w-1}{w} \frac{1}{\gamma +1}(\frac{\gamma +w^{-1}-1}{\gamma +1})^{\gamma }. \end{array}\right. } \end{aligned}$$
(30)

To the differential Eq. (30), we must add the boundary conditions. They correspond to the two ends

$$\begin{aligned} x = 0 \text { and }x = \bar{x} = \alpha ^{-\frac{1}{n+1}} t^{\frac{1}{n+1}} \end{aligned}$$

of (25). We have at

$$\begin{aligned} x = 0 : \text {X} = 0, \end{aligned}$$
(31)

and at

$$\begin{aligned} {\left\{ \begin{array}{ll} x = \bar{x} = \alpha ^{-\frac{1}{n+1}} t^{\frac{1}{n+1}} : \text {X} = x, \text {X}_x = B \\ \text { with } B = \frac{\gamma +w^{-1}-1}{\gamma +1}. \end{array}\right. } \end{aligned}$$
(32)

Note that (31) and the first equation of (32) express the fit with respect to position at the two ends of (25), while the second equation of (32) expresses the fit with respect to specific volume at the free end of (25). The last obtains from (22), (28), remembering the \(\text {X} = x\) at that place. The fit with respect to mass velocity at the same place could also be expressed, from (23) with \(V = \text {X}_t\). But this condition must be (and is) a consequence of the corresponding one for specific volume, since the conservation of matter has been safeguarded throughout our procedure.

Thus the differential Eq. (30), with boundary conditions (31), (32), describes the motion of the gas in (25), i.e., behind the detonation wave.

12. The similitude already referred to in (6) suggests that the connection of

$$\begin{aligned} \frac{x}{\bar{x}} \text { and }\frac{\text {X}}{\bar{x}} \end{aligned}$$

be independent of t. This means that we assume

$$\begin{aligned} \text {X} = \text {X}(x, t) = \alpha ^{-\frac{1}{n+1}} t^{\frac{1}{n+1}} f\left( \frac{x}{\alpha ^{-\frac{1}{n+1}}t^{\frac{1}{n+1}}}\right) : \end{aligned}$$
(33)

It is convenient to introduce

$$\begin{aligned} z = \frac{x}{\bar{x}} = \frac{x}{\alpha ^{-\frac{1}{n+1}} t^{\frac{1}{n+1}}}. \end{aligned}$$
(34)

Now the differential Eq. (30) becomes

$$\begin{aligned} {\left\{ \begin{array}{ll} z^2f_{zz}(z) + nzf_z(z) - nf(z) \\ = -Az^{-(q-1)}f(z)^{q-1}(z^{-2n+\gamma (q-1)} f(z)^{-\gamma (q-1)} f_z(z)^{-\gamma })_z. \end{array}\right. } \end{aligned}$$
(35)

And the boundary conditions (31), (32) become

$$\begin{aligned} z= & {} 0 : f(z) = 0 \end{aligned}$$
(36)
$$\begin{aligned} z= & {} 1 : f(z) = 1, f_z(z) = B. \end{aligned}$$
(37)

In (30), (32) we expressed A, B in terms of w. It is now convenient to express A, w in terms of B:

$$\begin{aligned} A= & {} (1-B) \cdot B^{\gamma }, \end{aligned}$$
(38)
$$\begin{aligned} w= & {} \frac{1}{(\gamma +1)B - (\gamma -1)}. \end{aligned}$$
(39)

13. Instead of making a direct attempt to integrate (35) with (36), (37), we first study the conditions at (36) somewhat more closely.

Throughout (25), the pressure can be expressed by (29) and (24), (28). The same computations by which the differential equation was derived in (11) give

$$\begin{aligned} p = \frac{A}{v_0} ((n+1)\alpha x^n)^{-2} x^{\gamma (q-1)} \text {X}^{-\gamma (q-1)} \text {X}_x^{-\gamma }, \end{aligned}$$

hence by (33), (34)

$$\begin{aligned} p = \frac{A}{(n+1)^2v_0} \cdot \alpha ^{\frac{2}{n+1}} t^{-\frac{2n}{n+1}} \cdot z^{-2n+\gamma (q-1)} f(z)^{-\gamma (q-1)} f_z(z)^{-\gamma }. \end{aligned}$$
(40)

A simple discussion of (35), which will not be given here, shows that \(p \rightarrow 0\) and \(p \rightarrow \infty \) are both impossible. Hence p is asymptotically like \(z^0\) for \(z \rightarrow \infty \). By (40) the same is true for \(z^{-2n+\gamma (q-1)} f(z)^{-\gamma (q-1)} f_z(z)^{-\gamma }\). Hence \(f(z)^{q-1}f_z(z)\) is asymptotically like \(z^{q-1-\frac{2n}{\gamma }}\). Now (36) necessitates \(q-1-\frac{2n}{\gamma } > -1\), i.e., \(n < \frac{\gamma q}{2}\), but this follows from (18). Further, (36) permits to infer from the above, that \(f(z)^q\) is asymptotically like \(z^{q-\frac{2n}{\gamma }}\). Hence f(z) is asymptotically like \(z^{1-\frac{2n}{\gamma q}}\).

So we can replace (36) by this stronger requirement (C undetermined but \(> 0\)):

$$\begin{aligned} f(z) = Cz^{1-\frac{2n}{\gamma q}} + \cdots \text { for }z \rightarrow \infty . \end{aligned}$$
(41)

3.10 Sec. 10–13 discussion

Here, von Neumann uses standard methods to solve a first-order in time nonlinear equation for \(\bar{x}(t)\), the position of the detonation front at time t. Pressure is bounded, so \(p \rightarrow 0\) and \(p \rightarrow \infty \) are both physically unrealizable. By (35), this means that p behaves asymptotically like \(z^0\) for \(z \rightarrow \infty \).

3.11 Sec. 14 transcription

14. It is convenient to put

$$\begin{aligned} f(z) = \mu g(vz). \end{aligned}$$
(42)

This leaves the form of (35) and (41) unaffected, except that it multiplies A and C by

$$\begin{aligned} \mu ^{-((\gamma -1)q+2)} v^{-((\gamma -1)q-2n)} \text { and }\mu ^{-1} v^{-(1-\frac{2n}{\gamma })}. \end{aligned}$$

Thus we can choose \(\mu \), v so as to make both these coefficients equal to 1. Since C is undetermined, we may use v instead of C as the undetermined quantity. In this way,

$$\begin{aligned} \mu= & {} A^{\frac{1}{(\gamma -1)q+2}} v^{-\frac{(\gamma -1)q-2n}{(\gamma -1)q+2}}, \end{aligned}$$
(43)
$$\begin{aligned} C= & {} A^{\frac{1}{(\gamma -1)q+2}} v^{\frac{2(n+1)}{(\gamma -1)q+2} - \frac{2n}{\gamma }}, \end{aligned}$$
(44)

ensue.

In this way, the differential Eq. (35) becomes

$$\begin{aligned} {\left\{ \begin{array}{ll} z^2g_{zz}(z) + nzg_z(z) - ng(z) \\ = -z^{-(q-1)} g(z)^{q-1} (z^{-2n+\gamma (q-1)} g(z)^{-\gamma (q-1)} g_z(z)^{-\gamma })_z, \end{array}\right. } \end{aligned}$$
(45)

and the boundary condition (36), in its stronger form (41), becomes

$$\begin{aligned} g(z) = z^{1-\frac{2n}{\gamma q}} + \cdots \quad \text { for }z \rightarrow 0. \end{aligned}$$
(46)

The remaining boundary condition (37) becomes now:

$$\begin{aligned} {\left\{ \begin{array}{ll} z = v: \quad g(z) = A^{-\frac{1}{(\gamma -1)q+2}} \cdot v^{\frac{(\gamma -1)q-2n}{(\gamma -1)q+2}}, \\ \quad \quad \quad \quad g_z(z) = A^{-\frac{1}{(\gamma -1)q+2}}B \cdot v^{-\frac{2(n+1)}{(\gamma -1)q+2}}. \end{array}\right. } \end{aligned}$$
(47)

Since v is arbitrary, we can formulate (37) like this: There must exist a z for which

$$\begin{aligned} g(z)= & {} A^{-\frac{1}{(\gamma -1)q+2}} \cdot z^{\frac{(\gamma -1)q-2n}{(\gamma -1)q+2}}, \end{aligned}$$
(48)
$$\begin{aligned} g_z(z)= & {} A^{-\frac{1}{(\gamma -1)q+2}}B \cdot z^{-\frac{2(n+1)}{(\gamma -1)q+2}}. \end{aligned}$$
(49)

This z then determines v by

$$\begin{aligned} v = z. \end{aligned}$$
(50)

We formulate (48), (49):

$$\begin{aligned} B= & {} \frac{zg_z(z)}{g(z)}, \end{aligned}$$
(51)
$$\begin{aligned} A= & {} z^{(\gamma -1)q-2n} \cdot g(z)^{-((\gamma -1)q+2)}. \end{aligned}$$
(52)

As we saw at the end of (12), B is an undetermined quantity, too. So we can interpret (51) as defining B, and eliminate A, B from (51), (52) by means of (38). In this way, the condition

$$\begin{aligned} {\left\{ \begin{array}{ll} z^{-(\gamma -1)(q-1)+(2n+1)} \cdot g(z)^{(\gamma -1)(q-1)} \cdot g_z(z)^{\gamma } \cdot \\ \quad \quad \cdot (g(z) - zg_z(z)) = 1 \end{array}\right. } \end{aligned}$$
(53)

obtains.

We can also substitute (51) in (39). This gives

$$\begin{aligned} w = \frac{g(z)}{(\gamma +1)zg_z(z)-(\gamma -1)g(z)}. \end{aligned}$$
(54)

3.12 Sec. 14 discussion

The manuscript’s multipliers of A and C appear to contain a minor error. When we substitute \(\mu g(vz)\) into (35), we obtain

$$\begin{aligned}&z^2\mu v^2g_{zz}(vz)+nz\mu g_z(vz)-n\mu g(vz)\\&\quad =-Az^{-(q-1)}\mu ^{q-1}g(vz)^{q-1}(z^{-2n+\gamma (q-1)}\mu ^{-\gamma (q-1)}g(vz)\mu ^{-\gamma }v^{-\gamma }g_z(vz)^{-\gamma })_z \end{aligned}$$

Substituting in \(\frac{z}{v}\) for z, we obtain

$$\begin{aligned}&\frac{z}{v}^2\mu v^2g_{zz}(z)+n\frac{z}{v}\mu g_z(z)-n\mu g(z)=\\&\quad -A\frac{z}{v}^{-(q-1)}\mu ^{q-1}g(z)^{q-1}(\frac{z}{v}^{-2n+\gamma (q-1)}\mu ^{-\gamma (q-1)}g(z)\mu ^{-\gamma }v^{-\gamma }g_z(z)^{-\gamma })_{\frac{z}{v}} \end{aligned}$$

Pulling u and v terms out does then yield the multiplier of A given in the manuscript, leaving the form of (35) unaffected.

However, when we substitute \(\mu g(vz)\) into (41), we obtain

$$\begin{aligned} \mu g(vz)= & {} Cz^{1-\frac{2n}{\gamma q}}\\ g(vz)= & {} \mu ^{-1}Cz^{1-\frac{2n}{\gamma q}}. \end{aligned}$$

Again, substitute in \(\frac{z}{v}\) for z; then

$$\begin{aligned} g(z)= & {} \mu ^{-1}C\frac{z}{v}^{1-\frac{2n}{\gamma q}}\\ g(z)= & {} \mu ^{-1}v^{-(1-\frac{2n}{\gamma q})}Cz^{1-\frac{2n}{\gamma q}}. \end{aligned}$$

The multiplier given by von Neumann lacks the q in the power of v; the q is again missing in (44). However, the boundary condition given in (46) is unaffected by this minor error. The rest of the section follows as stated by von Neumann.

3.13 Sec. 15–17 transcription

15. Summing up:

After \(\gamma \), q and n are chosen, the differential Eq. (45) must be integrated, beginning at \(z = 0\) with (46). The solution must be continued up to the point z where (53) holds. This z then determines w by (54). This is the process by which (17) is obtained.

16. The process of integrating (45) can be simplified by putting

$$\begin{aligned} z= & {} e^s \end{aligned}$$
(55)
$$\begin{aligned} g(z)= & {} e ^{Ps} \mu (s) \end{aligned}$$
(56)

where

$$\begin{aligned} P = \frac{(\gamma -1)q-2n}{(\gamma -1)q+2}. \end{aligned}$$
(57)

Then (45) becomes

$$\begin{aligned} {\left\{ \begin{array}{ll} (\frac{d}{ds} + (P-1))(\frac{d}{ds} + P)\mu (s) + n(\frac{d}{ds} + P)\mu (s) - n\mu (s) \\ = -\mu (s)^{q-1}(\frac{d}{ds} + Q)(\mu (s)^{-\gamma (q-1)}((\frac{d}{ds} + P)\mu (s))^{-\gamma }), \end{array}\right. } \end{aligned}$$
(58)

where

$$\begin{aligned} Q = -2n + \gamma q(1-P) = 2\frac{(\gamma +n)q - 2n}{(\gamma -1)q + 2}. \end{aligned}$$
(59)

Now put

$$\begin{aligned} w(s) = (\frac{d}{ds} + P)\mu (s), \end{aligned}$$
(60)

so that

$$\begin{aligned} g_z(z) = e^{(P-1)s} w(s). \end{aligned}$$
(61)

Then

$$\begin{aligned} \frac{d}{ds} = \frac{d\mu }{ds} \frac{d}{d\mu } = (w-P\mu )\frac{d}{d\mu }, \end{aligned}$$
(62)

and so (58) becomes

$$\begin{aligned}&((w - P\mu )\frac{d}{d\mu } + (P-1))w + nw - n\mu \\&\quad = -\mu ^{q-1}((w-P\mu )\frac{d}{d\mu } + Q)(\mu ^{-\gamma (q-1)}w^{-\gamma }), \end{aligned}$$

i.e.,

$$\begin{aligned} {\left\{ \begin{array}{ll} \frac{dw}{d\mu } \\ = \frac{1}{w-P\mu } \frac{[(1-P-n)w+n\mu ] + [\gamma (q-1)w-(\gamma (q-1)P+Q)\mu ]\mu ^{-((\gamma -1)(q-1)+1)} w^{-\gamma }}{1 - \gamma \mu ^{-((\gamma -1)(q-1)+1)} w^{-(\gamma +1)}}. \end{array}\right. } \end{aligned}$$
(63)

17. (62) gives further

$$\begin{aligned} s = \int {\frac{d\mu }{w-P\mu } + \text { Constant }}, \end{aligned}$$
(64)

and \(z \rightarrow 0\) clearly means \(s \rightarrow -\infty \) by 55. (46) becomes, by (55), (56)

$$\begin{aligned} \mu (s) = e^{\frac{Q}{\gamma q}s} + \cdots \quad \text { for }s \rightarrow -\infty . \end{aligned}$$
(65)

Similarly, differentiation of (46) gives, by (55), (60)

$$\begin{aligned} w(s) = (1-\frac{2n}{\gamma q}) e^{\frac{Q}{\gamma q} s} + \cdots \quad \text { for }s \rightarrow \infty . \end{aligned}$$
(66)

Before we discuss (65), (66) any further, let us also reformulate (53) and (54). (53) becomes

$$\begin{aligned} \mu ^{(\gamma -1)(q-1)} \cdot w^{\gamma } \cdot (\mu -w) = 1, \end{aligned}$$

i.e.,

$$\begin{aligned} \mu - w = \mu ^{-(\gamma -1)(q-1)} w^{-\gamma }. \end{aligned}$$
(67)

(54) becomes

$$\begin{aligned} w = \frac{\mu }{(\gamma +1)w - (\gamma -1)\mu }. \end{aligned}$$
(68)

3.14 Sec. 15–17 discussion

The goal of these sections is to solve for the Lagrangian coordinate X(x, t), which describes the position at time t of a gas particle that was initially at position x. Von Neumann obtains the solution through a series of variable transformations, substitutions, asymptotic arguments, and integrations.

3.15 Sec. 18–19 transcription

18. Let us now return to (65), (66).

By (18) the numerator of Q in (59) is

$$\begin{aligned} (\gamma +n)q-2n = \gamma q + (q-2)n {\left\{ \begin{array}{ll} \ge \gamma q> 0 \\ \text { for }q = 1, 2, \\ \end{array}\right. } {\left\{ \begin{array}{ll} \\ = 3\gamma - n \ge 3\gamma - \frac{3}{2}> \\ > 0 \\ \text { for }q = 3. \end{array}\right. } \end{aligned}$$

Consequently in any event

$$\begin{aligned} Q > 0. \end{aligned}$$
(69)

Thus (65), (66) imply that on that boundary point

$$\begin{aligned} u, w \rightarrow 0, \frac{w}{u} \rightarrow 1 - \frac{2n}{\gamma q}. \end{aligned}$$
(70)

Conversely: Combining (70) with (64) gives

$$\begin{aligned}&s \sim \int {\frac{du}{(1-\frac{2u}{\gamma q}-P)u}} = \int {\frac{du}{\frac{Q}{\gamma q} u}} \\&\quad = \frac{\gamma q}{Q} \ln {u} + \text { Constant }. \end{aligned}$$

We can adjust this constant so that

$$\begin{aligned} s \sim \frac{\gamma q}{Q} \ln {u} \end{aligned}$$

results, hence \(s \rightarrow -\infty \), \(u \sim e^{\frac{Q}{\gamma q} s}\), and so by (70) \(w \sim (1 - \frac{2n}{\gamma q})e^{\frac{Q}{\gamma q}s}\). Thus (65), (66) are valid.

In other words: (65), (66) may be replaced by (70).

19. Summing up: (Cf. (15))

After \(\gamma \), q and n are chosen, the differential Eq. (63) must be integrated, beginning with (70).

The solution must be continued up to the point u, w where (67) holds. These u, w then determine w by (68).

This is the process by which (17) is obtained.

3.16 Sec. 18–19 discussion

Although \(Q>0\) in (59), von Neumann’s statement in the unnumbered equation between (68) and (69) appears incorrect. Our own analysis is as follows.

When \(q=1\), the numerator of Q is greatest when \(n\rightarrow 0\) and smallest when \(n\rightarrow q/2\). In fact, Q can be less than 0 for \(\gamma <1/2\); however, the quantity \(\gamma \) is always greater than 1 by definition. Recall that \(\gamma \) is the ratio of the heat capacity of an ideal gas at constant pressure to the heat capacity at constant volume. This quantity can also be expressed as \(1+2/f\), where f is the degrees of freedom of a molecule, a positive value. Hence \(\gamma > 1\), and \(Q<0\) is unrealizable. When \(q=2\), the numerator of Q is always \(\gamma q\), and \(Q>0\) for all n. When \(q=3\), the numerator of Q is smallest when \(n\rightarrow 0\) and greatest when \(n\rightarrow q/2\). Similarly to the case when \(q=1\), Q can be less than zero for \(\gamma <1/2\).

In these final sections, von Neumann shows that the differential equation \(dw/d\mu \) found in (63) can be integrated, and by the inequalities above,

$$\begin{aligned} w \sim \left( 1 - \frac{2n}{\gamma q}\right) e^{\frac{Q}{\gamma q}s}. \end{aligned}$$

Von Neumann’s goal was to compute the functional relationship between w and n, and he succeeds with the expression above.

4 Discussion and conclusions

In his published work on shock waves (von Neumann 1963j) that ultimately led to ZND theory, von Neumann shows “when the so-called Chapman-Jouguet hypothesis is true, and what formulae are to be used when it is not true.” He states later in von Neumann (1963j) that “it is hoped that this will connect the present theory with the difficult questions of initiating a detonation of primers and boosters,” and in von Neumann (1963k) he states that he will publish a subsequent report on boosters. We could not locate any such work specific to boosters, and we offer two hypothesis for the apparent absence of this publication.

First, it is possible that von Neumann did finish a manuscript dealing with the topic of boosters, but that it has not been declassified. The two earliest works in Table  1, that is, von Neumann (1963j, (1963k), were originally restricted documents. However, both were declassified in the early 1950s, and we did not find other von Neumann documents that have been declassified since that time.

We offer a second, and we believe more likely, explanation for the absence of a published report on boosters, namely, that (von Neumann unpublished) was meant to be this work, but was eventually abandoned due to its scientific limitations. It is difficult to place the exact time of the writing of von Neumann (unpublished). On the one hand, we might assume that it is the paper on boosters von Neumann hinted at in 1942 and 1943 (von Neumann 1963j, k) and that he was working on it simultaneously to or immediately following these works. On the other hand, von Neumann (1963j, (1963k) take into account important aspects of detonations that von Neumann does not include in von Neumann (unpublished), perhaps suggesting that (von Neumann unpublished) came first. For example, the 1943 publication (von Neumann 1963k) specifically states the differences in modeling that must be accounted for in one, two or three dimensions; these are not addressed in von Neumann (unpublished). Additionally, although von Neumann states repeatedly in von Neumann (unpublished) that the CJ hypothesis is generally true, by 1942 he had already shown when it failed as a predictive method. He is careful to note this failure in subsequent works (von Neumann 1963b), and yet no mention of it is made in von Neumann (unpublished). Moreover, the majority of von Neumann’s works on the theory of shock waves following (von Neumann 1963j, k) (excluding those concerned with shock wave reflection) have the following characteristics:

  • incorporation of thermodynamic information that is absent in von Neumann (unpublished),

  • substantial usage of numerical methods rather than analytical ones,

  • formulation for detonation in a particular dimension, and

  • points or homogeneous spheres as the source of detonation.

In most of his shock wave papers, von Neumann emphasizes the infancy of the theory and the difficulty in finding theories that are generalizable and physically justifiable. The simplicity of the models in his published works perhaps suggests why von Neumann never published the more ambitious manuscript we have analyzed here. While von Neumann’s step-by-step derivations in the unpublished manuscript are mostly accurate, several fundamental assumptions of the model cannot be justified and the nuance of his discussion on the topic in his published works makes clear that he was aware of these difficulties. For example, his liberal use of the CJ hypothesis, the simplification of the discontinuity between the booster and the main charge in terms of energetic output to a monotonic decreasing function, the lack of specific thermodynamic laws in the energetic output equations, the generalization to multiple dimensions, and the exclusion of \(p_0\) without explanation are much less sophisticated than the work done in his published papers.

Von Neumann was diagnosed with cancer in 1955 and subsequently placed under security for fear that he would “reveal military secrets while heavily medicated,” passing away in 1957 (Macrae 2019). To date, there is no generalized theory regarding boostered detonations like the one von Neumann hoped to create in von Neumann (unpublished). In the 1960s, experimental observation made clear that the one dimensional theories could not, in general, capture the complex structures of three dimensional detonations. Detonations have asymmetrical dynamics specific to the material used. Those materials and their chemical properties result in complications that prohibit the creation of a generalizable theory. Instead, investigations typically came to rely on experiments and numerical simulations. From that point of view, the handwritten manuscript is not an undiscovered model that will propel forward the modern understanding of explosive shock waves; however, it is an intriguing analysis performed by a mathematical powerhouse. The manuscript (von Neumann unpublished) sheds light on the scientific hurdles faced in the development of detonation theory during a time of pressure to understand explosives as they were being rapidly developed and produced during World War II.