Question

Submarine Explosion
-----------------------------

A large mass of incompressible, inviscid fluid contains a spherical bubble obeying Boyle's Law:

p V = constant

At great distances from the bubble, the pressure is zero.
Neglecting body forces, show that the radius R(t) of the bubble at time t satisfies the equation:

(d/dt) (R^2 (d/dt)(R) ) - (1/2) R ((d/dt)(R))^2 = k/R^2

for a constant k

Answer 1

Solution 1.
------------

Here spherical symmetry applies and so:

phi(r, t) = (F(t)/r) + G(t)

Then we consider our boundary condition.
A unit normal to the boundary between the bubble and the fluid is e_r, and so:

grad(phi) dot n = (d/d r)(phi)

The expanding bubble represents a moving boundary, given by R = R(t), with velocity:

U = (d/dt)( R(t)) e_r

Applying the boundary condition at r=R, we have:

u(R,t) = (d/dt)(phi(R,t)) = - F(t)/R(t)^2 = (d/dt)(R(t))

so that F(t) = - R^2 (d/dt)(R) and

phi(r, t) = - (R(t)^2)/r (d/dt)(R(t)) + G(t), v(r, t) = R(t)^2/(r^2) (d/dt)(R(t)) e_r


Bernoulli's equation for unsteady irrotational flow under zero body forces takes the form:

p(r, t)/rho + (1/2)( R(t)^4 / r^4) (d/dt)(R(t)) - (1/r)(d/dt)( R(t)^2) (d/dt) (R(t)) + G'(t) = H(t)

Writing this as r goes to infinity, where the fluid is at zero pressure, leads to:
G'(t) = H(t).

Thus Bernoulli's equation reduces to:

p(r, t) / rho + (1/2) (R(t)^4 / r^4) ((d/dt)(R(t)))^2 - (1/r) (d/dt)( R(t)^2 (d/dt)(R(t)) ) = 0,

to hold everwhere.
In particular when r == R, i.e. at the boundary, we have:

p(R, t) / rho + (1/2)*((d/dt)(R(t)))^2 - (1/R(t)) (d/dt)(R(t)^2 (d/dt)(R(t)) ) = 0

Now at r == R, Boyle's law holds: pV = C, a constant:

p(R, t) = C/( (4/3) pi R^3 ),

so that

p(R, t)/ rho = k/ R^3,

for some constant k, and:

k / R(t)^3 = (1/R(t)) (d/dt) (R(t)^2 (d/dt)(R(t))) - (1/2) ( (d/dt)(R(t)) )^2.

Multiplying across by R, we obtain the required nonlinear equation, which we may solve numerically.

equation1 := (3/2) R(t) ( (d/dt)(R(t)) )^2 + R(t)^2 ((d/dt)^2(R(t)) ) - 10/R(t)^2 = 0

-----------------------
Maple code:
-----------------------

> restart: with(plots) : with(Detools):
k := 10
eq1 := diff( R(t)^2* diff(R(t), t) ,t) - 1/2*R(t)* (diff(R(t),t))^2 - k/R(t)^2 = 0;
p := dsolve( {eq1, R(0)=1, D(R)(0) = 0, R(t), type=numeric, range = 0..10} ):
odeplot(p);