Perturbing Values and determing an SFO
2 次查看(过去 30 天)
显示 更早的评论
I am using the following code
K1 = 0.005;
K2 = 0.005;
K3 = 0.005;
K4 = 0.005;
Kc = 0.5;
Kd = 0.02;
kdd = 0.01;
V2 = 1.5;
V4 = 0.5;
vd = 0.25;
VM1 = 3;
VM3 = 1.;
Synth = 0.025;
% Original Equations:
% PhsKnsKnt = [Synth-vd*X*(C/(Kd+C)) - kdd*C; M*VM3*((1-X)/(K3+(1-X))) - V4*
(X/(K4+X)); VM1*(C/(Kc+C))*((1-M)/(K1+(1-M))) - V2*(M/(K2+M))]; % Anonymous Function with substitutions:
PhsKnsKnt = @(T, P) [Synth - vd*P(2)*(P(1)/(Kd+P(1))) - kdd*P(1); P(3)*VM3*((1- P(2))/(K3+(1-P(2)))) - V4*(P(2)/(K4+P(2))); VM1*(P(1)/(Kc+P(1)))*((1-P(3))/(K1+(1- P(3)))) - V2*(P(3)/(K2+P(3)))];
[T P] = ode45(PhsKnsKnt, [0 100], [0.01; 0.01; 0.01]);
C = P(:,1);
X = P(:,2);
M = P(:,3);
figure(1)
plot(T, P)
legend('[C] (µ\itM\rm)', '[X] (µ\itM\rm)', '[M] (µ\itM\rm)', 'Location','NorthWest')
xlabel('Time (s)')
ylabel('Concentration')
grid
I am trying to perturb each of the values from K1- VM3 by 5%. For example for kd i tried to following : Kd = 0.02 -(0.02*0.05); I have changed the other constants accordingly by 5%. My issue is that my "C" value when running this code has yet to change even when i perturb different values. I am trying to figure out SFO to rank variables using the following:
%change in C / % change of the constant (being 5%).
0 个评论
采纳的回答
Star Strider
2014-3-3
编辑:Star Strider
2014-3-3
Actually, it does change. It just doesn’t change very much.
Run this:
Tvct = linspace(0,33);
Kdv = [0.95 1 1.05]*Kd;
for k1 = 1:3
Kd = Kdv(k1);
PhsKnsKnt = @(T, P) [Synth - vd*P(2)*(P(1)/(Kd+P(1))) - kdd*P(1); P(3)*VM3*((1-P(2))/(K3+(1-P(2)))) - V4*(P(2)/(K4+P(2))); VM1*(P(1)/(Kc+P(1)))*((1-P(3))/(K1+(1-P(3)))) - V2*(P(3)/(K2+P(3)))];
[T P] = ode45(PhsKnsKnt, [Tvct], [0.01; 0.01; 0.01]);
CXM(:,:,k1) = P;
end
Csq = squeeze(CXM(:,1,:));
Csqr = [Csq(:,1)./Csq(:,2) Csq(:,3)./Csq(:,2)];
Csqd = [Csq(:,1)-Csq(:,2) Csq(:,3)-Csq(:,2)];
I put in Tvct to be sure the function was evaluated at the same time points (since the ode solvers are adaptive). I created Csq, a matrix of the three runs for C(t) with values for Kd varying from 0.95 to 1.05, and then calculated the ratios of them in Csqr. Csqd are the differences, since that may be what you want.
(I put up a transient earlier post about sensitivity functions. That works, but only on actual objective functions. I forgot for a minute that this is a DE. Doesn’t apply to DEs.)
2 个评论
Star Strider
2014-3-3
编辑:Star Strider
2014-3-3
I restricted the time to one oscillation here because it’s easier to see the details that way. Kdv of 0.95*Kd = Kd-(Kd*0.5), similarly for 1.05.
The initial conditions of the DE in your original post were all 0.01. They have to be specified for each variable in the MATLAB ode solvers, thus the [3 x 1] vector.
P.S. — Change the solver from ode45 to ode15s. It’s a bit better at these sorts of DEs.
P.P.S. — Plot this to see the effect:
figure(2)
plot(T,Csqr)
legend('[C] (µ\itM\rm) K_d = 0.95xK_d', '[C] (µ\itM\rm) K_d = 1.05xK_d', 'Location','NorthWest')
grid
更多回答(1 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!