initial commit / add all books

21 files changed, 714 insertions, 0 deletions

diff --git a/1076/CH13/EX13.1/13_1.sce b/1076/CH13/EX13.1/13_1.sce
new file mode 100755
index 000000000..397468d9b
--- /dev/null
+++ b/1076/CH13/EX13.1/13_1.sce
@@ -0,0 +1,32 @@
+clear;

+clc

+

+Xd=.7

+pf=.8

+pfa=acos(pf)

+V=1

+I0=1* exp( %i * pfa *-1)

+E0=V+ (%i * Xd * I0)

+E=round(abs(E0)*100)/100

+d0=atand(imag(E0)/real(E0))

+E0=E * exp(%i * d0 * %pi/180)

+Pe0=E*V*sind(d0)/Xd

+Qe0=(E*V*cosd(d0)/Xd)-(V*V/Xd)

+

+mprintf("\n(a)\nPe= %.1f Qe=%.1f E= %.2f load angle=%.1f",Pe0, Qe0, E, d0);

+

+e1=E0

+E1=abs(e1)

+Pe1=1.2* Pe0;

+d1=asind(Pe1* Xd/ (V*E1))

+Qe1=(E1*V*cosd(d1)/Xd)-(V*V/Xd)

+

+mprintf("\n(b)\nPe= %.2f Qe=%.2f E= %.2f load angle=%.1f",Pe1, Qe1, E1, d1);

+

+e2=1.2 * E0

+E2=abs(e2)

+Pe2=Pe0;

+d2=asind(Pe2* Xd/ (V*E2))

+Qe2=(E2*V*cosd(d2)/Xd)-(V*V/Xd)

+

+mprintf("\n(c)\nPe= %.1f Qe=%.2f E= %.2f load angle=%.1f",Pe2, Qe2, E2, d2);

diff --git a/1076/CH13/EX13.10/13_10.sce b/1076/CH13/EX13.10/13_10.sce
new file mode 100755
index 000000000..5f73e6c0f
--- /dev/null
+++ b/1076/CH13/EX13.10/13_10.sce
@@ -0,0 +1,33 @@
+clear

+clc

+

+Pm=1

+Pe1=.25

+d1=round(asin(Pe1)*1000)/1000

+

+Pe2=.5

+d2=round(asin(Pe2)*1000)/1000

+

+d3=0;

+

+A1=((Pe2 * d2)+ ( 1 * cos(d2)))- ((Pe2 * d1)+ ( 1 * cos(d1)))

+

+E=10

+x=.811

+e=1e-3

+a=Pe2

+b=1

+c=A1 -((Pe2 * d2)+ ( 1 * cos(d2)))

+

+while (E>e)

+//for (i=1:4)

+    f=(a*x) +  (b * cos(x)) + c

+    df=a - (b*sin(x))

+    x1=x-(f/df)

+    disp(x1,f,df)

+    E=abs(x1-x)

+    x=x1

+end

+d3=x1 * 180/%pi

+

+mprintf("\ndelta 3 = %.1f deg",d3)

diff --git a/1076/CH13/EX13.11/13_11.sce b/1076/CH13/EX13.11/13_11.sce
new file mode 100755
index 000000000..2f35907df
--- /dev/null
+++ b/1076/CH13/EX13.11/13_11.sce
@@ -0,0 +1,19 @@
+clear

+clc

+

+G1=50

+H1=8

+

+G2=100

+H2=4

+Gb=100

+

+Ha=(H1*G1/Gb) + (H2*G2/Gb) 

+mprintf("\n(a)Ha= %d MJ/MVA", Ha)

+

+Hb=Ha*2

+mprintf("\n(b)Hb= %d MJ/MVA", Hb)

+

+He= (Ha*Hb)/(Ha+Hb)

+mprintf("\n(c)He= %.3f MJ/MVA", He)

+

diff --git a/1076/CH13/EX13.12/13_12.sce b/1076/CH13/EX13.12/13_12.sce
new file mode 100755
index 000000000..89ef99c85
--- /dev/null
+++ b/1076/CH13/EX13.12/13_12.sce
@@ -0,0 +1,15 @@
+clear

+clc

+

+G=100

+f=50

+H=5

