# Solving nonlinear differential to emulate Newtonian gravity

Tags:
1. Sep 24, 2014

### ellipsis

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?)

2. Sep 24, 2014

### ellipsis

(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?)

3. Sep 24, 2014

### Integral

Staff Emeritus
GIGO

4. Sep 24, 2014

### lpetrich

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.

5. Sep 24, 2014

### 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.

6. Sep 24, 2014

### Staff: Mentor

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

7. Sep 24, 2014

### ellipsis

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.

8. Sep 24, 2014

### ellipsis

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:

. . .

Last edited: Sep 24, 2014
9. Aug 16, 2015

### tomdodd4598

Could you possibly roughly explain how you got that parametric equation? I got stuck here the first time I tried:
And have not made any progress since.