'Undefined Function or Variable' Error when Function is Defined
显示 更早的评论
My model contains two functions, cp(t) and ct(t). In my original (working) code (see below), cp(t) and ct(t) are each defined by one functional form. Now, I would like to study the behavior of the dependent variables for different functional forms of cp(t) and ct(t). How can I do this?
My Attempt: I created two matrices in the ODE solver; in the first matrix, each column represents a function for cp(t), and in the second matrix, each column represents a function for ct(t). When I run the ODE solver, I get the error "Undefined Function or Variable" for cp(t).
Thank you!
-----------------------------------------------
ORIGINAL, WORKING CODE
ODE File:
function dy = short_ODE(t,y)
b = 21.5;
d = 0.167;
qe = 0.0542;
q1 = 0.000257;
q2 = 0.0008;
dy(1,1) = 10*b*y(9)*(cp(t)/(cp(t)+ct(t)))-ne(t)*y(1)-qe*y(1);
dy(2,1) = ne(t)*y(1)-n1(t)*y(2)-f1(t)*y(2)-q1*(y(2)^2)/cp(t);
dy(3,1) = n1(t)*y(2)-n2(t)*y(3)-f2(t)*y(3)-q2*(y(3)^2)/cp(t);
dy(4,1) = n2(t)*y(3)-np(t)*y(4)-fp(t)*y(4);
dy(5,1) = np(t)*y(4)/2-(1+d)*y(5)+y(9);
dy(6,1) = -(1+d)*y(6)+y(5);
dy(7,1) = -(1+d)*y(7)+y(6);
dy(8,1) = -(1+d)*y(8)+y(7);
dy(9,1) = -(1+d)*y(9)+y(8);
function y = T(t)
T0 = 21.5;
y = T0 + 0.6346*cos(t*2*pi/365) + 0.7731*sin(t*2*pi/365);
end
% Nested cp(t)
function y = cp(t)
y = 135000-35000*(1.446*cos(t*2*pi/365)+0.7109*sin(t*2*pi/365)...
+1.347*cos(2*t*2*pi/365)+1.408*sin(2*t*2*pi/365)...
-0.9942*cos(3*t*2*pi/365)-0.2396*sin(3*t*2*pi/365));
end
% Nested ct(t)
function y = ct(t)
epsilon = 1e-10;
y = 0*t + epsilon;
end
% Nested ne
function y = ne(t)
y = 0.693/(1.10316852+12.55262944/(1+(T(t)/15.21894679)^6.468691289));
end
% Nested n1
function y = n1(t)
y = 0.693/(2.711278535+18.35202441/(1+(T(t)/15.61912844)^5.845513817));
end
% Nested f1
function y = f1(t)
y = 0.666733-0.056977*T(t)+0.00125224*T(t)^2;
end
% Nested n2
function y = n2(t)
y = 0.693/(4.244561306+136.3457233/(1+(T(t)/10.3599642)^4.783523801));
end
% Nested f2
function y = f2(t)
y = 0.666733-0.056977*T(t)+0.00125224*T(t)^2;
end
% Nested np
function y = np(t)
y = 0.693/(8.583819474+18.15593517/(1+(T(t)/21.5533767)^7.852462638));
end
% Nested fp
function y = fp(t)
y = 0.666733-0.056977*T(t)+0.00125224*T(t)^2;
end
end
ODE Solver:
clear all;
tf = 1000;
dt = 1;
t = 0:dt:tf;
y0 = [3954200.8835438923;2310300.5386356956;2007000.7121853685;
103780.6956163843;1289400.7821495022/100/5;
1289400.7821495022/100/5;1289400.7821495022/100/5;
1289400.7821495022/100/5;1289400.7821495022/100/5];
[t,y] = ode15s(@short_ODE,t,y0);
A = y(:,5)+y(:,6)+y(:,7)+y(:,8)+y(:,9);
plot(t,A,'b');
xlabel('Time (in days)');
ylabel('Population');
h_xlabel = get(gca,'XLabel');
set(h_xlabel,'FontSize',14);
h_ylabel = get(gca,'YLabel');
set(h_ylabel,'FontSize',14);
set(gca,'fontName','Helvetica');
MODIFIED CODE WITH ERROR
ODE Solver:
clear;
tf = 1000;
dt = 1;
t = 1:dt:tf;
y0 = [3954200.8835438923;2310300.5386356956;2007000.7121853685;
103780.6956163843;1289400.7821495022/100/5;
1289400.7821495022/100/5;1289400.7821495022/100/5;
1289400.7821495022/100/5;1289400.7821495022/100/5];
m = min(c(t));
d = m/10;
i = 0;
for p = 0:d:m % This is where I created matrices cp(:,i) and ct(:,i)
i = i + 1;
a(1:length(c(t))) = p;
cp(:,i) = transpose(a); %#ok<*SAGROW>
ct(:,i) = transpose(c(t)) - transpose(a);
end
for j = 1:11
calle1 = @(t,y) e1_ODE(t,y,j);
[t,y] = ode15s(calle1,t,y0);
A(:,j) = y(:,5)+y(:,6)+y(:,7)+y(:,8)+y(:,9);
plot(t,A,'b');
xlabel('Time (in days)');
ylabel('Population');
h_xlabel = get(gca,'XLabel');
set(h_xlabel,'FontSize',14);
h_ylabel = get(gca,'YLabel');
set(h_ylabel,'FontSize',14);
set(gca,'fontName','Helvetica');
hold on
end
% Nested c(t)
function y = c(t)
y = 135000-35000*(1.446*cos(t*2*pi/365)+0.7109*sin(t*2*pi/365)...
+1.347*cos(2*t*2*pi/365)+1.408*sin(2*t*2*pi/365)...
-0.9942*cos(3*t*2*pi/365)-0.2396*sin(3*t*2*pi/365));
end
ODE File:
function dy = e1_ODE(t,y,j) % Here, I have index as an input, which wasn't needed
% for the working code
b = 21.5;
d = 0.167;
qe = 0.0542;
q1 = 0.000257;
q2 = 0.0008;
epsilon = 1e-10;
dy(1,1) = 10*b*y(9)*(cp(t,j)/(cp(t,j)+ct(t,j)))-ne(t)*y(1)-qe*y(1);
dy(2,1) = ne(t)*y(1)-n1(t)*y(2)-f1(t)*y(2)-q1*(y(2)^2)/cp(t,j)+epsilon;
dy(3,1) = n1(t)*y(2)-n2(t)*y(3)-f2(t)*y(3)-q2*(y(3)^2)/cp(t,j)+epsilon;
dy(4,1) = n2(t)*y(3)-np(t)*y(4)-fp(t)*y(4);
dy(5,1) = np(t)*y(4)/2-(1+d)*y(5)+y(9);
dy(6,1) = -(1+d)*y(6)+y(5);
dy(7,1) = -(1+d)*y(7)+y(6);
dy(8,1) = -(1+d)*y(8)+y(7);
dy(9,1) = -(1+d)*y(9)+y(8);
function y = T(t)
T0 = 21.5;
y = T0 + 0.6346*cos(t*2*pi/365) + 0.7731*sin(t*2*pi/365);
end
% Nested ne
function y = ne(t)
y = 0.693/(1.10316852+12.55262944/(1+(T(t)/15.21894679)^6.468691289));
end
% Nested n1
function y = n1(t)
y = 0.693/(2.711278535+18.35202441/(1+(T(t)/15.61912844)^5.845513817));
end
% Nested f1
function y = f1(t)
y = 0.666733-0.056977*T(t)+0.00125224*T(t)^2;
end
% Nested n2
function y = n2(t)
y = 0.693/(4.244561306+136.3457233/(1+(T(t)/10.3599642)^4.783523801));
end
% Nested f2
function y = f2(t)
y = 0.666733-0.056977*T(t)+0.00125224*T(t)^2;
end
% Nested np
function y = np(t)
y = 0.693/(8.583819474+18.15593517/(1+(T(t)/21.5533767)^7.852462638));
end
% Nested fp
function y = fp(t)
y = 0.666733-0.056977*T(t)+0.00125224*T(t)^2;
end
end
采纳的回答
更多回答(0 个)
类别
在 帮助中心 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!