Как рассчитать падение с изменяющимся ускорением свободного падения?

Первое уравнение Вы записали правильно (надо еще минус добавить), а проинтергировали - неправильно (слева вторая производная по времени, а вы интегрируете по расстоянию).
Имеем дифференциальное уравнение второго порядка:
r'' = -GM / r^2
Делаем замену p(r) = r'
Тогда r'' = dp/dr * dr/dt = p'p (производная сложной функции)
p'p = -GM / r^2
p * dp = -GM / r^2 * dr
p^2 / 2 = GM / r + C1
Переходим обратно к r':
(r')^2 = GM / r + C1
Надо найти константу интегрирования С1:
Пусть в момент времени t = 0 скорость тела была равна 0, а высота равнялась r0
0 = GM / r0 + C1
C1 = -GM / r0
r' = sqrt(GM * (1 / r - 1 / r0)
dr / sqrt(GM * (1 / r - 1/ r0) = dt
dr / sqrt(r0 / r - 1) = dt * sqrt(GM / r0)
Делаем замену u = sqrt(r0 / r - 1)
u^2 = r0 / r - 1
r = r0 / (u^2 + 1)
dr = -r0 * 2u * du / (u^2 + 1)^2
Уравнение приобретает вид:
-r0 * 2du / (u^2 + 1)^2 = dt * sqrt(GM / r0)
Делаем еще одну замену:
u = tg z
Уравнение приобретает вид:
-2 * r0 * (cos z)^2 * dz = dt * sqrt(GM / r0)
В свою очередь (cos z)^2 = (1 + cos (2z)) / 2
Подставляем, интегрируем по z, переходим обратно к u, а от него к r. Находим вторую константу интегрирования из начальных условий, получаем ответ.
P. S. Такое движение можно представить как движение по вырожденному эллипсу и рассчитать все с помощью третьего закона Кеплера, но это зависит от конкретной задачи.

Без двойного интегрирования не обойтись

Интересно ...Самому бы узнать ...- Спасибо за вопрос ...))))

Нужно вычислить ускорение и время падения, а следом уже и скорость.