Can ODE45 take in arrayed values?

조회 수: 11 (최근 30일)
nashyshan
nashyshan 2015년 6월 2일
댓글: nashyshan 2015년 6월 2일
Hi, I am getting an error in my code. The error is " Error using odearguments ( line 93) @(T,IP)V4.*(C+(1-ALPHA).*K4)./(C+K4)-(IR*IP) returns a vector of length 10000, but the length of initial conditions vector is 1. The vector returned by @(T,IP)V4.*(C+(1-ALPHA).*K4)./(C+K4)-(IR*IP) and the initial conditions vector must have the same number of elements ". The error is generated from the last ODE45 function.
This is the code
clear all; clc;
%--------------------------------------------------------------------------
% The simulation is based on the model described by DeYoung and Keizer in
% the paper titled " A single-pool inositol 1,4,5-trisphosphate-receptor-
% based model for agonist-stimulated oscillations in Ca2+ concentration"
%--------------------------------------------------------------------------
%--------------------------------------------------------------------------
% Initial conditions
%--------------------------------------------------------------------------
Ca_ER=10*10^-6; Ca_cyto=1.7*10^-6;
Ir=1; alpha=.5; p_open3=0.15;
%--------------------------------------------------------------------------
% Constants in micromolar
%--------------------------------------------------------------------------
c0=4*10^-6; c1=.185;
v1=6; v2=.11; v3=.09*10^-6; v4=1.2;
k3=.1*10^-6; k4=1.1*10^-6;
%--------------------------------------------------------------------------
% Receptor Binding Constants in micromolar per second
%--------------------------------------------------------------------------
a1=400*10^-6; a2=0.2*10^-6; a3=400*10^-6; a4=0.2*10^-6; a5=20*10^-6;
%--------------------------------------------------------------------------
% Receptor Dissociation Constants in micromolar
%--------------------------------------------------------------------------
d1=0.13*10^-6; d2=1.049*10^-6; d3=.9434*10^-9; d4=.1445*10^-9;
d5=.08234*10^-9;
%--------------------------------------------------------------------------
% ODE describing Ca2+ concentrations in the cyctosol. Refer Ca2+
% oscillations
%--------------------------------------------------------------------------
dc=@(t,c) (c1.*(v1.*(p_open3)+v2).*(Ca_ER)-c)-v3.*((c).^2)./(c.^2+(k3).^2);
[t,c]=ode45(dc,linspace(0, 100, 10000),.15*10^-6);
plot(t,c);
%--------------------------------------------------------------------------
% Obtaining Ca_ER from the conservation condition. Refer Ca2+ oscillations
%--------------------------------------------------------------------------
Ca_ER=(c0-c)./c1;
figure(2);
plot(t,Ca_ER);
%--------------------------------------------------------------------------
%ODE describing IP3 production by Ca2+ feedback. Refer equation 4
%--------------------------------------------------------------------------
dIP3= @(t,ip) v4.*(c+(1-alpha).*k4)./(c+k4)-(Ir*ip);
[t,ip3]=ode45(dIP3,linspace(0, 100),.2*10^-6);
plot(t,ip3);

채택된 답변

Torsten
Torsten 2015년 6월 2일
dy=@(t,y) [(c1.*(v1.*(p_open3)+v2).*(Ca_ER)-y(1))-v3.*((y(1)).^2)./(y(1).^2+(k3).^2); v4.*(y(1)+(1-alpha).*k4)./(y(1)+k4)-(Ir*y(2))];
[t,y]=ode45(dy,linspace(0, 100, 10000),[.15*10^-6 .2*10^-6]);
plot(t,y(:,1),t,y(:,2));
Best wishes
Torsten.

추가 답변 (1개)

Torsten
Torsten 2015년 6월 2일
You should solve the two equations for c and ip simultaneously, not in two steps.
Then the actual c at each time instant is always available for the calculation of ip.
Best wishes
Torsten.
  댓글 수: 1
nashyshan
nashyshan 2015년 6월 2일
@Torsten How do I solve them simultaneously? Can you suggest?

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Numeric Solvers에 대해 자세히 알아보기

태그

Community Treasure Hunt

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

Start Hunting!

Translated by