I need some manipulation with ODE45. Please help me!

Question

Tsogerdene 2011-8-3

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/12954-i-need-some-manipulation-with-ode45-please-help-me

I just wondering that why this code working wrong. Please help me! I can't manipulate the ODE function.

options=odeset('RelTol',10^-4);
[T,Y]=ode23('turshilt23',[0 300],[0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5 0.5],options);
plot(T,Y);
function dy=turshilt23(t,y,flag)
%here is some code
if(k==0)
  y(1)=y(1)/2;
end
dy(1)=mu*y(1)*(1-(y(1)/m_max));          
%here is some code
end

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Arnaud Miege 2011-8-3

2
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/12954-i-need-some-manipulation-with-ode45-please-help-me#answer_17719

Please format your code.

Also, it seems that in your function turshilt23, you haven't defined what k, mu or m_max are. Finally you compute dy(1), but you call the ode solver with an initial condition vector of length 9. You need to compute dy(2), dy(3), ... dy(9) in your function as well.

Have a look at the documentation for the ode solvers, there are various examples you can inspire yourself from.

HTH,

Arnaud

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

Answer 2

Tsogerdene 2011-8-3

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/12954-i-need-some-manipulation-with-ode45-please-help-me#answer_17725

Thank you for your answer. I wrote all of code.

 function dy=turshilt23(t,y,flag)
        global Global_CycB1;
        global Global_CycB2;
        global k;
        global y_1;
        dy=zeros(9,1);
        % size of protein
        mu=0.005;
        m_max=10;
        dy(1)=mu*y(1)*(1-(y(1)/m_max));           
        if Global_CycB2 >=0.099 & Global_CycB2<0.1
            if Global_CycB2>Global_CycB1
                      y(1)=y(1)*(1/2);
            end
        end     
        %CycBT
        k1=0.04;
        k2_1=0.04;
        k2_2=1;
        k2_3=1;
        dy(2)=k1-(k2_1+k2_2*y(3)+k2_3*y(5))*y(2);
        Keq=1000;
        beta=y(2)+y(7)+Keq^(-1);
        Trimer=(2*y(2)*y(7))/(beta+(beta^2-4*y(2)*y(7))^(1/2));
        CycB=y(2)-Trimer;
        Global_CycB1=Global_CycB2;
        Global_CycB2=CycB;
        %Cdh1
        k3_1=1;
        k3_2=10;
        k4_1=2;
        k4=35;
        J3=0.04;
        J4=0.04; 
        dy(3)= (((k3_1+k3_2*y(5))*(1-y(3)))/(J3+1-y(3)))-(((k4*y(1)*CycB+k4_1*y(8))*y(3))/(J4+y(3))); 
        %Cdc20T
        k5_1=0.005;
        k5_2=0.2;
        k6=0.1;
        J5=0.3;
        n=4;
        dy(4)= k5_1+k5_2*((y(1)*CycB)^n/(J5^n+(y(1)*CycB)^n))-k6*y(4);
        %Cdc20A
        k6=0.1;
        k7=1;
        k8=0.5;
        J7=1*(10^(-3));
        J8=1*(10^(-3));
        Mad=1;
        dy(5)= ((k7*y(6)*(y(4)-y(5)))/(J7+y(4)-y(5)))-((k8*Mad*y(5))/(J8+y(5)))-k6*y(5);
        %IEP
        k9=0.1;
        k10=0.02;
        dy(6)= k9*y(1)*CycB*(1-y(6))-k10*y(6);
        %CKI
        k11=1;
        k12_1=0.2;
        k12_2=50;
        k12_3=100;
        dy(7) = k11-(k12_1+k12_2*y(8) +k12_3*y(1)*CycB)*y(7);
        %SK
        k13_1=0;
        k13_2=1;
        k14=1;
        dy(8)= k13_1+k13_2*y(9)-k14*y(8);
        %TF
        k15_1=1.5;
        k15_2=0.05;
        k16_1=1;
        k16_2=3;
        J15=0.01;
        J16=0.01;
        dy(9)= (((k15_1*y(1)+k15_2*y(8))*(1-y(9)))/(J15+1-y(9)))-(((k16_1+k16_2*y(1)*CycB)*y(9))/(J16+y(9)));

I defined all of them. But doesn't working. Please help me!

2 个评论
显示无隐藏无

Arnaud Miege 2011-8-3

It works for me. What error message do you get?

Some comments:

* Your function is a function of t, y and flag but t and flag are not used. Consider replacing those with ~

* Your function defines k and y_1 as global variables, but those aren't used in the function. You can remove those two lines.

The only changes I made were therefore:

function dy=turshilt23(~,y,~)

global Global_CycB1;

global Global_CycB2;

% global k;

% global y_1;

I used 0.1 for the other two global variables Global_CycB1 and Global_CycB2. The solver solved OK.

Walter Roberson 2011-8-3

Note: Using ~ in that way is not available in 2008 and earlier.

请先登录，再进行评论。

Answer 3

Tsogerdene 2011-8-8

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/12954-i-need-some-manipulation-with-ode45-please-help-me#answer_18176

Thank you for helping. I just asking that y(1)=y(1)*(1/2)? I would like to divide actual value of y(1) when some conditions are satisfied. I used y(1)=y(1)*(1/2), but nothing is changed.