Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.136868e-13) at time t.

Question

Sveta am 3 Dez. 2023

0
Verknüpfen

Direkter Link zu dieser Frage

https://de.mathworks.com/matlabcentral/answers/2055674-unable-to-meet-integration-tolerances-without-reducing-the-step-size-below-the-smallest-value-allow

Kommentiert: Sveta am 5 Dez. 2023

for tf=1:1:8

tspan = t0:step:tf*dose_length;

lt=length(tspan);

sx=size(xx1,1);

x0=xx1(sx,:);

% t0=t0+dose_length;

options = odeset('RelTol',10,'AbsTol',[10 10 10 10 ]);

[tt,xx1]= ode15s (@odebladder_exp,tspan,x0,options,muB,muE,alpha,pA,pT,pR,pK,pD,gRate,0, 0,K,g_T);

t0=t0+dose_length;

sx=size(xx1,1);

x1((tf-1)*lt+1:lt*tf,:)=xx1;

t1((tf-1)*lt+1:tf*lt)=tt;

end

5 Kommentare
3 ältere Kommentare anzeigen3 ältere Kommentare ausblenden

Sveta am 4 Dez. 2023

Bearbeitet: Torsten am 4 Dez. 2023

In MATLAB Online öffnen

thanks a lot for your help!

the program works when I change only one parameter in ParamBCG file:alpha = 0.7*10^(-4) => alpha = 3.7*10^(-4)

Wait for your answer.

The best,

Svetlana.

% model for Bladder :

echo off

close all

clear

% 1 - bacilla, 2 - immune respfonse, 3 - tumor with bacilla, 4 - tumor without bacilla

% mu1 - death of bacilla

% mu2 - death of factor

% mu3 - death of cells of immune response (ther experiments show that the average lifetime of a lymphocyte is 24.25 days,

% giving a death rate of 1/24.25 = 0.412 Philis )

% beta - proportion rate of infection of tumor cells with bacilla

% par1 - death bacilla by macrophag

% par2 - tumor cells stimulate growth of effector cells

%

%global mu1 mu3 beta alpha par1 par2 par3 gRate b_incj count ;

ParamBCG;

MAINT = 1;

global muE muB;

global pK pR pA pT;

global r;

global alpha pD K g_T IFN

% TU_i = 7.8*10^6;%100;%2*10^7;

DOSE = 2.2*10^6;

t0=0;

xx1 = [0 20 0 8*10^6 ];

%x1=zeros(12100,4);

step = 0.1;

dose_length = 7;

xxt=xx1(4);

YINDEX = 500;

urot_State(1:YINDEX) = struct('growth_rate',0,'init_size',0,'result',zeros(2100),'t1',zeros(2100),'color','r', 'lnwdth',1.5);

sd=0;%stable desease

pr=0;nr=0;cr=0;

T0_ini= 10^6;%

T0_fin= 7.6*10^7;

T0_step= 10; %WAS 5; %determines how many elements in the vector

T0=linspace(T0_ini,T0_fin,T0_step);

DOSE = 2.2*10^6;%was 2.2

r0_ini= 0.016;%MICHAEL 0.008

r0_fin= 0.024;

% r0_ini= 0.0025;

% r0_fin= 0.0038;%%3.7.22- was 0.0015-0.0025 ; was 98

r0_step = 5; %WAS10;%detemines how many elelments in the vector

r0=linspace(r0_ini,r0_fin,r0_step);

kk=size(r0);

f_r=randperm(kk(2));

kk=size(T0);

f_T = randperm(kk(2));

k=1;

for i_T=1:T0_step;

for i_r=1:r0_step;

gNum = f_r(i_r);

gRate = r0(gNum);

t0=0; tf=7;

tspan = [t0 tf];

gNumT=f_T(i_T);

t_0 = T0(gNumT)

% xx1 = [b_incj APC_i AT_i AB_i EB_i ET_i IL2_i TI_i t_0];

xx1 = [0 20 0 t_0];

%before treatment

for tf=1:1:8

tspan = t0:step:tf*dose_length;

lt=length(tspan);

sx=size(xx1,1);

x0=xx1(sx,:);

% t0=t0+dose_length;

options = odeset('RelTol',0.001,'AbsTol',[0.001 0.001 0.001 0.001 ]);

[tt,xx1]= ode15s (@odebladder_exp,tspan,x0,options,muB,muE,alpha,pA,pT,pR,pK,pD,gRate,0, 0,K,g_T);

