Adding a piecewise time dependent term in system of differential equation

8 ビュー (過去 30 日間)
Problem Statement: I have a system of differential equation and out of which few has time dependent coefficients and in my case i have piecewise time dependence. I am unable to add the time dependence and solve the equation.
Simplified Differential Equation:
My time span runs from 0 to 360 days and my aim is to have certain nonzero value for for half a day and zero for other half.
My unsuccessful attempt:
Function definition: MWE_fn.m
function rk1 = MWE_fn(t,y)
r2=0.5;b2=1;d_A=0.05; c=0.5;
rk1(1)=r2*y(1)*(1-b2*y(1))-c*y(1)*y(2);
rk1(2)= A0(t)- d_A*y(2);
rk1=rk1(:);
end
function fa=A0(t)
if rem(t,1)==0
fa=0.04
else
fa=0
end
end
Function Body: MWE_body.m
timerange= 0:0.5:360;
IC= [0.1,0];%initial conditions
[t,y] =ode45(@(t,y) MWE_fn(t,y),timerange, IC);
plot(t,y(:,1),'b');
hold on
plot(t,y(:,2),'r');
I was aiming to write some kind of conditional statement for the function to take those values but that didn't go correct. The thing i had in mind was that for half a day we have time at integer values and for other half i have time at 0.5 , so i used a rem function to check that, but i was wrong and it didn't work.
Anyhelp would be appreciated.
Thank You.

採用された回答

Thiago Henrique Gomes Lobato
Thiago Henrique Gomes Lobato 2019 年 10 月 27 日
編集済み: Thiago Henrique Gomes Lobato 2019 年 10 月 27 日
The conditional statement idea is right, the problem was that you considered that the ode45 would evaluate your function only at the time steps you gave, which is not right. To perform the integration many other time intervals are needed, so you just have to adjust your conditional:
function fa=A0(t)
if t>round(t)
fa=0.04;
else
fa=0;
end
end
In this case, for every value x.y where y<5, fa = 0.04, which represents basically all your half days
  16 件のコメント
Muhammad Sinan
Muhammad Sinan 2021 年 8 月 13 日
for example, the ode is dxdt = -x; and
if 1<t & t<=1.6
x0 = 0;
elseif 1.6<t & t<=3.2
.
.
end
then how will this condition defined?
Walter Roberson
Walter Roberson 2021 年 8 月 13 日
Using a modified form to illustrate holding over values:
x01_1 = 0;
x02_1 = .1;
x01_2 = pi/2;
x01 = [x01_1, x02_1];
[t1, x1] = ode45(@(t,x) [x(2)^2-x(1); x(2)-1/20], [1 1.6], x01);
x02 = x1(end,:); x02(1) = x02_1;
[t2, x2] = ode45(@(t,x) [x(2)^2-x(1); x(2)-1/20], [1.6 3.2], x02);
t = [t1;t2];
x = [x1; x2];
plot(t, x); legend({'x(1)', 'x(2)'});

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

その他の回答 (1 件)

Walter Roberson
Walter Roberson 2019 年 10 月 27 日
You should examine the odeball example to see how to set up events to terminate integration every time there is a change that is not continuous differentiable at least twice. ode45 and related are for continuous systems only.
Or in your case where the determination can be done based just on time, you can work it without events by looping over time intervals.
When you use conditional statements with ode45 and kin then you need to arrange so that the condition is always true or always false within the boundaries of the current integration.
  1 件のコメント
Vira Roy
Vira Roy 2019 年 10 月 28 日
Well, very special thanks to Walter and Thiago for helping me out.

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

カテゴリ

Help Center および File ExchangeLoops and Conditional Statements についてさらに検索

製品


リリース

R2016a

Community Treasure Hunt

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

Start Hunting!

Translated by