+dL=50

+t=.6

+

+

+J = G*H*1e3;

+dJ=dL*1e3*t

+f2=sqrt((J-dJ)/J)*f

+fd=(f-f2)/f;

+mprintf("Freq deviation = %.3f percent", fd*1e2)

diff --git a/1076/CH13/EX13.13/13_13.sce b/1076/CH13/EX13.13/13_13.sce
new file mode 100755
index 000000000..f4b7a7119
--- /dev/null
+++ b/1076/CH13/EX13.13/13_13.sce
@@ -0,0 +1,13 @@
+clear

+clc

+

+Pi=1

+Pma=1.75

+Pmb=.4

+Pmc=1.25

+

+d0=asin(Pi/Pma)

+dm=%pi - asin(Pi/Pmc)

+

+dcc=acosd(((Pi*(dm-d0))- (Pmb*cos(d0))+ (Pmc*cos(dm)))/(Pmc-Pmb))

+mprintf("Critical Clearing angle = %.1f deg", dcc)

diff --git a/1076/CH13/EX13.14/13_14.jpeg b/1076/CH13/EX13.14/13_14.jpeg
new file mode 100755
index 000000000..0d5fe40f1
--- /dev/null
+++ b/1076/CH13/EX13.14/13_14.jpeg
diff --git a/1076/CH13/EX13.14/13_14.sce b/1076/CH13/EX13.14/13_14.sce
new file mode 100755
index 000000000..121c1f7b9
--- /dev/null
+++ b/1076/CH13/EX13.14/13_14.sce
@@ -0,0 +1,150 @@
+clear

+clc

+clf

+

+Sb = 50;

+S=50;

+V =1;

+Xd = 0.2; 

+X1 =0.4; 

+X2 = 0.4;

+H = 2.7;

+E=1.05 ;

+G=1;

+

+M = G*H/(180*50);

+

+pe0 = (E*V/X1);

+d0=asind(S/(Sb*pe0 ));

+Pe0=pe0 * sind(d0);

+

+pe1 = (E*V/(X1+X2+Xd));

+

+pe2 = (E*V/(X1+Xd));

+

+dt=.05

+c_1=dt*dt/M

+

+for i=1:14

+    if i==1 then 

+        m_t(i)=0;

+        m_Pm(i)=Pe0

+        m_sind(i)=sind(d0)

+        m_Pe(i)=S/Sb

+        m_Pa(i)=0

+        m_cPe(i)=c_1 * m_Pa(i)

+        m_dd(i)=0

+        m_d(i)=d0

+    else if i==2 then 

+        m_t(i)=0;

+        m_Pm(i)=pe1

+        m_d(i)=d0

+        m_sind(i)=sind(m_d(i))

+        m_Pe(i)=m_sind(i)*m_Pm(i)

+        m_Pa(i)=(1 - m_Pe(i) + m_Pa(i-1))/2

+        m_cPe(i)=c_1 * m_Pa(i)

+        m_dd(i)=0

+    else

+        m_t(i)=m_t(i-1) +dt;

+        m_Pm(i)=pe1

+        m_dd(i)=m_dd(i-1) + m_cPe(i-1)

+        m_d(i)=m_d(i-1)+m_dd(i)

+        m_sind(i)=sind(m_d(i))

+        m_Pe(i)=m_Pm(i) * m_sind(i)

+        m_Pa(i)=(1 - m_Pe(i))

+        m_cPe(i)=c_1 * m_Pa(i)

+    end

+    end

+    

+end

+res1(:,1)=m_t(:)

+res1(:,2)=m_Pm(:)

+res1(:,3)=m_sind(:)

+res1(:,4)=m_Pe(:)

+res1(:,5)=m_Pa(:)

+res1(:,6)=m_cPe(:)

+res1(:,7)=m_dd(:)

+res1(:,8)=m_d(:)

+res1=round(res1*1000)/1000

+i=1

