Solving nonlinear differential to emulate Newtonian gravity

ellipsis
Messages
158
Reaction score
24
This is part of a personal project... I've recently become addicted to modeling various physical systems from scratch, such that I find explicit solutions of position as a function of time, and graph em.

But I've hit a glass ceiling trying to find an analytic solution to the 1-dimensional problem free-fall problem where acceleration due to gravity is the inverse square of the distance:
$$
a = \frac{-1}{x^2}
$$
To keep it as simple as possible, I'm assuming the initial velocity is zero, and the initial position is positive. With a little differential manipulation I found:
$$
v = -\sqrt{\frac{2}{x}-\frac{2}{x_i}}
$$
I'm 99% confident that's correct where ##0<x\leq x_i##, and it seems intuitive. It follows from that, then:

$$
t = \int \frac{-1}{\sqrt{\frac{2}{x}-\frac{2}{x_i}}}\,dx
$$
Now this is a nasty integral, but it does have an explicit solution. You just cannot explicitly solve for x, but I'm perfectly happy with t in terms of x rather than x in terms of t (Because I can just convert it to a parametric form).

The only problem is: Wolfram_Alpha and MATLAB disagree on what the integral is, and I have no idea what to do with the constant, besides. Also, Wolfram_Alpha returns an expression with complex values, which I don't understand or want.

