MATLAB Answers

fzero error in fa value

21 views (last 30 days)
Matias Calzada Gianero
Matias Calzada Gianero on 6 Dec 2020
Answered: Alan Weiss on 7 Dec 2020
Hi Guys!
I´m trying to implement a code that is a function which 3 parameters as arguments. It´s a model for Comparmental model study. I get the following error:
Error in fzero (line 226)
fa = localFirstFcnEval(FunFcn,FunFcnIn,a,varargin{:});
And as fzero fails, that leads to the next error inherent to the function
Error in cortinez (line 46)
ke0 = fzero(@getKe0, ike0, options);
The only arguments that I had to put on the function are age, weight and a exponent. I will let the code for anyone that maybe can find the error. But It´s an internal error, so I´m a bit lost. Thanks in advance!
function [sys] = cortinez(Age,Weight,stepSize)
% Implements the model of Cortinez using table S2a for the pooled data
% Cortínez LI, Anderson BJ, Penna A, Olivares L, Muñoz HR, Holford NH, Struys MM, Sepulveda P.
% Influence of obesity on propofol pharmacokinetics: derivation of a pharmacokinetic model.
% Br J Anaesth. 2010;105:448-56.
% Standard volumes (Liters for 70 kg patient)
V1std = 4.48;
V2std = 21.2;
V3std = 237;
%Standard clearances (Liters/min for 70 kg patient)
CL1std = 1.92;
Q2std = 1.45;
Q3std = 0.86;
%Age correction factors
SLV2 = -0.0164;
SLQ2 = -0.0153;
wa = Weight/70;
wap = wa^0.75;
aa = Age - 50;
V1 = V1std * wa;
V2 = V2std * wa * exp(SLV2*aa);
V3 = V3std * wa;
volume = [V1 V2 V3];
CL1 = CL1std* wap;
CL2 = Q2std* wap * exp(SLQ2*aa);
CL3 = Q3std*wap;
clearance = [CL1 CL2 CL3];
% Solve for the value of ke0 which yields a time to peak effect of 1.4
% minutes. We search in the range of .1 to 2.0
tpeak = 1.6;
ike0 = [.1;2];
t=0:.0005:3;
u=zeros(size(t));
u(1)=1;
options = optimset('Display','off'); % Turn off Display
ke0 = fzero(@getKe0, ike0, options);
%disp(ke0);
function [error] = getKe0(tke0)
[sys] = mam2ss(clearance,volume,tke0);
[y,t] = lsim(sys,u,t);
s = lsiminfo(y,t,0); % final value is 0
error = tpeak - s.MaxTime;
end
sys = mam2ss(clearance,volume, ke0);
if (nargin==3)
sys=c2d(sys, stepSize,'zoh');
end
end
function [sys] = mam2ss(clearance,volume, ke0)
% returns a state space model from vectors of clearance
% and compartment volumes. All states are drug amount
% in the compartment; to get concentration, divide by compartment
% volume in the C matrix. Note that this code can handle an
% arbitrary number of compartments, but assumes the first compartment is
% the central compartment (mam is short for mamillary).
% The effect site compartment is represented as a compartment of negligible
% volume adjoined to the central compartment.
n = length(clearance);
if (nargin == 2)
k1 = clearance./volume(1);
k2 = clearance(2:n)./volume(2:n);
A = [-sum(k1) k2;transpose(k1(2:n)) -diag(k2)];
B = [1;zeros(n-1,1)]; % Direct addition of drug into the central compartment
C = [1/volume(1) zeros(1,n-1)]; % Observation of the central compartment concentration.
else
EFFECT_VOL_FACTOR=10000; % ratio of central compartment vol to effect
k1 = [clearance./volume(1) ke0/EFFECT_VOL_FACTOR];
k2 = [clearance(2:n)./volume(2:n) ke0];
A = [-sum(k1) k2;transpose(k1(2:n+1)) -diag(k2)];
B = [1;zeros(n,1)]; % Direct addition of drug into the central compartment
C = [zeros(1,n) EFFECT_VOL_FACTOR/volume(1)]; % Observation of the effect site compartment concentration.
end
C = C./1000; % All compartment volumes are in L, convert to ml
D=0;
sys = ss(A,B,C,D);
end

  0 Comments

Sign in to comment.

Answers (1)

Alan Weiss
Alan Weiss on 7 Dec 2020
Next time, please give the full error thrown by the function. As it is, I don't know what the error is.
I guess that the error is that your function @getKe0 does not differ in sign at the two endpoints ike0. But that is just a guess because I don't know what the error is.
Alan Weiss
MATLAB mathematical toolbox documentation

  0 Comments

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!