Will Flannery
The freefall wiki entry wiki Freefall has an analytic solution for freefall distance in a gravitational field, but ... it doesn't seem to work ... at least i can't get it to work ... here is my MATLAB program to test it ...
clear
G=6.7e-11; % gravitational constant m^3/(kg*s^2)
mEarth = 5.9742e24; % mass of Earth in kg
rEarth = 6.378e6; % radius of Earth in m

% drop a 1kg ball from 100m above the surface of the earth
mu = G*(mEarth + 1);
y0 = rEarth + 100;

t = 1;
x = (3/2*(pi/2-t*sqrt(2*mu/y0^3)))*2/3;
y1 = y0*(x - x^2/5 - 3*x^3/175 - 23*x^4/7875 - 1894*x^5/3931875 - 3293*x^6/21896875 - 2418092*x^7/62077640625);

y0-y1

ans =

9.5639e+04

So, the apple fell 95000 m in the 1st second ... (or I made a mistake)
I tried to check the references, the first has a lot of formulas but not the one above, and the second is behind paywall (The Physics Teacher) and the abstract does not look promising.

Is there an analytic formula anywhere ?

Last edited:

Haborix
Your definition for x has a typo in the exponent at the end. As an aside, I recall a lot of these series having terrible convergence properties. You might try doing the calculation truncating the inverted series at different powers to see if the answer varies wildly.

Will Flannery
Thanks, corrected ...
x = (3/2*(pi/2-t*sqrt(2*mu/y0^3)))^2/3;

However, now the apple is falling upwards ...

y0-y1 = ans =

-2.2703e+04

The last term in the series ... 2418092*x^7/62077640625
ans =

0.0028

but ... *y0 = 1.8177e+04

So, it's a long way from converging ... However, as t increases, x decreases and convergence improves, so when t = 500 the last term (including multiplication by y0) is 0.1935 ! So, is it converging to the true value ?

Last edited:
2022 Award
Try ^(2/3). a^2/3 is interpreted as ##a^2/3##, but you want ##a^{2/3}##.

Last edited:
• davenn
there are a couple of old threads here where an analytical solution for position vs time is derived. It is messier than one might imagine.

Will Flannery
Try ^(2/3). a^2/3 is interpreted as ##a^2/3##, but you want ##a^{2/3}##.
By Jove ! ... and I've only been Matlabbing for 100 years !

G=6.7e-11; % gravitational constant m^3/(kg*s^2)
mEarth = 5.9742e24; % mass of Earth in kg
rEarth = 6.378e6; % radius of Earth in m

mu = G*(mEarth + 1);
y0 = rEarth;

t = 1;
x = (3/2*(pi/2-t*sqrt(2*mu/y0^3)))^(2/3);
y1 = y0*(x - x^2/5 - 3*x^3/175 - 23*x^4/7875 - 1894*x^5/3931875 ...
- 3293*x^6/21896875 - 2418092*x^7/62077640625);

y0-y1 = -2.9134e+04
y0*2418092*x^7/62077640625 = 1.3488e+04

so it's still not converging on the high end ...
for t = 500
y0 -y1 = 1.3217e+06
y0*2418092*x^7/62077640625 = 296.6925

So, it's better but still not so good.
Conclusion - there are not enough terms in this thing to tell if it works or not (or there is yet another bug ... )

Last edited:
Will Flannery
Just use LaTeX, and one can read you math!

https://www.physicsforums.com/help/latexhelp/
That's not math, that's MATLAB.

In any case, now the question is, why does every book ignore the problem of an analytical solution for the 1-D gravity problem (as noted in the ref. in the wiki article, the guy looked at 100 physics books) and what is the solution?

From the reference in the wiki article, the approach to a solution is thru Kepler's eq. Note that when eccentricity e approaches 1 in Kepler's eq, the orbit does not approach a parabola, Kepler's original eq. only works for ellipses, and as e -> 1 the orbit gets flatter and flatter and when e = 1 you get a straight line !, I just plotted one, So, apparently we can use either a numeric or analytic solution to Kepler's eq, with e = 1, to plot a 1-d trajectory. And so the questions for the day is ... how?

So, if we drop a ball from a height of 100m, a = 100+rEarth and we can plot the trajectory ...
(1.5 hours later) ... that seems to work ...
clear
G=6.7e-11; % gravitational constant m^3/(kg*s^2)
mEarth = 5.9742e24; % mass of Earth in kg
rEarth = 6.378e6; % radius of Earth in m

a = 100+rEarth;
T = sqrt(a^3 /(G*mEarth/(4*pi^2))); % from Kepler's 3rd Law
e = 1;

N=100;
for i = 1:4
E = 2*pi*(i-1)/N + pi; % Eccentric anomaly, start at apogee
M = E - e*sin(E); % Mean anomaly
t(i) = M*T/(2*pi)-T/2; % M(t) = 2*pi*t/T
x(i) = a*cos(E);
end

plot(t,-x-rEarth) % invert graph

• weirdoguy
Gold Member
2022 Award
Why are you writing when you don't want to be read? SCNR.

Will Flannery
Why are you writing when you don't want to be read? SCNR.
? That doesn't make sense, I want it to be read... even debugged if necessary ...

Generally, things become much clearer when you actually have to program them and check if they work ... For example, it took me over an hour, and many corrections, to get that program to work.

2022 Award
If it's code, put it in code tags.
Matlab:
X=1 % It's even got a syntax highlighter

• • Nugatory and vanhees71
If it's code, put it in code tags.
Code:
SCNR

• 