Can ODE45 take in arrayed values?

1 view (last 30 days)
nashyshan
nashyshan on 2 Jun 2015
Commented: nashyshan on 2 Jun 2015
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);

Accepted Answer

Torsten
Torsten on 2 Jun 2015
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.

More Answers (1)

Torsten
Torsten on 2 Jun 2015
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 Comment
nashyshan
nashyshan on 2 Jun 2015
@Torsten How do I solve them simultaneously? Can you suggest?

Sign in to comment.

Tags

Community Treasure Hunt

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

Start Hunting!