yankeefan11
yankeefan11

Reputation: 485

Numerical solution of coupled dif eq

So I have a set of 3 differentialequations I want to solve. They can be seen in my code. My problem, is that I want to combine these codes so that I can have a for loop with respect to R (as will be shown.

What I have:

T2 = 1;
[T,Y] = ode45(@ball, [0 5*T2] ,[0 0 -10]);
figure
plot(T,Y(:,1),'-r',T,Y(:,2),'-g',T,Y(:,3),'-b')
legend('x(t)','y(t)','z(t)')
xlabel('Time (in units of T2)')
title(['Plot for RT2 = ',num2str(R)])

Where the @ball is

`function dr = ball(t,b)

T2 = 1;
T1 = T2/2;
d  = 0;
R  = 0.2;

dr = zeros(3,1);
dr(1) = (-1/T2)*b(1)-d*b(2);
dr(2) = (-1/T2)*b(2) + d*b(1) + R*b(3);
dr(3) = (-1/T1)*b(3) - R*b(2) ;

end`

What I want is a single program that will do this, but allow me to include a for loop so I can vary R and make a couple subplots. IS this possible?

Upvotes: 0

Views: 36

Answers (2)

Nras
Nras

Reputation: 4311

Without the usage of an anonymous function (which is a decent way to get your result done), you can also pass the argument directly in the ode45-call. After the initial condition the next argument are the options, which can be left empty. After the options, additional Parameters can be submitted:

function main

T2 = 1;
opt = []; % // no further options
R = 0.2; % // the parameter R to give to the function ball
[T,Y] = ode45(@ball, [0 5*T2] ,[0 0 -10], opt, R); %% // added opt and R as parameter
figure
plot(T,Y(:,1),'-r',T,Y(:,2),'-g',T,Y(:,3),'-b')
legend('x(t)','y(t)','z(t)')
xlabel('Time (in units of T2)')
title(['Plot for RT2 = ',num2str(R)])

end


function dr = ball(t,b, R)

T2 = 1;
T1 = T2/2;
d  = 0;
% R  = 0.2; % // not needed anymore

dr = zeros(3,1);
dr(1) = (-1/T2)*b(1)-d*b(2);
dr(2) = (-1/T2)*b(2) + d*b(1) + R*b(3);
dr(3) = (-1/T1)*b(3) - R*b(2) ;

end

Upvotes: 0

David
David

Reputation: 8459

You can use an anonymous function for this.

Change ball.m to remove the hard-coded R and replace it with an input argument:

function dr = ball(t,b,R)

T2 = 1;
T1 = T2/2;
d  = 0;
%// etc.

and then replace your ode45 call with this:

R=0.4;
[T,Y] = ode45(@(t,b) ball(t,b,R), [0 5*T2] ,[0 0 -10]);

where @(t,b) ball(t,b,R) is a function with inputs t and b that calls ball.m with the value of R specified on the previous line. So you can construct the for loop as follows:

for R=0.2:.02:1 %// or whatever range you want
    [T,Y] = ode45(@(t,b) ball(t,b,R), [0 5*T2] ,[0 0 -10]);
    %// etc.
end

Upvotes: 1

Related Questions