If I have a sigmoid equation depend on time and I want to input to Runge Kutta Orde 4Th equation. So, the graph show the sigmoid from Runge_kutta . How to call the function?

2 views (last 30 days)
This is my sigmoid equation depend on time ( j is looping for time)
%input for time
t(1)=0;
dt=0.1; %time interval
t=0:dt:100; %time span
T=zeros(length(t)+1,1); %empty array for t
M1=zeros(length(t)+1,1); %empty array for M1
for j = 1:length(t)-1
T(j+1)=T(j)+dt;
M1(j+1)= M1(j)+1./(1+exp(-T(j))); %this is sigmoid equation
end
figure
plot (T,M1,'r','Linewidth',3)
xlabel('time')
ylabel('M1')
That equation I want to input into Runge kutte orde 4 and show the sigmoid graph from RK 4
Thank you

Accepted Answer

Sam Chak
Sam Chak on 29 May 2023
Edited: Sam Chak on 29 May 2023
I'm unsure if I understand the issue correctly. Do you expect a plot like this?
tspan = linspace(0, 100, 10001);
M0 = 0; % initial value
[t, M] = ode45(@odefcn, tspan, M0);
% Solution Plot
plot(t, M, 'linewidth', 1.5),
grid on, xlabel('t'), ylabel('M(t)')
% System
function Mdot = odefcn(t, M)
Mdot = 1/(1 + exp(- t));
end
  2 Comments
Sam Chak
Sam Chak on 29 May 2023
The integration method in your code looks like the forward Euler method and it produce similar result like the RK4.
t(1) = 0;
dt = 1e-2; % time interval
t = 0:dt:100; % time span
T = zeros(length(t) + 1, 1); % empty array for t
M1 = zeros(length(t) + 1, 1); % empty array for M1
for j = 1:length(t)-1
T(j+1) = T(j) + dt;
M1(j+1) = M1(j) + dt./(1 + exp(- T(j)));
end
figure
plot(t, M1(1:end-1), '--'), grid on
xlabel('time')
ylabel('M1')

Sign in to comment.

More Answers (1)

Sam Chak
Sam Chak on 30 May 2023
The RK4 algorithm by Cleve Moler, is available in File Exchange.
t = linspace(0, 100, 1001);
M0 = 0; % initial value
Mout = ode4(@odefcn, 0, 0.1, 100, M0);
% Solution Plot
plot(t, Mout, 'linewidth', 1.5),
grid on, xlabel('t'), ylabel('M(t)')
% System
function Mdot = odefcn(t, M)
Mdot = 1/(1 + exp(- t));
end
% Classical Runge-Kutta 4th-order ODE solver
function yout = ode4(F, t0, h, tfinal, y0)
y = y0;
yout = y;
for t = t0 : h : tfinal-h
s1 = F(t,y);
s2 = F(t+h/2, y+h*s1/2);
s3 = F(t+h/2, y+h*s2/2);
s4 = F(t+h, y+h*s3);
y = y + h*(s1 + 2*s2 + 2*s3 + s4)/6;
yout = [yout; y];
end
end

Categories

Find more on Programming 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!