t0=t0+dose_length;

sx=size(xx1,1);

x1((tf-1)*lt+1:lt*tf,:)=xx1;

t1((tf-1)*lt+1:tf*lt)=tt;

end

b_incj = 2.2*10^6; %DOSE ;

for tf=9:1:14 %was 14

tspan = t0:step:tf*dose_length;

lt=length(tspan);

sx=size(xx1,1);

xx1(sx,1)=xx1(sx,1)+b_incj;

x1((tf-1)*lt+1:lt*tf,:)=xx1;

disp(size(xx1));

x0=xx1(sx,:);

disp(size(xx1));

%tic

[tt,xx1]= ode15s (@odebladder_exp,tspan,x0,options,muB,muE,alpha,pA,pT,pR,pK,pD,gRate,0, IFN,K,g_T);

%toc

t0=t0+dose_length;

disp(size(t1));

disp(size(tt));

t1((tf-1)*lt+1:tf*lt)=tt;

disp(size(xx1));

end

%without treatment

for tf=15:1:27

tspan = t0:step:tf*dose_length;

lt=length(tspan);

sx=size(xx1,1);

x0=xx1(sx,:);

%tic

[tt,xx1]= ode15s (@odebladder_exp,tspan,x0,options,muB,muE,alpha,pA,pT,pR,pK,pD,gRate,0, 0,K,g_T);

%toc

t0=t0+dose_length;

sx=size(xx1,1);

x1((tf-1)*lt+1:lt*tf,:)=xx1;

t1((tf-1)*lt+1:tf*lt)=tt;

end

%maintanance treatment

b_incj=0;

kmaint=0;

for tf=28:1:140

tspan = t0:step:tf*dose_length;

lt=length(tspan);

sx=size(xx1,1);

x0=xx1(sx,:);

%tic

[tt,xx1]= ode15s (@odebladder_exp,tspan,x0,options,muB,muE,alpha,pA,pT,pR,pK,pD,gRate,0, IFN,K,g_T);

%toc

t0=t0+dose_length;

sx=size(xx1,1);

if MAINT ==1

if (tf==28 || tf==40 || tf==64 || tf==88 ||tf==112 || tf==136 || kmaint>0 && kmaint<3)

% if (tf==27 || tf==42 || tf==60 || tf==76 || tf==94 || tf==108 || kmaint>0 && kmaint<3)

b_incj=DOSE ;

kmaint=kmaint+1;

else

b_incj = 0;

kmaint=0;

end

xx1(sx,1)=xx1(sx,1)+b_incj;

disp(size(xx1));

x1((tf-1)*lt+1:lt*tf,:)=xx1;

t1((tf-1)*lt+1:tf*lt)=tt;

end

%after treatment

b_incj = 0;

for tf=141:1:170

tspan = t0:step:tf*dose_length;

lt=length(tspan);

sx=size(xx1,1);

x0=xx1(sx,:);

[tt,xx1]= ode15s (@odebladder_exp,tspan,x0,options,muB,muE,alpha,pA,pT,pR,pK,pD,gRate,0, IFN,K,g_T);

t0=t0+dose_length;

sx=size(xx1,1);

x1((tf-1)*lt+1:lt*tf,:)=xx1;

t1((tf-1)*lt+1:tf*lt)=tt;

end

zz = 1:length(t1);

t1=t1(1:length(t1));

adTw1=x1(zz,4);

if adTw1(end)<0.1

cr=cr+1;%complete response

col='-b';

lnwdth = 2.3;

else

if adTw1(end)<=(0.7*adTw1(1))%0.5

lnwdth = 4;

col='--g';

pr=pr+1;%partial response

else

if (adTw1(end)>0.7*adTw1(1))&&(adTw1(end)<=1.2*adTw1(1))

sd=sd+1;%stable desease

col=':k';

lnwdth = 1.5;

else

col='-r';

lnwdth = 0.8;

nr=nr+1; %progressive desease

end

hold on

figure (1);

%subplot(2,1,2);

plot( t1,adTw1,col,'LineWidth',lnwdth);

urot_State(k).color = col;

urot_State(k).growth_rate=gRate;

urot_State(k).init_size=T0(gNumT);

urot_State(k).result=x1;

urot_State(k).t1=t1;

urot_State(k).lnwdth=lnwdth;

k=k+1;

end;

end

t_0 = 26000000

