Hi all I have some problem when I run this script. It has two scripts, and the main script is running very slow. I tried to vectorize all the loop I can and the for loop outside is the only loop I can't vectorize.
script one:function [Vstep]=iteratingvol(~)
%--------------engine speces---------------- b=3.7; %bore s=3.9; %stroke l=5.7; %connecting rod length a=(1/2)*s; %crank radius
%x is piston location; at TDC, x=0; at BDC, x=s theta=-10:1:90;
Vstep=(a+l-(a.*cosd(theta)+(l^2-a^2.*(sind(theta)).^2).^0.5))*(b/2)^2*pi;
Vstep=fliplr(Vstep);
main script:P2=177.8034; T2=1185.4; V2=6; n2=0.000048544;
kcomb=1.4; kcomp=1.4; mrc=1.05; dT=5035; %-----------------------slice model set up--------------------------------- slicenm=100;
P=repmat(P2,1,slicenm); T=repmat(T2,1,slicenm);
V=repmat((V2/slicenm),1,slicenm);
n=repmat((n2/slicenm),1,slicenm);
Vstep=iteratingvol();
%------------Iteration terms----------- deltap=1; %incremental pressure for iterations
for i=1:slicenm
%selected slice combust and increase temperature T(1,i)=T2+dT; n(1,i)=(n2/slicenm)*mrc; P(1,i)=n(1,i)*T(1,i)*R*12/V(1,i);
%expand to the same pressure as unburn slices if i==slicenm V(1,i)=(P(1,i)*V(1,i)^kcomb/P(1,i-1))^(1/kcomb); P(1,i)=P(1,i-1); else V(1,i)=(P(1,i)*V(1,i)^kcomb/P(1,i+1))^(1/kcomb); P(1,i)=P(1,i+1); end
T(1,i)=P(1,i)*V(1,i)/(n(1,i)*R*T(1,i)*12);
%compress slices back to orginal volume before combustion
while sum(V(1,:)) > Vstep1(1,i)
V=((P.*V.^kcomp)./(P+deltap)).^(1/kcomp);
P=P+deltap;
T=P.*V./(n.*R.*12);
end
Tr=T;
Vr=V;
Pr=P;
%compress slices into smaller volume due to piston motion
while sum(V(1,:)) > Vstep1(1,i+1)
V=((P.*V.^kcomp)./(P+deltap)).^(1/kcomp);
P=P+deltap;
T=P.*V./(n.*R.*12);
end
end
The script one is used to define term 'Vstep' which is going to use in main script. If 'Vstep' is defined purely by typing each element one by one, then the main script is running fine. But if term 'Vstep' is calculated by another script, then the main script is running super slow.
Let me know if you guys have any thought, I will appreciate that very much.

Antworten (1)

shanghua chen
shanghua chen am 24 Okt. 2017

0 Stimmen

I may retype the scripts.
script one:
function [Vstep]=iteratingvol(~)
%--------------engine speces---------------- b=3.7; %bore s=3.9; %stroke l=5.7; %connecting rod length a=(1/2)*s; %crank radius
%x is piston location; at TDC, x=0; at BDC, x=s theta=-10:1:90;
Vstep=(a+l-(a.*cosd(theta)+(l^2-a^2.*(sind(theta)).^2).^0.5))*(b/2)^2*pi;
Vstep=fliplr(Vstep);
main script:
P2=177.8034; T2=1185.4; V2=6; n2=0.000048544;
kcomb=1.4; kcomp=1.4; mrc=1.05; dT=5035; %-----------------------slice model set up--------------------------------- slicenm=100;
P=repmat(P2,1,slicenm); T=repmat(T2,1,slicenm);
V=repmat((V2/slicenm),1,slicenm);
n=repmat((n2/slicenm),1,slicenm);
Vstep=iteratingvol();
%------------Iteration terms----------- deltap=1; %incremental pressure for iterations
for i=1:slicenm
%selected slice combust and increase temperature T(1,i)=T2+dT; n(1,i)=(n2/slicenm)*mrc; P(1,i)=n(1,i)*T(1,i)*R*12/V(1,i);
%expand to the same pressure as unburn slices if i==slicenm V(1,i)=(P(1,i)*V(1,i)^kcomb/P(1,i-1))^(1/kcomb); P(1,i)=P(1,i-1); else V(1,i)=(P(1,i)*V(1,i)^kcomb/P(1,i+1))^(1/kcomb); P(1,i)=P(1,i+1); end
T(1,i)=P(1,i)*V(1,i)/(n(1,i)*R*T(1,i)*12);
%compress slices back to orginal volume before combustion
while sum(V(1,:)) > Vstep1(1,i)
V=((P.*V.^kcomp)./(P+deltap)).^(1/kcomp);
P=P+deltap;
T=P.*V./(n.*R.*12);
end
Tr=T;
Vr=V;
Pr=P;
%compress slices into smaller volume due to piston motion
while sum(V(1,:)) > Vstep1(1,i+1)
V=((P.*V.^kcomp)./(P+deltap)).^(1/kcomp);
P=P+deltap;
T=P.*V./(n.*R.*12);
end
end

1 Kommentar

Jan
Jan am 24 Okt. 2017
No, do not "retype"it and post in again in the section for answers. Remove this pseudo answer and use the "{} Code" button to format your code. See http://www.mathworks.com/matlabcentral/answers/13205-tutorial-how-to-format-your-question-with-markup#answer_18099

Melden Sie sich an, um zu kommentieren.

Kategorien

Gefragt:

am 24 Okt. 2017

Kommentiert:

Jan
am 24 Okt. 2017

Community Treasure Hunt

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

Start Hunting!

Translated by