If anybody knows how to do this, any help would be vastly appreciated... I've only formally taken up to Calculus II, and I've been struggling with this for weeks now. (I fear I simply don't have the required math knowledge yet... e.g. what is a Lagrangian?)
 
Physics news on Phys.org
(X-post from Physics)

This is part of a personal project... I've recently become addicted to modeling various physical systems from scratch, such that I find explicit solutions of position as a function of time, and graph em.

But I've hit a glass ceiling trying to find an analytic solution to the 1-dimensional problem free-fall problem where acceleration due to gravity is the inverse square of the distance:
$$
a = \frac{-1}{x^2}
$$
To keep it as simple as possible, I'm assuming the initial velocity is zero, and the initial position is positive. With a little differential manipulation I found:
$$
v = -\sqrt{\frac{2}{x}-\frac{2}{x_i}}
$$
I'm 99% confident that's correct where ##0<x\leq x_i##, and it seems intuitive. It follows from that, then:

$$
t = \int \frac{-1}{\sqrt{\frac{2}{x}-\frac{2}{x_i}}}\,dx
$$
Now this is a nasty integral, but it does have an explicit solution. You just cannot explicitly solve for x, but I'm perfectly happy with t in terms of x rather than x in terms of t (Because I can just convert it to a parametric form).

The only problem is: Wolfram_Alpha and MATLAB disagree on what the integral is, and I have no idea what to do with the constant, besides. Also, Wolfram_Alpha returns an expression with complex values, which I don't understand or want.

If anybody knows how to do this, any help would be vastly appreciated... I've only formally taken up to Calculus II, and I've been struggling with this for weeks now. (I fear I simply don't have the required math knowledge yet... e.g. what is a Lagrangian?)
 
Your units do not agree in your very first equation.
GIGO
 
Here goes:

Rearrange the integral for t:
$$ t = - \int \sqrt{ \frac{x x_i}{2(x_i - x)}} dx $$

This suggests that we could use a change of variables into a trigonometric function: ## x = x_i \sin^2 u ##:
$$ t = - \sqrt{2 (x_i)^3} \int \sin^2 u \, du $$
giving
$$ t = - \sqrt{2 (x_i)^3} \frac12 (u - \sin u \cos u) $$

A word of warning about this solution: it cancels to first order near u = 0, giving ## O(u^3) ## behavior. You could try doing the integral numerically as an alternative.
 
  • Like
Likes ellipsis
Is that relevant to the solution of the differential equation? I neglected to add in realistic parameters to avoid unnecessary complexity, but if you insist...

$$
a = -\frac{GM}{(R+h)^2}
$$

I detect a certain dubiousness on your part, Integral, but I think that's fair. I'll explicitly go through the steps to find velocity in terms of current height, initial velocity, and initial height. So you know I'm not talking out of my butt here.

First, some obvious definitions:
$$
\begin{align}
a = \frac{dv}{dt}\\
v = \frac{dh}{dt}
\end{align}
$$

Via the chain rule:
$$
\begin{align}
a = \frac{dv}{dt} = \frac{dv}{dh}\frac{dh}{dt} = v\frac{dv}{dh}\\

v\frac{dv}{dh} = -\frac{GM}{(R+h)^2}\\
\\
v~dv = -\frac{GM}{(R+h)^2}~dh
\end{align}
$$

Now we integrate:
$$\begin{align}
\int v\,dv = \int -\frac{GM}{(R+h)^2}\,dh\\
\frac{1}{2}v^2 = \frac{GM}{R+h}+C_1
v^2 = \frac{2GM}{R+h}+C_1
\end{align}$$

(Now we must briefly solve for ##C_1##)
$$\begin{align}
v_i^2 = \frac{2GM}{R+h_i}+C_1\\
C_1 = -\frac{2GM}{R+h_i}+v_i^2
\end{align}$$

And replace ##C_1## in the original:
$$\begin{align}
v^2 = \frac{2GM}{R+h}-\frac{2GM}{R+h_i}+v_i^2
\end{align}$$

Finally, we obtain:
$$\begin{align}
v = \pm\sqrt{\frac{2GM}{R+h}-\frac{2GM}{R+h_i}+v_i^2}
\end{align}$$

All of the units should now check out. You'll have to use common sense to realize whether the velocity should positive (up) or negative (down).

In the case of Earth...

$$\begin{align}
v = \pm\sqrt{\frac{797200884}{6371+h}-\frac{797200884}{6371+h_i}+v_i^2}
\end{align}$$
A few sanity checks:
For initial height of 10, initial velocity 0, and final height 0, I get a final velocity of very nearly 14 m/s. Which is what the normal kinematic equations give.
 
That equation for a should be -k/x^2. Then you will have a ##\sqrt{2k}## in the equation, instead of a 2. Multiply numerator and denominator in the equation for t by ##\sqrt{xx_i}##, and then you can do the integration using trig substitution.

Chet
 
Thank you so much. You are a wizard of mathematics, a scholar of the arcane. Although I didn't end up using your trigonometric substitution trick (though, that might make for a more efficient simulation), your initial simplification allowed me to convince W|A to give me the right integral.

When I finish fixing it up I will post it here for posterity.
 
Here's the final solution I found, in parametric form:
$$\begin{cases}
x = \sqrt{\frac{x_i}{2}}\Bigg(\sqrt{Tx_i-T^2} - x_i\sin^{-1}{\Big(\sqrt{\frac{T}{x_i}}\Big)} + \frac{x_i\pi}{2}\Bigg)\\
y = T
\end{cases}$$
where x is x and y is t. When graphing, it is best to do it from ##T = 0 \to x_i##

I could really not have imagined a better result... I have spent months convinced that an analytical solution to this DE was impossible. In celebration, here's how it should look when graphed:

newton.gif
. . .

numerical.png
 
Last edited:
ellipsis said:
Here's the final solution I found, in parametric form:
$$\begin{cases}
x = \sqrt{\frac{x_i}{2}}\Bigg(\sqrt{Tx_i-T^2} - x_i\sin^{-1}{\Big(\sqrt{\frac{T}{x_i}}\Big)} + \frac{x_i\pi}{2}\Bigg)\\
y = T
\end{cases}$$
where x is x and y is t. When graphing, it is best to do it from ##T = 0 \to x_i##

Could you possibly roughly explain how you got that parametric equation? I got stuck here the first time I tried:
ellipsis said:
Finally, we obtain:
\begin{align} v = \pm\sqrt{\frac{2GM}{R+h}-\frac{2GM}{R+h_i}+v_i^2} \end{align}

And have not made any progress since.
 

Similar threads

Back
Top