+head=['' '' 't'  '' '' 'Pm'  '' '' 'sin d'  '' '' 'Pe'  '' '''Pa'  '' '''8.33Pa'  '' '''d delta' 'delta']

+disp(res1, head, "(a)")

+plot(m_t, m_d)

+title('Swing Curve(Blue: Sustained fault; Red: Fault cleared in .1 sec)');

+xlabel('seconds');

+ylabel('degrees');

+

+

+//(b)

+while i<15

+    if i==1 then 

+        m_t2(i)=0;

+        m_Pm2(i)=Pe0

+        m_sind2(i)=sind(d0)

+        m_Pe2(i)=S/Sb

+        m_Pa2(i)=0

+        m_cPe2(i)=c_1 * m_Pa2(i)

+        m_dd2(i)=0

+        m_d2(i)=d0

+    else if i==2 then 

+            m_t2(i)=0;

+            m_Pm2(i)=pe1

+            m_d2(i)=d0

+            m_sind2(i)=sind(m_d2(i))

+            m_Pe2(i)=m_sind2(i)*m_Pm2(i)

+            m_Pa2(i)=(1 - m_Pe2(i) + m_Pa2(i-1))/2

+            m_cPe2(i)=c_1 * m_Pa2(i)

+            m_dd2(i)=0

+        else 

+            m_t2(i)=m_t2(i-1) +dt;

+            if m_t2(i) == .1 then

+                m_Pm2(i)=pe1

+                m_dd2(i)=m_dd2(i-1) + m_cPe2(i-1)

+                m_d2(i)=m_d2(i-1)+m_dd2(i)

+                m_sind2(i)=sind(m_d2(i))

+                m_Pe2(i)=m_Pm2(i) * m_sind2(i)

+                m_Pa2(i)=(1 - m_Pe2(i))

+                m_cPe2(i)=c_1 * m_Pa2(i)

+                

+                i=i+1

+                m_t2(i)=m_t2(i-1)

+                m_Pm2(i)=pe2

+                m_dd2(i)=m_dd2(i-1)

+                m_d2(i)=m_d2(i-1)

+                m_sind2(i)=sind(m_d2(i))

+                m_Pe2(i)=m_Pm2(i) * m_sind2(i)

+                m_Pa2(i)=(1 - m_Pe2(i) + m_Pa2(i-1))/2

+                m_cPe2(i)=c_1 * m_Pa2(i)

+            else 

+                m_Pm2(i)=m_Pm2(i-1)

+                m_dd2(i)=m_dd2(i-1) + m_cPe2(i-1)

+                m_d2(i)=m_d2(i-1)+m_dd2(i)

+                m_sind2(i)=sind(m_d2(i))

+                m_Pe2(i)=m_Pm2(i) * m_sind2(i)

+                m_Pa2(i)=(1 - m_Pe2(i))

+                m_cPe2(i)=c_1 * m_Pa2(i)

+            end

+        end

+    end

+    i=i+1

+end

+res2(:,1)=m_t2(:)

+res2(:,2)=m_Pm2(:)

+res2(:,3)=m_sind2(:)

+res2(:,4)=m_Pe2(:)

+res2(:,5)=m_Pa2(:)

+res2(:,6)=m_cPe2(:)

+res2(:,7)=m_dd2(:)

+res2(:,8)=m_d2(:)

+res2=round(res2*1000)/1000

+disp(res2,head, "(b)")

+

+plot(m_t2, m_d2,'r')

+//(c)

+D0=d0 * %pi/180

+Pi=1

+Dm=%pi -  asin(Pi/pe2)

+

+dcc=acosd(((Pi * (Dm-D0))-(pe1*cos(D0))+(pe2*cos(Dm)))/(pe2 -pe1))

+tcc=.395

+mprintf("\n\n(c) dcc= %.1f deg; clearing time=%.3f sec", dcc,tcc)

diff --git a/1076/CH13/EX13.15/13_15.sce b/1076/CH13/EX13.15/13_15.sce
new file mode 100755
index 000000000..c6e1a3453
--- /dev/null
+++ b/1076/CH13/EX13.15/13_15.sce
@@ -0,0 +1,25 @@
+clear

+clc

+

+P=50

+pe=100

+

+

+E=10

+x=1

+e=1e-3

+

+d0=30/180 * %pi

+while (E>e)

+    f=cos(d0)+cos(x) - ((%pi - d0-x)*sin(x))

+    df= (-%pi+d0+x)*cos(x)

+    x1=x-(f/df)

+

+    E=abs(x1-x)

+    x=x1

+end

+d1=x1 * 180/%pi

+

+P1=sind(d1) *pe

+Pr=P1-P

+mprintf("\nPower Required = %.2f MW",Pr)

diff --git a/1076/CH13/EX13.16/13_16.sce b/1076/CH13/EX13.16/13_16.sce
new file mode 100755
index 000000000..fe4f6386e
--- /dev/null
+++ b/1076/CH13/EX13.16/13_16.sce
@@ -0,0 +1,28 @@
+clear;

+clc

+

+E=1.5

+V=1

+xg=1

+xm=1

+

+xt1=.1

+xt2=.1

+x1=.4

+x2=.4

+xl=(x1*x2)/(x1+x2)

+X1=xl+xg+xt1

+X2=xt2+xm

+xL=.8

+xC=-.8

+

+

+X3=xC

+X=((X1*X2) + (X2*X3) + (X1*X3))/X3

+P1=E*V/X

+mprintf("\n(a) Steady State limit = %.3f", P1)

+

+X3=xL

+X=((X1*X2) + (X2*X3) + (X1*X3))/X3

+P2=E*V/X

+mprintf("\n(b) Steady State limit = %.3f", P2)

diff --git a/1076/CH13/EX13.17/13_17.sce b/1076/CH13/EX13.17/13_17.sce
new file mode 100755
index 000000000..84698061d
--- /dev/null
+++ b/1076/CH13/EX13.17/13_17.sce
@@ -0,0 +1,14 @@
+clear

+clc

+

+f=50

+w=2*%pi*f

+I=8800

+pf=.85

+J=.5*I*w*w*1e-6

+P=60

+MVA=P/pf

+H=J/MVA

+M=J/(180*f)

+

+mprintf("Inertia const= %.3f MJ/MVA, Momentum= %.5f MJ-s/elec deg",H,M)

diff --git a/1076/CH13/EX13.18/13_18.sce b/1076/CH13/EX13.18/13_18.sce
new file mode 100755
index 000000000..2320cb4cc
--- /dev/null
+++ b/1076/CH13/EX13.18/13_18.sce
@@ -0,0 +1,30 @@
+clear

+clc

+

+f=50

+w=2*%pi*f

+I1=25000

+pf1=.8

+J1=.5*I1*w*w*1e-6

+P1=45

+G1=P1/pf1

+H1=J1/G1

+M1=J1/(180*f)

+

+mprintf("machine 1 \nInertia const= %.2f MJ/MVA, Momentum= %.3f MJ-s/elec deg",H1,M1)

+

+I2=9000

+pf2=.75

+J2=.5*I2*w*w*1e-6

+P2=60

+G2=P2/pf2

+H2=J2/G2

+M2=J2/(180*f)

+

+mprintf("\n\nmachine 2 \nInertia const= %.2f MJ/MVA, Momentum= %.5f MJ-s/elec deg",H2,M2)

+

+M=(M1*M2)/(M1+M2)

+Gb=100

+H=M*180 * f /Gb

+mprintf("\n\n Equivalent constant at 100MVA base \nInertia const= %.4f MJ/MVA, Momentum= %.5f MJ-s/elec deg",H,fix(M*100000)/100000)

+

diff --git a/1076/CH13/EX13.19/13_19.sce b/1076/CH13/EX13.19/13_19.sce
new file mode 100755
index 000000000..59de3cf68
--- /dev/null
+++ b/1076/CH13/EX13.19/13_19.sce
@@ -0,0 +1,17 @@
+clear

+clc

+

+Pm=2

+Pi=1

+H=6

+G=1

+f=50

+p=Pi/Pm

+M=G*H/(%pi*f)

+d0=asin(p)

+

+dcc=acos(((p*(%pi - (2*d0)))- (Pi*cos(d0)))/(Pm-Pi))

+mprintf("Critical Clearing angle = %.4f rad\n\n", dcc)

+

+tcc=sqrt(2*M*(dcc-d0)/Pi)

+mprintf("Critical Clearing time = %.3f sec = %.2f cycles", tcc , tcc*50)

diff --git a/1076/CH13/EX13.2/13_2.sce b/1076/CH13/EX13.2/13_2.sce
new file mode 100755
index 000000000..d328c9b92
--- /dev/null
+++ b/1076/CH13/EX13.2/13_2.sce
@@ -0,0 +1,30 @@
+clear;

+clc

+

+P=4

+f=50

+G=200

+H=6

+J=G*H

+mprintf("\n(a)stored energy = %.0f MJ",J)

+

+Pl=120

+Pl1=160

+Pa=Pl-Pl1

+M=J/(180*f)

+a=Pa/M

+mprintf("\n(b)acceleration = %.0f elec deg/ sec sq",a)

+

+c=5

+t=c/f

+dd=.5*a*t*t

+N=120*f/P

+a=a*60/(180*P)

+Nn=N+(a*t)

+mprintf("\n(c)change in power angle = %.1f elec deg, change in speed = %.1f rpm",dd,Nn)

+

+G2=150

+H2=4

+Gb=100

+Heq=((G*H)+(H2*G2))/Gb

+mprintf("\n(d)Equivalent inertia constant = %.0f MJ/MVA",Heq)

diff --git a/1076/CH13/EX13.3/13_3.sce b/1076/CH13/EX13.3/13_3.sce
new file mode 100755
index 000000000..b175fc862
--- /dev/null
+++ b/1076/CH13/EX13.3/13_3.sce
@@ -0,0 +1,45 @@
+clear;

+clc

+

+f=50;

+H=9;

+x=.6

+P=.7

+pf=.8

+pfa=acos(pf)

+V=1

+D=.14

+dd=10 * %pi/180

+

+S=P/pf * exp(%i * pfa)

+I=conj(S/V)

+E=V+ (%i * x *I)

+d0=atand(imag(E)/real(E))

+Pr=abs(E) * V * cosd(d0)/x

+mprintf("\nSynchronising Power Coefficient = %.3f pu",Pr)

+

+w1=sqrt(Pr *f *%pi/H)

+mprintf("\nUndamped Speed Of oscillations = %.2f rad/sec",w1)

+w1=round(w1*100)/100

+

+z=(D/2) * sqrt(%pi * f/(H*Pr))

+mprintf("\nDamping Ratio = %.4f ",z)

+

+wd=w1 *sqrt(1-(z*z))

+Wd=wd / (2*%pi)

+mprintf("\ndamped angular frequency Of oscillations = %.3f rad/sec = %.3f Hz",wd, Wd)

+

+z=round(z*1e4)/1e4

+c_1=10/sqrt(1-z^2)

+c_2=z*w1

+

+mprintf("\n\nd = %.2f + %.3f exp(-%.3f t) sin(%.3f t + %.1f deg)", d0, c_1, c_2,wd, acosd(z))

+

+c_3=w1 * 10 /(360 * sqrt(1-z^2))

+mprintf("\n\nf = %.0f - %.4f exp(-%.3f t) sin(%.3f t)", f, c_3, c_2,wd)

+

+

+

+

+

+

diff --git a/1076/CH13/EX13.4/13_4.sce b/1076/CH13/EX13.4/13_4.sce
new file mode 100755
index 000000000..196b14e6a
--- /dev/null
+++ b/1076/CH13/EX13.4/13_4.sce
@@ -0,0 +1,36 @@
+clear;

+clc

+dP=.1

+f=50;

+H=9;

+x=.6

+P=.7

+pf=.8

+pfa=acos(pf)

+V=1

+D=.14

+dd=10 * %pi/180

+

+S=P/pf * exp(%i * pfa)

+I=conj(S/V)

+E=V+ (%i * x *I)

+d0=atand(imag(E)/real(E))

+Pr=abs(E) * V * cosd(d0)/x

+

+w1=sqrt(Pr *f *%pi/H)

+w1=round(w1*100)/100

+

+z=(D/2) * sqrt(%pi * f/(H*Pr))

+

+wd=w1 *sqrt(1-(z*z))

+Wd=wd / (2*%pi)

+

+z=round(z*1e4)/1e4

+c_1=1/sqrt(1-z^2)

+c_2=z*w1

+c_3=180 * f * dP/(9 * w1^2)

+

+mprintf("\n\nd = %.2f + %.3f(1- (%.3f exp(-%.3f t) sin(%.3f t + %.1f deg)))", d0, c_3,c_1, c_2,wd, acosd(z))

+

+c_4= dP * f /(w1* 9 * sqrt(1-z^2))

+mprintf("\n\nf = %.0f + %.4f exp(-%.3f t) sin(%.3f t)", f, c_4, c_2,wd)

diff --git a/1076/CH13/EX13.5/13_5.sce b/1076/CH13/EX13.5/13_5.sce
new file mode 100755
index 000000000..9eae00b26
--- /dev/null
+++ b/1076/CH13/EX13.5/13_5.sce
@@ -0,0 +1,32 @@
+clear;

+clc

+

+E=1.1

+V=1

+xg=.7

+xt=.1

+x1=.2

+X1=x1+xg+xt

+X2=.2

+xL=.5

+xC=-.5

+

+

+X=X1+X2

+P1=E*V/X

+mprintf("\n(a) Steady State limit = %.4f", P1)

+

+X3=xL

+X=((X1*X2) + (X2*X3) + (X1*X3))/X3

+P2=E*V/X

+mprintf("\n(b) Steady State limit = %.4f", P2)

+

+X3=xC

+X=((X1*X2) + (X2*X3) + (X1*X3))/X3

+P3=E*V/X

+mprintf("\n(c) Steady State limit = %.3f", P3)

+

+X=X1+xC+X2

+P4=E*V/X

+mprintf("\n(d) Steady State limit = %.2f", P4)

+

diff --git a/1076/CH13/EX13.6/13_6.sce b/1076/CH13/EX13.6/13_6.sce
new file mode 100755
index 000000000..6fd103a43
--- /dev/null
+++ b/1076/CH13/EX13.6/13_6.sce
@@ -0,0 +1,18 @@
+clear;

+clc

+

+f=50

+xg=1

+xl=.5

+E=1.1

+V=1

+H=5

+p=.5

+

+X=xl+xg

+d0=asin(p)

+Pr=E*V*cos(d0)/X

+M=H/(%pi*f)

+wn=sqrt(Pr/M)/(2*%pi)

+

+mprintf("Freg of oscillation = %.2f Hz", wn)

diff --git a/1076/CH13/EX13.7/13_7.sce b/1076/CH13/EX13.7/13_7.sce
new file mode 100755
index 000000000..f20820d43
--- /dev/null
+++ b/1076/CH13/EX13.7/13_7.sce
@@ -0,0 +1,45 @@
+clear

+clc

+

+xd=.2

+x1=.4

+x2=.4

+Pi=1.5

+E=1.2

+V=1

+

+

+X=xd  +((x1*x2)/(x1+x2))

+pe=E*V/X

+d0 = asin(Pi/pe)

+

+X2=xd+x1

+pe2=E*V/X2

+d1 = asin(Pi/pe2)

+dm=%pi- d1

+

+A1=((Pi * d1)+ ( pe2 * cos(d1)))- ((Pi * d0)+ ( pe2 * cos(d0)))

+A2=((Pi * dm)+ ( pe2 * cos(dm)))- ((Pi * d1)+ ( pe2 * cos(d1)))

+

+if abs(A1)<abs(A2) then

+    mprintf("STABLE\n\n")

+else

+    mprintf("UNSTABLE\n\n")

+end

+

+E=10

+x=2

+e=1e-3

+a=Pi

+b=pe2

+c=-A1 +(Pi * d1 ) + (cos (d1) * pe2)

+

+while (E>e)

+    f=(a*x) +  (b * cos(x)) - c

+    df=1.5 - (2*sin(x))

+    x1=x-(f/df)

+    E=abs(x1-x)

+    x=x1

+end

+d2=x1 * 180/%pi

+mprintf("\ndelta 2 = %.2f deg",d2)

diff --git a/1076/CH13/EX13.8/13_8.sce b/1076/CH13/EX13.8/13_8.sce
new file mode 100755
index 000000000..8ac137760
--- /dev/null
+++ b/1076/CH13/EX13.8/13_8.sce
@@ -0,0 +1,79 @@
+clear

+clc

+

+xd=.2

+x1=.4

+x2=.4

+Pi=1.5

+E=1.2

+V=1

+

+Xs1=xd  +((x1*x2)/(x1+x2))

+pe=E*V/Xs1

+d0 = asin(Pi/pe)

+dc= (%pi/2)-d0

+dc=round(dc*1e3)/1e3

+

+X1=x1;

+X2=x2/2

+X3=x2/2

+Xs2=((X1*X2) + (X2*X3) + (X1*X3))/X3

+pe2=E*V/Xs2

+

+

+

+Xs3=xd+ x1

+pe3=E*V/Xs3

+d2 = asin(Pi/pe3)

+dm=%pi- d2

+

+

+if pe2<Pi then

+    mprintf("UNSTABLE for sustained fault\n\n")

+else

+    mprintf("STABLE for sustained fault\n\n")

+end

+

+

+A1=((Pi * dc)+ ( pe2 * cos(dc)))- ((Pi * d0)+ ( pe2 * cos(d0)))

+A2=((Pi * dm)+ ( pe3 * cos(dm)))- ((Pi * dc)+ ( pe3 * cos(dc)))

+

+if abs(A1)<abs(A2) then

+    mprintf("STABLE system\n\n")

+else

+    mprintf("UNSTABLE  system\n\n")

+end

+

+

+E=10

+x=2

+e=1e-3

+a=Pi

+b=pe3

+c=-A1 +(Pi * dc ) + (cos (dc) * pe3)

+

+while (E>e)

+    f=(a*x) +  (b * cos(x)) - c

+    df=1.5 - (2*sin(x))

+    x1=x-(f/df)

+    E=abs(x1-x)

+    x=x1

+end

+d2=x1 * 180/%pi

+mprintf("\ndelta 2 = %.2f deg",d2)

+

+Pmb=pe2

+Pmc=pe3

+d0=round(d0*1000)/1000

+dcc=acosd(((Pi*(dm-d0))- (Pmb*cos(d0))+ (Pmc*cos(dm)))/(Pmc-Pmb))

+mprintf("\nCritical Clearing angle = %.1f deg", dcc)

+

+

+

+

+

+

+

+

+

+

diff --git a/1076/CH13/EX13.9/13_9.sce b/1076/CH13/EX13.9/13_9.sce
new file mode 100755
index 000000000..67f411255
--- /dev/null
+++ b/1076/CH13/EX13.9/13_9.sce
@@ -0,0 +1,23 @@
+clear

+clc

+

+xd=.2

+x1=.4

+x2=.4

+Pi=1.5

+E=1.2

+V=1

+

+Xs1=xd  +((x1*x2)/(x1+x2))

+pe=E*V/Xs1

+d0 = asin(Pi/pe)

+

+Xs3=xd+ x1

+pe3=E*V/Xs3

+d2 = asin(Pi/pe3)

+dm=%pi- d2

+

+Pmb=0

+Pmc=pe3

+dcc=acosd(((Pi*(dm-d0))- (Pmb*cos(d0))+ (Pmc*cos(dm)))/(Pmc-Pmb))

+mprintf("Critical Clearing angle = %.2f deg", dcc)

diff --git a/1076/CH13/EX19.20/13_20.sce b/1076/CH13/EX19.20/13_20.sce
new file mode 100755
index 000000000..0e03cfbd4
--- /dev/null
+++ b/1076/CH13/EX19.20/13_20.sce
@@ -0,0 +1,30 @@
+clear

+clc

+

+f=50

+G=20

+V=13.2

+H=9

+nP=4

+

+J=G*H

+mprintf("\n(a)Stored Energy = %.0f MJ", J)

+disp("The unit is incorrectly mentioned as "'mJ'' in the textbook.");

+

+Pi= 25*.735

+P=15

+Pa=Pi-P

+M=G*H/(180*f)

+a=Pa/M

+mprintf("\n(b)Acceleration = %.2f elec deg/sec sq",a )

+

+

+c=15

+t=c/f

+dd=.5*a*t*t

+mprintf("\n(c)change in angle = %.2f deg",dd )

+

+A=a * 60 / (180*nP)

+Ns=120*f/nP

+N=Ns+(A*t)

+mprintf("\n(d)New speed = %.3f rpm",fix(N*1000)/1000 )