run takes too long!

3 vues (au cours des 30 derniers jours)
arian hoseini
arian hoseini le 21 Jan 2023
Commenté : arian hoseini le 21 Jan 2023
i have this code that it runs and answers are correct but it takes so long...why?is there any solution?
clc
clear;
% percent of the relay failure rate
SE=[0.99 0.9 0.8 0.5 0];
ME=[0 0 0 0 0];
eta=0.1;
for k=1:5
syms Tc
% permanent fault (lpf)
l1=2/8760;
% temporary faults (ltf)
l2=15/8760;
% relay failure (lp)
l3=0.01/8760;
% mal-trip failure rate (lext)
l4=0.00001/8760;
% reclosing switching rate (lrct)
l5=10800;
% main switching rate (lmct)
l6=30857.14;
% back-up switching rate (lbct)
l7=8640;
% relay instantaneous mal-trip (lrt-op)
l8=0.001/Tc;
% common cause failure rate of C and P (lcp)
l9=0.000001;
% manual switching rate to restore C (litc)
l10=0.5;
% manual switching rate to restore X (litx)
l11=0.5;
% failure rate of self-cheking system (lsc)
l12=0.002/8760;
% failure rato of monitoring circuit (lmn)
% l13=0.0005/8760;
% routine inseption rates (lrt)
l19=1/Tc;
% number of relay tested by self-cheking (msc)
m1=720;
% number of relay tested by inspected routine test (mrt)
m2=1;
% repair rates (mp)
m3=1;
% repaire rate of C (mc)
m4=1;
% portion of relay failure rate by self-cheking (lpsc)
l14=(1-eta)*l3*SE(k);
% portion of relay failure rate by monitoring (lpmn)
l15=0.0005/8760;
% portion of relay failure rate not detected (lprt)
l16=((1-eta)*l3*(1-SE(k)-ME(k)));
% portion of relay potential mal-trip failure rate (lpt-sc)
l17=eta*l3*SE(k);
% portion of relay potential mal-trip failure rate (lpt-rt)
l18=eta*l3*(1-SE(k)-ME(k));
A(1,1)=1-(l19+l1+l2+l12+l16+l14+l15+l18+l17+l4+l9);
A(1,2)=l1;A(1,4)=l2;A(1,6)=l12;A(1,7)=l19;A(1,8)=l16;A(1,9)=l14;A(1,10)=l15;
A(1,11)=l18;A(1,12)=l17;A(1,13)=l4;A(1,15)=l9;
A(2,2)=1-l6;
A(2,3)=l6;
A(3,3)=1-(l3+m4);
A(3,1)=m4;
A(3,17)=l3;
A(4,4)=1-l6;
A(4,5)=l6;
A(5,5)=1-l5;
A(5,1)=l5;
A(6,6)=1-(l1+l2+m1);
A(6,1)=m1;A(6,15)=(l1+l2);
A(7,7)=1-(m2+l1+l8);
A(7,1)=m2;A(7,13)=l8;A(7,15)=l1;
A(8,8)=1-(l19+l1);
A(8,10)=l19;A(8,15)=l1;
A(9,9)=1-(l12+l1+l2);
A(9,10)=l12;A(9,15)=l1+l2;
A(10,10)=1-(l1+l2+m3);
A(10,1)=m3;A(10,15)=l1+l2;
A(11,11)=1-(l19+l1+l3);
A(11,10)=l19;A(11,13)=l3;A(11,15)=l1;
A(12,12)=1-(l12+l3+l1+l2);
A(12,10)=l12;A(12,13)=l3;A(12,15)=l1+l2;
A(13,13)=1-l6;
A(13,14)=l6;
A(14,14)=1-(l10+l10);
A(14,1)=l10;A(14,10)=l10;
A(15,15)=1-l7;
A(15,16)=l7;
A(16,16)=1-l11;
A(16,17)=l11;
A(17,17)=1-(m3+m4);
A(17,3)=m3;A(17,10)=m4;
A(17,:)=1;
B=[0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
1];
% format long
P=inv(A)*B;
% R=VPA(P);
% disp(P)
% disp('P1=')
% disp(P(1,1))
% disp('P2=')
% disp(P(2,1))
% disp('P3=')
% disp(P(3,1))
% disp('P4=')
% disp(P(4,1))
% disp('P5=')
% disp(P(5,1))
% disp('P6=')
% disp(P(6,1))
% disp('P7=')
% disp(P(7,1))
% disp('P8=')
% disp(P(8,1))
% disp('P9=')
% disp(P(9,1))
% disp('P10=')
% disp(P(10,1))
% disp('P11=')
% disp(P(11,1))
% disp('P12=')
% disp(P(12,1))
% disp('P13=')
% disp(P(13,1))
% disp('P14=')
% disp(P(14,1))
% disp('P15=')
% disp(P(15,1))
% disp('P16=')
% disp(P(16,1))
% disp('P17=')
% disp(P(17,1))
PI(k)=P(1);
PII(k)=P(2,1)+P(3,1)+P(4,1)+P(5,1);
PIII(k)=P(6,1)+P(7,1)+P(8,1)+P(9,1)+P(10,1)+P(11,1)+P(12,1);
PIV(k)=P(15,1)+P(16,1)+P(17,1);
PV(k)=P(13,1)+P(14,1);
% PI=double(PI);
end
for c=1:1:500
%%
% Pi1(c)=subs(PI(1),Tc,c);
% Pi1(c)=double(Pi1(c));
%
% Pi2(c)=subs(PI(2),Tc,c);
% Pi2(c)=double(Pi2(c));
%
% Pi3(c)=subs(PI(3),Tc,c);
% Pi3(c)=double(Pi3(c));
%
% Pi4(c)=subs(PI(4),Tc,c);
% Pi4(c)=double(Pi4(c));
%
% Pi5(c)=subs(PI(5),Tc,c);
% Pi5(c)=double(Pi5(c));
%%
Piii1(c)=subs(PIII(1),Tc,c);
Piii1(c)=double(Piii1(c));
Piii2(c)=subs(PIII(2),Tc,c);
Piii2(c)=double(Piii2(c));
Piii3(c)=subs(PIII(3),Tc,c);
Piii3(c)=double(Piii3(c));
Piii4(c)=subs(PIII(4),Tc,c);
Piii4(c)=double(Piii4(c));
Piii5(c)=subs(PIII(5),Tc,c);
Piii5(c)=double(Piii5(c));
%%
Piv1(c)=subs(PIV(1),Tc,c);
Piv1(c)=double(Piv1(c));
Piv2(c)=subs(PIV(2),Tc,c);
Piv2(c)=double(Piv2(c));
Piv3(c)=subs(PIV(3),Tc,c);
Piv3(c)=double(Piv3(c));
Piv4(c)=subs(PIV(4),Tc,c);
Piv4(c)=double(Piv4(c));
Piv5(c)=subs(PIV(5),Tc,c);
Piv5(c)=double(Piv5(c));
% Pv(c)=subs(PV,Tc,c);
% Pv(c)=double(Pv(c));
end
% figure(1)
% plot(Pi1(150:end));
% hold on
% plot(Pi2(150:end));
% hold on
% plot(Pi3(150:end));
% hold on
% plot(Pi4(150:end));
% hold on
% plot(Pi5(150:end));
figure(2)
% plot(Pii(150:end));
% hold on
% figure(3)
plot(Piii1(150:end));
hold on
plot(Piii2(150:end));
hold on
plot(Piii3(150:end));
hold on
plot(Piii4(150:end));
hold on
plot(Piii5(150:end));
figure(3)
plot(Piv1(150:end));
hold on
plot(Piv2(150:end));
hold on
plot(Piv3(150:end));
hold on
plot(Piv4(150:end));
hold on
plot(Piv5(150:end));
% figure(5)
% plot(Pv(150:end));
% hold on

