Code covered by the BSD License  

Highlights from
Two-dimensional, Low-thrust Earth-to-Mars Trajectory Analysis with MATLAB

Two-dimensional, Low-thrust Earth-to-Mars Trajectory Analysis with MATLAB

by

 

27 Jun 2013 (Updated )

Determines optimal, two-dimensional low-thrust Earth-to-Mars interplanetary trajectories.

ilt_shoot (x)
function [f, g] = ilt_shoot (x)

% objective function and equality constraints

% simple shooting method

% inputs

%  x(1) = current value of transfer time or throttle setting
%  x(2, nsegments) = current values of steering angle alpha

% outputs

%  f = vector of equality constraints and
%      objective function evaluated at x

% Orbital Mechanics with MATLAB

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

global iopt tfactor nsegments tof alpha_wrk xinitial throttle

if (iopt == 1)
    
    % minimize transfer time - fixed throttle setting
    
    throttle = 1.0;
    
else
    
    % maximize final mass - current throttle setting
    
    throttle = x(1);
    
end

% compute duration of each time interval (non-dimensional)

if (iopt == 1)
    
    % minimize transfer time option
    
    deltat = x(1) / nsegments;
    
    tof = x(1);
    
else
    
    % maximize final mass option
    
    deltat = tof / nsegments;
    
end

% specify number of differential equations

neq = 6;

% truncation error tolerance

tetol = 1.0e-10;

% initialize initial time

ti = -deltat;

% set initial conditions

yi(1) = xinitial(1);

yi(2) = xinitial(2);

yi(3) = xinitial(3);

yi(4) = xinitial(4);

yi(5) = 0.0;

yi(6) = 0.0;

% step size guess (non-dimensional time)

h = (1200.0 / 86400.0) / tfactor;

% integrate for all segments

for i = 1:1:nsegments
    
    alpha_wrk = x(i + 1);
    
    % increment initial and final times
    
    ti = ti + deltat;
    
    tf = ti + deltat;
    
    % integrate from current ti to tf
    
    yfinal = rkf78('ilt_eqm', neq, ti, tf, h, tetol, yi);
    
    % reset integration vector
    
    yi = yfinal;
    
    % check for end of simulation
    
    if (tf >= tof)
        
        break;
        
    end
    
end

% objective function (minimize accumulated delta-v)

f(1) = yfinal(6);

% compute equality constraints (final state boundary conditions)

f(2) = yfinal(1);

f(3) = yfinal(2);

f(4) = yfinal(3);

% transpose

f = f';

% no derivatives

g = [];

Contact us