MATLAB: Selecting last value of tspan ode45

for looplast valueode45tspan

I'm trying to extract results only for the last value of the tspan in ode45. However, I end up with 2 results for tfinal… I don't understand why this is happening.
Here's the code:
PS. The ode45 is ran in a for loop because I want to solve the set of diff eqs for 130 different values of ABP (which is a parameter in the eqs).
clear all 
clc
%%=========== ODE45 ==================   
ABP= linspace(40,170,131);
     for i=1:1:length(ABP)  %change in ABP at every loop 
         sol = ode45(@first,[0 100], [0 0 0 0 0.205 97.6 12 49.67],[],ABP(i)); %(function, timespan, initial condition for xq,xc,xm1,xm,Ca,P1,V_sa,P2)
     end
%%========= FUNCTION ==================
function dvdt = first(t,v,ABP)
xq= v(1); %variables
xc= v(2);
xm= v(3);
xm1= v(4);
Ca=v(5);
P1= v(6);
V_sa= v(7);
P2 = v(8);
   % -----Constants -----
  R_la= 0.4;
  R_sa_b= 5.03;
  R_sv= 1.32;
  R_lv= 0.56;
  P_v= 6;
  V_la=1;
  V_sa_b= 12;
  P_ic= 10;
  k_ven= 0.186;
  P_v1= -2.25;
  V_vn= 28; 
  tau_q= 20; 
  Pa_co2_b= 40; 
  tau_co2= 40; 
  tau1= 2; 
  tau2= 4;
  tau_g= 1; 
  C_a_p=2.87; 
  C_a_n= 0.164;
  g=0.2; 
  E=0.4; 
  K= 0.15; 
  V0= 0.02; 
  q_b=12.5;
  G_q=3; 
  Pa_co2=40; 
  Ca_b= 0.205;
  eps=1;
  u=0.5; 
  Pa_b=100;
  ka=3.68;
deltaCa_p=2.87;
deltaCa_n=0.164;
% =========== System of diff eq =========================================
dxq= (-xq+G_q*( ( ( (P1-P2)/(R_sv+0.5*(R_sa_b*(V_sa_b/V_sa)^2))) -q_b) /q_b) )/tau_q;
dxc=(-xc +0.3+3*tanh(Pa_co2/Pa_co2_b -1.1))/tau_co2;
dxm=(eps*u-tau2*xm1-xm)/tau1^2;
dxm1=xm1;
sum_xqxcxm=xm+xc-xq; %sum of feedback mechanisms
if t==100 % Last value of tspan!
    disp(sum_xqxcxm) %should get 14 displayed values (because i=14) but I get 28!
end
   if (ABP==40)||(ABP==41) ||(ABP==42) ||(ABP==43) ||(ABP==44) ||(ABP==45) ||(ABP==46) ||(ABP==47)||(ABP==48)||(ABP==49)||(ABP==50)||(ABP==51)||(ABP==52)
      delta=deltaCa_p; %because sum_xqxcxm >0 for ABP=40
elseif (ABP==53)||(ABP==54) %add other values later
     delta=deltaCa_n; %sum<0
 end
dCa=0.5*delta*(1-tanh(2*sum_xqxcxm/delta)^2);
dP1= 1/Ca * ((ABP-P1)/(R_la+0.5*(R_sa_b*(V_sa_b/V_sa)^2)) - (P1-P2)/(R_sv+0.5*(R_sa_b*(V_sa_b/V_sa)^2)) -dCa*(P1-P_ic));
dV_sa= dCa*(P1-P_ic) + Ca*dP1;
dP2=1/(1/(k_ven*(P2-P_ic-P_v1))) *((P1-P2)/(R_sv+0.5*(R_sa_b*(V_sa_b/V_sa)^2)) -(P2-P_v)/R_lv);
dvdt=[dxq;dxc;dxm;dxm1;dCa;dP1;dV_sa;dP2];
%%==== print to file====
Rsa= R_sa_b*(V_sa_b/V_sa)^2;
CBF= (P1-P2)/(Rsa*0.5 + R_sv);
CBF_ch= (CBF-q_b)/q_b *100;
if t==100 %save only final solution
    fileID=fopen('results.txt','a');
     fprintf(fileID,' %-5s %-2s %-5s %-5s %-5s %-5s %-5s %-5s %-5s %-5s\n','xq','xc','xm','xm1','Ca','P1','V_sa','P2', 'CBF','CBFchange'); %how to write this only as a header and not repeat every time?
fprintf(fileID,' %4.3f %4.3f %4.3f %4.3f %4.3f %4.3f %4.3f %4.3f %4.3f %4.3f\n\n',xq,xc,xm,xm1,Ca,P1,V_sa,P2,CBF,CBF_ch)
fclose(fileID);
end
end

Best Answer

Function "first" is called several times for the same t-values.
Better call the function "first" once more after ode45 has finished (i.e. after sol = ode45(...)). Or use OutputFcn.
Best wishes
Torsten.
Related Question