I'm trying to model a rocket launch with MathematicaMathematica but I've run into a little problem since I don't know how to turn the thrust off. I'm using Newton's Law of Universal Gravitation plus an added term for the rocket's thrust and trying to find its position as a function of time. Here is the code so far (it is at a very basic stage at the moment with constant mass and no drag added, plus the thrust never stops, which is my biggest concern right now):
(*Gravitational Constant*)
G = 6.672*10^-11
(*Mass of Earth and rocket*)
M = AstronomicalData["Earth"AstronomicalData["Earth", "Mass"]"Mass"]
m = 2800000
(*Rocket thrust*)
T = 34020000
(*Radius of Earth*)
r = AstronomicalData["Earth"AstronomicalData["Earth", "Radius"]"Radius"]
(*Numerical solution modelling the gravitation interation between the \
Earth and a launching rocket*)
(*NOTE: Rocket mass will change over time; also, add in drag*)
soln = NDSolve[{
x''[t] == -((G M x[t])/(x[t]^2 + y[t]^2 + z[t]^2)^(3/2)),
y''[t] == -((G M y[t])/(x[t]^2 + y[t]^2 + z[t]^2)^(3/2)) + 0.25 T/m,
z''[t] == -((G M z[t])/(x[t]^2 + y[t]^2 + z[t]^2)^(3/2)) + 0.75 T/m,
x[0] == 0, y[0] == 0, z[0] == r, x'[0] == 0, y'[0] == 0,
z'[0] == 0}, {x[t], y[t], z[t]}, {t, 0, 20000},
MaxSteps ->> 10000000, Method ->> "StiffnessSwitching"]"StiffnessSwitching"]
RocketPlot =
ParametricPlot3D[{x[t], y[t], z[t]} /. soln, {t, 0, 10000000}, AxesLabel ->> {x, y, z}]
Does anyone know if conditional statements can be used inside NDSolveNDSolve so that the thrust can be stopped at a certain position or time? Any help would be appreciated, cheers guys.