71 4

1 568

71 1

71 4

1 639

71 1

71 4

1 710

71 1

71 4

1 781

71 1

71 4

1 852

71 1

71 4

1 923

71 1

71 4

71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4

t_0 = 26000000

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4

t_0 = 26000000

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

1 12070

71 1

71 4

71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4 71 4

Warning: Failure at t=6.160957e+02. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.818989e-12) at time t.

1 4

Unable to perform assignment because the size of the left side is 71-by-4 and the size of the right side is 1-by-4.

text(350,120e7,[' CR=',num2str(cr)],'Color','b','FontSize',10);

text(350,100e7,[' PR=',num2str(pr)],'Color','g','FontSize',10);

text(350,80e7,[' SD=',num2str(sd)],'Color','k','FontSize',10);

text(350,60e7,[' NR=',num2str(nr)],'Color','r','FontSize',10);

disp([cr ,sd , pr , nr ])

xlabel('days');

ylabel('tumor cells');

title(' BCG induction&maintenance for 50 patients')

text (-45, 23*10^8, 'B ');

% num_k=randperm(500);

% figure(2);

% hold on

% for i_T=1:100;%length(num_k)/2;

% k=num_k(i_T);

% t1=urot_State(k).t1;

% x1=urot_State(k).result;

% col=char(urot_State(k).color);

% t1=t1(1:length(t1));

% zz = 1:length(t1);

% lnwdth = urot_State(k).lnwdth;

%

% adTw1=x1(zz,9);

% plot( t1,adTw1,col,'LineWidth',lnwdth);

% end;

% hold off

% zz = 1:length(t1);

% t1=t1(1:length(t1));

% adB=x1(zz,1);

% adE=x1(zz,2);

% adTb=x1(zz,3);

% adTw=x1(zz,4);

%

% %%plotyy

% hold on

% figure (1);

% subplot(2,1,1);

% %plot( t1,adTw,'-r','LineWidth',3.1);

% % plot( t1,adB,'-k','LineWidth',1.9);

%

% plot( t1,adB,'--k',t1,adTw,'-r','LineWidth',2.1);

% legend('BCG','T');

%

% subplot(2,1,2);

% %plot( t1,adTw,'-r','LineWidth',3.1);

% plot( t1,adB,'-k','LineWidth',1.9);

%

% xlabel('time(days)');

% ylabel('number of cells ');

function xdot = odebladder_exp(t,x,muB,muE,alpha,pA,pT,pR,pK,pD,r,b_injc,IFN,K,g_T)

% ode23s ode45- good for eps = 0

%global mu1 mu3 beta alpha par1 par2 par3 gRate b_incj count ;

%%%%xdot = zeros (4,1);

% 1 - bacilla, 2 - immune response, 3 - tumor with bacilla, 4 - tumor without bacilla

% this time we try without INF

%n=3;

if IFN == 1

Q1= 10; Q2 = 7.82;% Q1= 30; Q2 = 21.82; %Q1=3*10; Q2 = 3*7.3;

else

Q1=1; Q2 = 1;

end

xdot(1) = b_injc+ x(1)*(- pT*x(4) - pA*x(2) -muB) ;%- pT*x(4)

xdot(2) = ((pR*x(1)*x(2))+ (alpha-pD)*x(3)*x(2)+alpha*x(2)*x(3))*Q2 -x(2)*muE;

xdot(3) = pT*x(1)*x(4) - Q1*pK*x(2)*x(3);

%xdot(4) = r*x(4)-pT*x(1)*(1-x(4)/K) - Q1*pK*x(4)*x(2)/(1+5.2*10^3);% WORK!

xdot(4) = r*x(4)*(1-x(4)/K)-pT*x(1)*x(4) - Q1*pK*x(4)*x(2)/(1+g_T);%;*x(4)RIGHT

%xdot(4) = r*x(4)*(1-x(4)/K) - Q1*pK*x(4)*x(2)/(1+5.2*10^3) ;%;*x(4)-pT*x(1)*(1-x(4)/K)

xdot=xdot';

end

function ParamBCG()

global muE muB;

muE = 0.19;

% muE1=0.19;muE2=0.034;

muB = 0.1;

global pA pT pD pK pR %pE = pR

%pA= 2.16*10^(-6); this value do not allow eradication of the tumor

pA= 0.1925*2.16*10^(-6);

pT= 0.028*10^(-2);%0.028*10^(-3);%MICHAEL

