Spaceship Flight

clc; clear all; close all;

%%%%%%%%%%%%%%%%
%% Initialize %%
%%%%%%%%%%%%%%%%

ShipMass    = 75000; % kg
OrbDist     = 2000000; % meters
FuelDensity = 4500; % N*s/kg
FuelRate    = [500:20:8000]; % kg/s
N           = length(FuelRate);
dt          = 2; % seconds
g           = 9.8; % m/s^2

%%%%%%%%%%%%%%%%%%%%%%
%%%% Calculations %%%%
%%%%%%%%%%%%%%%%%%%%%%

for i=1:N

   %%%%%%%%%%%%%%%%%%%
   %% Non-Constants %%
   %%%%%%%%%%%%%%%%%%%

   FuelMass     = 5000000; % kg
   Position     = 0; % meters
   Velocity     = 0; % m/s
   Acceleration = 0; % m/s^2
   G            = 0; % N
   GForce(i)    = 0; % Stores value each trajectory
   OrbitTime(i) = 0; % Stores value each trajectory

   while Position < OrbDist

       TotalMass    = ShipMass + FuelMass;
       PropForce    = FuelRate(i).*FuelDensity;
       Acceleration = PropForce./TotalMass;
       Velocity     = Velocity + ((1/2).*Acceleration.*dt.^2);
       Position     = Position + Velocity.*dt;
       G            = Acceleration./g;

       if G > GForce(i)

           GForce(i) = G;

       end

       FuelMass     = FuelMass - FuelRate(i).*dt;
       OrbitTime(i) = OrbitTime(i)+dt;

   end
end

OrbitTime = OrbitTime./60 % Change time from seconds to minutes

figure(1)
plot (GForce, FuelRate, 'b')   % Shows what the Fuel Rate would be at
                               % different g-forces

figure(2)
plot (OrbitTime, GForce, 'g')  % Shows g-force based on time
                               % until orbit

figure(3)
plot (OrbitTime,FuelRate, 'r') % Shows Fuel Rate based on time
                               % until orbit