Réponse acceptée

Matt J
Matt J le 21 Jan 2023
Modifié(e) : Matt J le 21 Jan 2023
i have this code that it runs and answers are correct but it takes so long...why?is there any solution?
It is unnecessary to use sym variables here. It is much faster to use numeric variables:
clc
clear;
% percent of the relay failure rate
tic;
for Tc=1:500
out(Tc)=doIt(Tc);
end
toc
Elapsed time is 0.111667 seconds.
Piii=[out.PIII];
Piv=[out.PIV];
figure(2)
plot(Piii(:,150:end)'); legend
figure(3)
plot(Piv(:,150:end)'); legend
function out=doIt(Tc)
SE=[0.99 0.9 0.8 0.5 0];
ME=[0 0 0 0 0];
eta=0.1;
for k=1:5
% permanent fault (lpf)
l1=2/8760;
% temporary faults (ltf)
l2=15/8760;
% relay failure (lp)
l3=0.01/8760;
% mal-trip failure rate (lext)
l4=0.00001/8760;
% reclosing switching rate (lrct)
l5=10800;
% main switching rate (lmct)
l6=30857.14;
% back-up switching rate (lbct)
l7=8640;
% relay instantaneous mal-trip (lrt-op)
l8=0.001/Tc;
% common cause failure rate of C and P (lcp)
l9=0.000001;
% manual switching rate to restore C (litc)
l10=0.5;
% manual switching rate to restore X (litx)
l11=0.5;
% failure rate of self-cheking system (lsc)
l12=0.002/8760;
% failure rato of monitoring circuit (lmn)
% l13=0.0005/8760;
% routine inseption rates (lrt)
l19=1/Tc;
% number of relay tested by self-cheking (msc)
m1=720;
% number of relay tested by inspected routine test (mrt)
m2=1;
% repair rates (mp)
m3=1;
% repaire rate of C (mc)
m4=1;
% portion of relay failure rate by self-cheking (lpsc)
l14=(1-eta)*l3*SE(k);
% portion of relay failure rate by monitoring (lpmn)
l15=0.0005/8760;
% portion of relay failure rate not detected (lprt)
l16=((1-eta)*l3*(1-SE(k)-ME(k)));
% portion of relay potential mal-trip failure rate (lpt-sc)
l17=eta*l3*SE(k);
% portion of relay potential mal-trip failure rate (lpt-rt)
l18=eta*l3*(1-SE(k)-ME(k));
A(1,1)=1-(l19+l1+l2+l12+l16+l14+l15+l18+l17+l4+l9);
A(1,2)=l1;A(1,4)=l2;A(1,6)=l12;A(1,7)=l19;A(1,8)=l16;A(1,9)=l14;A(1,10)=l15;
A(1,11)=l18;A(1,12)=l17;A(1,13)=l4;A(1,15)=l9;
A(2,2)=1-l6;
A(2,3)=l6;
A(3,3)=1-(l3+m4);
A(3,1)=m4;
A(3,17)=l3;
A(4,4)=1-l6;
A(4,5)=l6;
A(5,5)=1-l5;
A(5,1)=l5;
A(6,6)=1-(l1+l2+m1);
A(6,1)=m1;A(6,15)=(l1+l2);
A(7,7)=1-(m2+l1+l8);
A(7,1)=m2;A(7,13)=l8;A(7,15)=l1;
A(8,8)=1-(l19+l1);
A(8,10)=l19;A(8,15)=l1;
A(9,9)=1-(l12+l1+l2);
A(9,10)=l12;A(9,15)=l1+l2;
A(10,10)=1-(l1+l2+m3);
A(10,1)=m3;A(10,15)=l1+l2;
A(11,11)=1-(l19+l1+l3);
A(11,10)=l19;A(11,13)=l3;A(11,15)=l1;
A(12,12)=1-(l12+l3+l1+l2);
A(12,10)=l12;A(12,13)=l3;A(12,15)=l1+l2;
A(13,13)=1-l6;
A(13,14)=l6;
A(14,14)=1-(l10+l10);
A(14,1)=l10;A(14,10)=l10;
A(15,15)=1-l7;
A(15,16)=l7;
A(16,16)=1-l11;
A(16,17)=l11;
A(17,17)=1-(m3+m4);
A(17,3)=m3;A(17,10)=m4;
A(17,:)=1;
B=[0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
1];
% format long
P=inv(A)*B;
out.PI(k,1)=P(1);
out.PII(k,1)=P(2,1)+P(3,1)+P(4,1)+P(5,1);
out.PIII(k,1)=P(6,1)+P(7,1)+P(8,1)+P(9,1)+P(10,1)+P(11,1)+P(12,1);
out.PIV(k,1)=P(15,1)+P(16,1)+P(17,1);
out.PV(k,1)=P(13,1)+P(14,1);
% PI=double(PI);
end
end
  2 commentaires
arian hoseini
arian hoseini le 21 Jan 2023
thank u adn another question if u see the plot data1,2,3,4,...are very close and i saw something in articles an papers that they zoom on that part but the zoom plot is on the main plot...can u help me with that.
Matt J
Matt J le 21 Jan 2023
Here is one way you can bring out the differences in the plots, but as you can see they are very tiny. I can't imagine what significance the differences could have:
for Tc=1:1:500
out(Tc)=doIt(Tc);
end
Piii=[out.PIII]; Piii=Piii(2:5,:)-Piii(1,:);
Piv=[out.PIV]; Piv=Piv(2:5,:)-Piv(1,:);
figure(2)
semilogy(Piii(:,150:end)'); legend
figure(3)
semilogy(Piv(:,150:end)'); legend

Connectez-vous pour commenter.

Plus de réponses (1)

the cyclist
the cyclist le 21 Jan 2023
Vectorizing the for loop gives a nice speedup. There are almost certainly other optimizations possible.
clc
clear;
% percent of the relay failure rate
SE=[0.99 0.9 0.8 0.5 0];
ME=[0 0 0 0 0];
eta=0.1;
for k=1:5
syms Tc
% permanent fault (lpf)
l1=2/8760;
% temporary faults (ltf)
l2=15/8760;
% relay failure (lp)
l3=0.01/8760;
% mal-trip failure rate (lext)
l4=0.00001/8760;
% reclosing switching rate (lrct)
l5=10800;
% main switching rate (lmct)
l6=30857.14;
% back-up switching rate (lbct)
l7=8640;
% relay instantaneous mal-trip (lrt-op)
l8=0.001/Tc;
% common cause failure rate of C and P (lcp)
l9=0.000001;
% manual switching rate to restore C (litc)
l10=0.5;
% manual switching rate to restore X (litx)
l11=0.5;
% failure rate of self-cheking system (lsc)
l12=0.002/8760;
% failure rato of monitoring circuit (lmn)
% l13=0.0005/8760;
% routine inseption rates (lrt)
l19=1/Tc;
% number of relay tested by self-cheking (msc)
m1=720;
% number of relay tested by inspected routine test (mrt)
m2=1;
% repair rates (mp)
m3=1;
% repaire rate of C (mc)
m4=1;
% portion of relay failure rate by self-cheking (lpsc)
l14=(1-eta)*l3*SE(k);
% portion of relay failure rate by monitoring (lpmn)
l15=0.0005/8760;
% portion of relay failure rate not detected (lprt)
l16=((1-eta)*l3*(1-SE(k)-ME(k)));
% portion of relay potential mal-trip failure rate (lpt-sc)
l17=eta*l3*SE(k);
% portion of relay potential mal-trip failure rate (lpt-rt)
l18=eta*l3*(1-SE(k)-ME(k));
A(1,1)=1-(l19+l1+l2+l12+l16+l14+l15+l18+l17+l4+l9);
A(1,2)=l1;A(1,4)=l2;A(1,6)=l12;A(1,7)=l19;A(1,8)=l16;A(1,9)=l14;A(1,10)=l15;
A(1,11)=l18;A(1,12)=l17;A(1,13)=l4;A(1,15)=l9;
A(2,2)=1-l6;
A(2,3)=l6;
A(3,3)=1-(l3+m4);
A(3,1)=m4;
A(3,17)=l3;
A(4,4)=1-l6;
A(4,5)=l6;
A(5,5)=1-l5;
A(5,1)=l5;
A(6,6)=1-(l1+l2+m1);
A(6,1)=m1;A(6,15)=(l1+l2);
A(7,7)=1-(m2+l1+l8);
A(7,1)=m2;A(7,13)=l8;A(7,15)=l1;
A(8,8)=1-(l19+l1);
A(8,10)=l19;A(8,15)=l1;
A(9,9)=1-(l12+l1+l2);
A(9,10)=l12;A(9,15)=l1+l2;
A(10,10)=1-(l1+l2+m3);
A(10,1)=m3;A(10,15)=l1+l2;
A(11,11)=1-(l19+l1+l3);
A(11,10)=l19;A(11,13)=l3;A(11,15)=l1;
A(12,12)=1-(l12+l3+l1+l2);
A(12,10)=l12;A(12,13)=l3;A(12,15)=l1+l2;
A(13,13)=1-l6;
A(13,14)=l6;
A(14,14)=1-(l10+l10);
A(14,1)=l10;A(14,10)=l10;
A(15,15)=1-l7;
A(15,16)=l7;
A(16,16)=1-l11;
A(16,17)=l11;
A(17,17)=1-(m3+m4);
A(17,3)=m3;A(17,10)=m4;
A(17,:)=1;
B=[0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
1];
% format long
P=inv(A)*B;
% R=VPA(P);
% disp(P)
% disp('P1=')
% disp(P(1,1))
% disp('P2=')
% disp(P(2,1))
% disp('P3=')
% disp(P(3,1))
% disp('P4=')
% disp(P(4,1))
% disp('P5=')
% disp(P(5,1))
% disp('P6=')
% disp(P(6,1))
% disp('P7=')
% disp(P(7,1))
% disp('P8=')
% disp(P(8,1))
% disp('P9=')
% disp(P(9,1))
% disp('P10=')
% disp(P(10,1))
% disp('P11=')
% disp(P(11,1))
% disp('P12=')
% disp(P(12,1))
% disp('P13=')
% disp(P(13,1))
% disp('P14=')
% disp(P(14,1))
% disp('P15=')
% disp(P(15,1))
% disp('P16=')
% disp(P(16,1))
% disp('P17=')
% disp(P(17,1))
PI(k)=P(1);
PII(k)=P(2,1)+P(3,1)+P(4,1)+P(5,1);
PIII(k)=P(6,1)+P(7,1)+P(8,1)+P(9,1)+P(10,1)+P(11,1)+P(12,1);
PIV(k)=P(15,1)+P(16,1)+P(17,1);
PV(k)=P(13,1)+P(14,1);
% PI=double(PI);
end
cvec = 1:500;
Piii1=double(subs(PIII(1),Tc,cvec));
Piii2=double(subs(PIII(2),Tc,cvec));
Piii3=double(subs(PIII(3),Tc,cvec));
Piii4=double(subs(PIII(4),Tc,cvec));
Piii5=double(subs(PIII(5),Tc,cvec));
Piv1=double(subs(PIV(1),Tc,cvec));
Piv2=double(subs(PIV(2),Tc,cvec));
Piv3=double(subs(PIV(3),Tc,cvec));
Piv4=double(subs(PIV(4),Tc,cvec));
Piv5=double(subs(PIV(5),Tc,cvec));
% figure(1)
% plot(Pi1(150:end));
% hold on
% plot(Pi2(150:end));
% hold on
% plot(Pi3(150:end));
% hold on
% plot(Pi4(150:end));
% hold on
% plot(Pi5(150:end));
figure(2)
% plot(Pii(150:end));
% hold on
% figure(3)
plot(Piii1(150:end));
hold on
plot(Piii2(150:end));
hold on
plot(Piii3(150:end));
hold on
plot(Piii4(150:end));
hold on
plot(Piii5(150:end));
figure(3)
plot(Piv1(150:end));
hold on
plot(Piv2(150:end));
hold on
plot(Piv3(150:end));
hold on
plot(Piv4(150:end));
hold on
plot(Piv5(150:end));
% figure(5)
% plot(Pv(150:end));
% hold on
  1 commentaire
arian hoseini
arian hoseini le 21 Jan 2023
thank u sir

Connectez-vous pour commenter.

Tags

Produits


Version

R2016b

Community Treasure Hunt

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

Start Hunting!

Translated by