Can someone explain why my graph is not generating any lines. I copied my code from this video https://ww​w.youtube.​com/watch?​v=vaMiCkrR​C6Q

15 views (last 30 days)
clear
clf
%All Constants Used
G = 6.67408e-11;
M = 5.9722e24;
m = 1000;
%Changing Variables
r = [6.67e6, 0 , 0]; %Radius
v = [0, sqrt(G*M/norm(r(1,:),2)),0 ]; %Hohmann Equation
t = 0.0; %Start time
dt = 1; %Steps/Change in Calculation
DeltaV = 2420; %Change in Velocity
TtoDeltaV = 5431;
%First Orbit
rplot = r;
while t <= TtoDeltaV
rlen = norm(r,2); %Radius Magnitude
accel = -G*M/rlen^2;
nextv = v + accel*dt*r/rlen;
nextr = r + dt*v;
t = t + dt;
rplot = [rplot, nextr];
v = nextv;
r = nextr;
end
hold on;
Orbit1 = plot(rplot(:,1), rplot(:,2),'b','LineWidth',1);
%Second Orbit
rplot = r;
theta = 0;
rp = r;
rplen = norm(rp,2);
nextv = v + DeltaV*[0 1 0];
v = nextv;
while theta <= 3.14159 %Position Vector
rlen = norm(r,2);
accel = -G*M/rlen^2;
nextv = v + dt*accel*r/rlen;
nextr = r + dt*v;
t = t + dt;
rplot = [rplot, nextr];
v = nextv;
r = nextr;
tempDot = dot(rp,r); %Dot product
costheta = tempDot/(rlen*rplen);
theta = acos(costheta);
end
disp("Sat reaches forbits at t = " + t + "s")
Sat reaches forbits at t = 24323s
Orbit2 = plot(rplot(:,1), rplot(:,2),'r','LineWidth',1);
%Orbit 3
rplot = r;
tstart = t; %Start or Orbit
OrbitP = (2*pi)*sqrt(norm(r,2).^3)/(G*M);
OrbitV = sqrt((G*M)/norm(r,2));
DeltaV = OrbitV - norm(v,2);
nextv = v + DeltaV*[0 -1 0]; %-1 to add the direction of going down
v = nextv;
while t < (tstart +OrbitP);
rlen = norm(r,2);
accel = -G*M/rlen^2;
nextv = v + dt*accel*r/rlen;
nextr = r + dt*v;
t = t + dt;
rplot = [rplot, nextr];
v = nextv;
r = nextr;
end
Orbit3 = plot(rplot(:,1),rplot(:,2),'k','LineWidth',1);
title('Hohmanns Transfer Orbit');

Answers (1)

Walter Roberson
Walter Roberson on 4 Dec 2023
clear
clf
%All Constants Used
G = 6.67408e-11;
M = 5.9722e24;
m = 1000;
%Changing Variables
r = [6.67e6, 0 , 0]; %Radius
v = [0, sqrt(G*M/norm(r(1,:),2)),0 ]; %Hohmann Equation
t = 0.0; %Start time
dt = 1; %Steps/Change in Calculation
DeltaV = 2420; %Change in Velocity
TtoDeltaV = 5431;
%First Orbit
rplot = r;
while t <= TtoDeltaV
rlen = norm(r,2); %Radius Magnitude
accel = -G*M/rlen^2;
nextv = v + accel*dt*r/rlen;
nextr = r + dt*v;
t = t + dt;
rplot = [rplot, nextr];
v = nextv;
r = nextr;
end
hold on;
Orbit1 = plot(rplot(:,1), rplot(:,2),'b','LineWidth',1);
%Second Orbit
rplot = r;
theta = 0;
rp = r;
rplen = norm(rp,2);
nextv = v + DeltaV*[0 1 0];
v = nextv;
while theta <= 3.14159 %Position Vector
rlen = norm(r,2);
accel = -G*M/rlen^2;
nextv = v + dt*accel*r/rlen;
nextr = r + dt*v;
t = t + dt;
rplot = [rplot, nextr];
v = nextv;
r = nextr;
tempDot = dot(rp,r); %Dot product
costheta = tempDot/(rlen*rplen);
theta = acos(costheta);
end
disp("Sat reaches forbits at t = " + t + "s")
Sat reaches forbits at t = 24323s
Orbit2 = plot(rplot(:,1), rplot(:,2),'r','LineWidth',1);
%Orbit 3
rplot = r;
tstart = t; %Start or Orbit
OrbitP = (2*pi)*sqrt(norm(r,2).^3)/(G*M);
OrbitV = sqrt((G*M)/norm(r,2));
DeltaV = OrbitV - norm(v,2);
nextv = v + DeltaV*[0 -1 0]; %-1 to add the direction of going down
v = nextv;
t
t = 24323
tstart
tstart = 24323
OrbitP
OrbitP = 0.0045
dt
dt = 1
OrbitP is less than dt. When you add t = t + dt then the result is going to be more than tstart + OrbitP
while t < (tstart +OrbitP);
rlen = norm(r,2);
accel = -G*M/rlen^2;
nextv = v + dt*accel*r/rlen;
nextr = r + dt*v;
t = t + dt;
rplot = [rplot, nextr];
v = nextv;
r = nextr;
end
so you are only going to have one point to plot.
You also have the problem that v is a 1 x 3 vector, so nextr will be a 1 x 3 vector, so your [rplot, nextr] is horzontal concatenating together 1 x 3 vectors each time -- but you expect to extract columns out of the result.
Orbit3 = plot(rplot(:,1),rplot(:,2),'k','LineWidth',1);
title('Hohmanns Transfer Orbit');

Categories

Find more on CubeSat and Satellites in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!