How do I solve these differential equations using a while loop?

[Roger Stafford]

If you want to allow delta t to approach zero as a limit, you can solve these equations using one of the ode functions. The first equation, for example, would have the form:

dU/dt = k1-k2*X./((X.^2+Y.^2+Z.^2).^(3/2))

On the other hand if you wish to solve them using delta t as a fixed nonzero value, then do so with a for-loop to provide the iteration, not a while-loop. Just carry out the operations you have given here within the for-loop at each step going from the n-th values to the n+1-st values.

Answer 1

[Abraham Boayue]

i = 1;
while i <= n-2
      i = i +1;
      % write all your code
      % here. This will produce 
      % the same results as the 
      % for loop.
end

[Abraham Boayue]

You change the for loop to the while loop above, it does the same operation as the for loop.

Answer 2

[Abraham Boayue]

clear variables
close all

% Define parameters
dt = dt;
t = t0:dt:tf;
n = length(t);
m  = m;
thx = thx; 
thy = thy;
thz = thz;
G   =   G;
Me  =  Me;

% Initializations       Initial conditions     Boundary conditions
u = zeros(1,n);          u(1) =  u0;           u(n) = un;
v = u;                   v(1) =  v0;           v (n)= vn;
w = v;                   w(1) =  w0;           w(n) = wn;
x = w;                   x(1) =  x0;           x(n) = xn;
y = x;                   y(1) =  y0;           y(n) = yn;
z = y;                   z(1) =  z0;           z(n) = zn;

for i = 2: n-1

     u(i+1) = u(i) + (thx/m - G*Me*(x(i)/(x(i)^2 +y(i)^2 +z(i)^2)^(3/2)))*dt;
     v(i+1) = v(i) + (thy/m - G*Me*(y(i)/(x(i)^2 +y(i)^2 +z(i)^2)^(3/2)))*dt;
     w(i+1) = w(i) + (thz/m - G*Me*(z(i)/(x(i)^2 +y(i)^2 +z(i)^2)^(3/2)))*dt;

     x(i)  = x(i)  + u(i+1)*dt;
     y(i)  = y(i)  + v(i+1)*dt;
     z(i)  = z(i)  + w(i+1)*dt;
  end

figure;
plot(t,u,'linewidth',2); 
hold on
plot(t,v,'linewidth',2);
plot(t,w,'linewidth',2);
plot(t,x,'linewidth',2);
plot(t,y,'linewidth',2);
plot(t,z,'linewidth',2);
a = ylabel('Pressure');
set(a,'Fontsize',14);
a = xlabel('x');
set(a,'Fontsize',14);
a=title(['Solution to system of ode - dt = ' num2str(dt)]);
legend('u', 'v','w','x','y','z')
xlim([0 1]);
set(a,'Fontsize',16);
grid;

[Abraham Boayue]

Hey Christopher Here is a model solution to your equations. I could not test the code since you did not provide any parameter. The coding was straight forward using a for loop as Roger suggested. I hope this will help you get started with your coding. Best of lucks!

[Christopher Maraj]

Thank you very much for your help. However this is for a school project and we are required to use a while loop. While what you provided have provided a lot of clarification, is there any way to adapt it using a while loop?

[Abraham Boayue]

Hey Chris, what are these equations by the way? Do you have access to the actual differential equations? Can you post the values of the parameters as well as the initial and boundary conditions? I would like to see what they represent by running the code.

[Christopher Maraj]

yes sure, here they are:

[Christopher Maraj]

The equations are for a project I'm working on about modelling orbits for 6 satellites.