pD = 1.03*10^(-10);

pK=0.4*10^(-6);%10^(-6);%(-7);

pR = 0.223e-8;%WAS0.223e-8;

global g_T r;

r = 0.0015;%0.0038;

g_T=5.2*10^3;

global alpha b_injc

alpha = 0.7*10^(-4);%3.7*10^(-6);Michael

b_injc = 2.2e6;

global K;

K = 10^(12);

global IFN;

global MAINT;

IFN = 0;%1;

MAINT = 0;

end

Sam Chak am 4 Dez. 2023

In MATLAB Online öffnen

@Sveta, I just want to check the stability of the original system with some random initial conditions. So far, I haven't encountered any integration failure messages in this relatively simple code. I prefer to include the fixed parameters inside the ODE function for easy reference. Which parameters do you intend to change during the integration?

tspan = [0 1000];

x0 = rand()*[1e7 2e4 3e4 4e4];

[t, x] = ode45(@odebladder_exp, tspan, x0);

plot(t, x), grid on

function xdot = odebladder_exp(t, x)

xdot = zeros (4, 1);

%% parameters (some are unused)

muB = 0.1;

muE = 0.19;

alpha = 0.7e-4;

pA = 0.1925*2.16e-6;

pT = 0.028e-2;

pR = 0.223e-8;

pK = 0.4e-6;

pD = 1.03e-10;

r = 0.0015;

b_injc = 2.2e6;

IFN = 0;

K = 1e12;

g_T = 5.2e3;

Q1 = 1;

Q2 = 1;

%% differential equations

xdot(1) = b_injc + x(1)*(- pT*x(4) - pA*x(2) - muB);

xdot(2) = (pR*x(1)*x(2) + (alpha - pD)*x(3)*x(2) + alpha*x(2)*x(3))*Q2 - x(2)*muE;

xdot(3) = pT*x(1)*x(4) - Q1*pK*x(2)*x(3);

xdot(4) = r*x(4)*(1 - x(4)/K) - pT*x(1)*x(4) - Q1*pK*x(4)*x(2)/(1 + g_T);

end

Melden Sie sich an, um zu kommentieren.

Melden Sie sich an, um diese Frage zu beantworten.

Answer 1

Torsten am 4 Dez. 2023

0
Verknüpfen

Direkter Link zu dieser Antwort

https://de.mathworks.com/matlabcentral/answers/2055674-unable-to-meet-integration-tolerances-without-reducing-the-step-size-below-the-smallest-value-allow#answer_1365299

Bearbeitet: Torsten am 4 Dez. 2023

In MATLAB Online öffnen

I had to add the line

global alpha pD K g_T IFN

at the beginning of your program to make the code work because the above global variables were not available.

In some phase of your code run, the integrator is no longer able to integrate your system of differential equations and gives up with the error message

Warning: Failure at t=6.160957e+02. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.818989e-12) at time t.

After each call of the integrator, you expect an array of size (71x4) as result and you want to save it into a different array. But since the integrator stops already at time 616,.. , the (71x4) array is not available and you make an assignment that would only have been possible if the integrator had succeeded. This gives the error message

Unable to perform assignment because the size of the left side is 71-by-4 and the size of the right side is 1-by-4.

and MATLAB stops.

I cannot help you in this respect since the error stems from the integration of your model equation that I don't know.

1 Kommentar
-1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden

Sveta am 5 Dez. 2023

thanks for your message. I sow matLAB has settled on integration. it is a problem of the function of Matlab - ode45 or ode15s. Maybe is it possible to change them to another function or to change set of parameters in options?

I try to solve it too much time and cannot to find something good...

Please, help me.

Melden Sie sich an, um zu kommentieren.

Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.136868e-13) at time t.

5 Kommentare
3 ältere Kommentare anzeigen3 ältere Kommentare ausblenden

Antworten (1)

1 Kommentar
-1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden

Siehe auch

Kategorien

Tags

Produkte

Version

Community Treasure Hunt

Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.136868e-13) at time t.

5 Kommentare 3 ältere Kommentare anzeigen3 ältere Kommentare ausblenden

Antworten (1)

1 Kommentar -1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden

Siehe auch

Kategorien

Tags

Produkte

Version

Community Treasure Hunt

5 Kommentare
3 ältere Kommentare anzeigen3 ältere Kommentare ausblenden

1 Kommentar
-1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden