I had this code working previously and made a few adjustments and I can not work out how I should resolve this and get it working again, the code is in 2 separate files. The equation at the bottom of the first file (CsBEV) is what I added in and I need to plot it
First file is as follows
function diff=ode_sys_2equations(t,var)
y1=var(1);
y2=var(2);
y3=var(3);
y4=var(4);
% Parameters for DOX
Dose_DOX = (170/1000); % g/m^2/s check this
Tiv_DOX = 7200; % s
t_half_DOX = 1194.8; % s
kclr_DOX = log(2)/t_half_DOX; % 1/s
ke_DOX = 5.8E-4; % 1/s
Vdist_DOX = 0.00225; % m^3/kg
Mass = 85; % kg
P_DOX = 3.0E-6; % m/s
% Parameters for Bevacizumab
Dose_BV = (850/1000); % g/s for 7200s at 5mg/kg every 2 weeks
Tiv_BV = 7200; % s
t_half_BV = 1728000; % s (half life = 20 days for a male)
kclr_BV = log(2)/t_half_BV
ka_BV = 0; % 1/s *** NOT REQUIRED
Vdist_BV = 0.00329; % m3 or 4.7E-5 m3/kg for a 70kg male
Mass = 85; % kg
% Paramaters
alpha=-4.1E-6;
beta=1.22E-5;
gamma=-8.1E-6;
kak=5E-7;
r_tumour=0.03; % m
Vt=(4/3)*pi*r_tumour^3; % m^3
a=0.35; % Dimensionless (taken from table 1 - paramaters for brain tumour/normal tissue)
Ve=Vt*a; % m^3
s_v=20000; % m^-1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if 0<=t & t<=Tiv_DOX
Rin_DOX=1; % mg/kg *** CHECK THIS ***
elseif 86400<=t & t<=93600
Rin_DOX=1;
elseif 172800<=t & t<=180000
Rin_DOX=1
elseif 259200<=t & t<=266400
Rin_DOX=1
elseif 345600<=t & t<=352800
Rin_DOX=1
else
Rin_DOX=0
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if 0<=t && t<=Tiv_BV
Rin_BV=1
elseif 1123200<=t && t<1130400
Rin_BV=1
else
Rin_BV=0
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Concn TMZ BLOOD TEST
diff(1,1)=(-kclr_DOX*y1)+(Rin_DOX*Dose_DOX); % y1 - Zhan *********
% Concn BV BLOOD TEST - ZHAN - y2
diff(2,1)=(-kclr_BV*y2)+(Rin_BV*Dose_BV); % y2 - Zhan ************
% Normalised Bevacizumab Concentration - Check the units for this equation
% or the equation itself and then everything should work ok.
% Anti-Angiogenic Effect
diff(3,1)=y3*(alpha+beta*y3+gamma*y3^2)-y3*kak*CsBEV; % Skye ********
% Concn of TMZ in Tumour
diff(4,1)=((P_DOX*y4*s_v*Vt)*(y1-y4)-(ke_DOX*y4*Ve))/Ve % Skye ***
end
CsBEV = y2/(Dose_BV/Vdist_BV);
Second FILE
clc
% Working with ODE 1, 2, & 3 (Concn Blood, Cancer Cell Effect &
% Endothelial Cell Effect)
range=[0:86400*30]; % 864000(s)=10days
% Initial conditions of system [Ctmz]
ICs=[0,0,1,0];
% Solve ODEs using ODE45
[tsol,varsol]=ode45(@WORKING_ODE_SYS, range, ICs);
% plot of Concentration of TMZ in Blood against time (g/m2/s)
figure(1)
plot(tsol/86400,varsol(:,1))
title('Concentration of DOX in Blood vs. Time')
xlabel('Time (hr)')
ylabel('Concentration of DOX in Blood (g/m^2)')
% plot of Concentration of BV in Blood against time (g)
figure(2)
plot(tsol/86400,varsol(:,2))
title('Concentration of Bevacizumab vs. Time')
xlabel('Concentration of Bevacizumab in Blood (g)')
ylabel('Time (hr)')
% plot of Anti-Angiogenic Effect
figure(3)
plot(tsol/86400,varsol(:,3))
title('Plot of Anti-Angiogenic Effect vs. Time')
xlabel('Time (hr)')
ylabel('Anti-Angiogenic Effect (%)')
% plot of Concnetration of TMZ in Tumour against time
figure(4)
plot(tsol/86400,varsol(:,4))
title('Concentration of DOX in Tumour vs. Time')
xlabel('Time (hr)')
ylabel('Concentration of DOX (g/m^2)')
figure(5)
plot(tsol,varsol(:,1))
hold on
plot(tsol,varsol(:,4))
title('Concentration of DOX in Blood & Tumour')
xlabel('Time (hr)')
ylabel('Concentration of DOX (g/m^2)')
hold off
