| secant1(inibod,endbod,inid,init,dt2,tol,tol2,miter,inig)
|
function [alpha1, iniv1, error, kiter] = secant1(inibod,endbod,inid,init,dt2,tol,tol2,miter,inig)
%once gravity 2 has been run, and given the entry parameters (starting body, time ,
%target body) , it computes the necessary angle to impact the target
%body with a starting velocity equal to the escape velocity
%IT MAY NOT HAVE A SOLUTION FOR CERTAIN INITIAL CONDITIONS
%inibod: starting body
%endbod: target body
%inid: distance from the initial body from where the rocket is launched
%dt2: time step size
%init: launching time
%inidir: launching direction [cos(alpha; sin(alpha)]
%tol and tol2: tolerances for numerical iterative processes
%miter: maximum number of iterations iterative processes
close all
dbstop if error
global g mass ndime
if ndime == 3
disp('Not ready for 3 dimensional case');
stop
end
alpha0 = inig;
alpha1 = alpha0-0.01;
x0 = alpha0;
x1 = alpha1;
ve = sqrt(2*g*mass(inibod)/inid);
[mdist] = mindist(inibod,inid,dt2,init,[cos(x0);sin(x0)],ve,tol2,miter,Inf);
fx0 = mdist(endbod);
[mdist] = mindist(inibod,inid,dt2,init,[cos(x1);sin(x1)],ve,tol2,miter,Inf);
fx1 = mdist(endbod);
error = 1;
kiter = 0;
while error > tol && kiter < miter
kiter = kiter +1;
grad1 = (fx1-fx0)/(x1-x0);
%grad1 = grad_bo([fx1,fx0,fx_1],[x1';x0';x_1']);
xk = x1 - fx1/grad1; %avanza en la direccin del gradiente
[mdist] = mindist(inibod,inid,dt2,init,[cos(xk);sin(xk)],ve,tol2,miter,Inf);
fxk = mdist(endbod);
error = abs(fxk);
x0 = x1;
fx0 = fx1;
x1 = xk;
fx1 = fxk;
end
if x1 > 2*pi
x1 = rem(x1,2*pi);
elseif x1 < 0
x1 = rem(x1,2*pi)+2*pi;
end
alpha1 = x1;
iniv1 = ve;
[mdist] = mindist(inibod,inid,dt2,init,[cos(x1);sin(x1)],ve,tol2,miter,1);
|
|