MATLAB Answers

No plot is being generated

4 ビュー (過去 30 日間)
Kyle Broder
Kyle Broder 2016 年 8 月 10 日
コメント済み: Kyle Broder 2016 年 8 月 10 日
I'm considering a system of differential equations given by dR/dt = U and dU/dt = (L-GM)/R^2. I'm trying to use RK4 to solve this system and generate some plots. My code is given below, however I'm not generating any plots. Is therefore something wrong with my RK4 code? Note that the code does run.
function [x,t] =nma_RK4_4( )
delT = 0.001; %time grid
t = linspace(0,5,round(5/delT)); %time
N = length(t);
x = zeros(2,N);
x(1,1) = 100000; %initial conditions
x(2,1) = 100000;
G = 6.67*power(10,-11);
M = 1.5*1.989*power(10,41);
L = x(1,1)*x(2,2);
theta = -x(2,1)/(power(x(1,1),2));
xvar = x(1,1)*cos(theta);
yvar = x(1,1)*sin(theta);
for i = 1:(N-1)
k1 = delT * f( t(i) , x(:,i));
k2 = delT * f( t(i)+1/2*delT , x(:,i)+1/2*k1);
k3 = delT * f( t(i)+1/2*delT , x(:,i)+1/2*k2);
k4 = delT * f( t(i)+delT , x(:,i)+k3);
x(:,i+1) = x(:,i)+1/6*(k1+2*k2+2*k3+k4);
end
function r=f(t,x)
r=[x(2)
(L-G*M)/(power(x(1),2))];
end
figure()
plot(xvar,yvar)
end

  0 件のコメント

サインインしてコメントする。

回答 (1 件)

Walter Roberson
Walter Roberson 2016 年 8 月 10 日
Your xvar and your yvar are scalars, so plot(xvar,yvar) is plotting only one point. The default when using plot() is to not use any markers, so the one point is not showing up. It would show up if you used plot(xvar, yvar, 'r*-')
Is there a reason you are calculating x but not using it to plot?

  1 件のコメント

Kyle Broder
Kyle Broder 2016 年 8 月 10 日
Thanks for your comment, you've made me realise an error in my code. After fixing the error, I'm still not getting a plot. Essentially, what I want to do is solve the above system of ODEs for R and then use that value obtained by using RK4 to then plot the result. I want to plot x against y however, and this is obtained from setting x = Rcos(theta) and y=Rsin(theta). As can be seen this is not working for me.
I don't want to plot just a point, I want something like this, https://arxiv.org/pdf/1503.05861.pdf, see the first graph on page 8. Note that the plot I've given as an example in that link does not look like what I have, given that it is a different system of ODEs, but the type of plot it what I want to generate.

サインインしてコメントする。

Community Treasure Hunt

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

Start Hunting!

Translated by