Help to speed up the code
显示 更早的评论
Hello everyone, I need to speed up this simulation code because it took around several hours, can somebody help me how to speed it up?
pc=3.16;
T=293;
rwc3s=0.47;
rwc2s=0.27;
rwc3a=0.10;
rwc4af=0.09;
mc=0.4;
mw=0.157;
ms=0.658;
mg=1.129;
wc=mw/mc;
vc=1/((wc*pc)+1);
ro=(3*vc/(4*pi))^1/3;
%% Proposed Model
yg=0.25;
yw=0.15;
RH=0.88;
b1=1231;
b2=7579;
B293=0.3*10^-10;
C=5*10^7;
ER=5364;
De293=((rwc3s*100)^2.024)*(3.2*10^-14);
kr293=(8.05*(10^-10))*((rwc3s*100+rwc3a*100)^0.975);
v=2.06;
cw=((RH-0.75)/0.25)^3;
pw=1;
%De=De293*exp(-b2*(1/T-1/293));
B=B293*exp(-b1*(1/T-1/293));
kr=kr293*exp(-ER*(1/T-1/293));
%% Determine the day
hr=15;
da=24*hr;
dt=1;
Nn=(da)*(1/dt);
time=1:Nn;
time2=1:Nn+1;
%% Simulation MC
n=100000; %this is the number of trial
rlower=0.5;
rupper=125;
r_initial=rlower+(rupper-rlower)*rand(1,n);
m=length(r_initial);
n1=length(time);
%% Pre-allocation speed
n2=length(time2);
Sr=zeros(1,n1);
cst=zeros(1,n1);
alfa=zeros(1,n1);
kd=zeros(1,n1);
De=zeros(1,n1);
rt_inc=zeros(1,n1);
Rt_inc=zeros(1,n1);
Alfa=zeros(m,n1);
Rt_inc1=zeros(m,n1);
rt_inc1=zeros(m,n1);
ro_init1=zeros(m,1);
alfag1=zeros(m,n1);
Vdp1=zeros(m,1);
vdp1=zeros(m,1);
for i = 1:m
ro_init=r_initial(i)*0.001;
L=((4*pi*(wc*pc/pw+1)/3)^(1/3))*ro_init;
for AA=1:Nn
if (AA==1)
rt_inc(1)=ro_init(1)-0.00001; %boundary condition
Rt_inc(1)=ro_init(1)+0.00001; %boundary condition
else
rt_inc(1)=ro_init(1)-0.00001;
Rt_inc(1)=ro_init(1)+0.00001;
end
end
for BB=1:Nn
if Rt_inc(BB)/L<=ro
Sr(BB)=0;
elseif (Rt_inc(BB)/L>=ro)&&(Rt_inc(BB)/L<0.5)
Sr(BB)=4*pi*(Rt_inc(BB)/L)^2;
elseif (0.5<=Rt_inc(BB)/L)&&(Rt_inc(BB)/L<0.5*(2^0.5))
Sr(BB)=(4*pi*(Rt_inc(BB)/L)^2)-(6*pi*(1-(0.5/(Rt_inc(BB)/L))));
elseif (Rt_inc(BB)/L>=0.5*(2^0.5)) && (Rt_inc(BB)/L<0.5*(3^0.5))
syms y x
fun = @(y,x) 8*(Rt_inc(BB)/L)./(sqrt((Rt_inc(BB)/L)^2-(x.^2)-(y.^2)));
ymin=sqrt((Rt_inc(BB)/L)^2-0.5);
xmin=@(x) sqrt((Rt_inc(BB)/L)^2-0.25-x.^2);
Sr(BB)=integral2(fun,ymin,0.5,xmin,0.5);
else
Sr(BB)=0;
end
cst(BB)=Sr(BB)/(4*pi*Rt_inc(BB)^2);
alfa(BB)=1-(rt_inc(BB)/ro_init)^3;
kd(BB)=(B/(alfa(BB)^1.5))+C*(Rt_inc(BB)-ro_init)^4;
De(BB)=De293*(log(1/alfa(BB)))^1.5;
rt_inc(BB+1)=rt_inc(BB)-(dt*((pw*cst(BB)*cw)/((yw+yg)*pc*rt_inc(BB)^2))*1/((1/(kd(BB)*rt_inc(BB)^2))+(((1/rt_inc(BB))-(1/Rt_inc(BB)))/De(BB))+(1/(kr*rt_inc(BB)^2))));
Rt_inc(BB+1)=((v-1)*((rt_inc(BB)^2)/(Rt_inc(BB)^2))*((rt_inc(BB)-rt_inc(BB+1))/dt))*dt+Rt_inc(BB);
end
Sr1(i,:)=Sr;
Alfa(i,:) = alfa ;
Rt_inc1(i,:) = Rt_inc(1:n1);
rt_inc1(i,:) = rt_inc(1:n1);
rt_inc2(i,:) = rt_inc;
ro_init1(i,:)= ro_init;
b=0.0517;
n=1.0145;% 0.6;
Vdp=(1-exp(-b*((2*ro_init*1000)^n)));
vdp=2*b*n*exp(-b*(2*ro_init*1000)^n)*(2*ro_init*1000)^(n-1);
N=6*(10^12)*mc*vdp/(pc*pi*(2*ro_init*1000)^3);
alfag=alfa*vdp;
alfag1(i,:)=alfag;
Vdp1(i,:)=Vdp;
vdp1(i,:)=vdp;
N1(i,:)=N;
end
Thankyou in advance.
Best Regads,
Kevin
采纳的回答
Durganshu
2020-10-19
0 个投票
Well, the code inside the loops can be edited and vectorized to obtain faster results. I would suggest you go through this documentation on vectorization that may help you:
2 个评论
Kevin Isakayoga
2020-10-19
Durganshu
2020-10-26
Sorry for a late reply. I'm actually quite busy. i'll try to speed up your code when I'm free.
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Simulink 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
