Skip to content

Instantly share code, notes, and snippets.

@codersidprogrammer

codersidprogrammer/gaus_seidel.m

Last active September 3, 2023 03:21
Show Gist options
  • Save codersidprogrammer/e8ba9132e02f9bc3b32d0a502757c252 to your computer and use it in GitHub Desktop.
Save codersidprogrammer/e8ba9132e02f9bc3b32d0a502757c252 to your computer and use it in GitHub Desktop.
Download ZIP
Cuplikan kode MATLAB yg digunankan untuk melakukan analisa pada sistem 3 bus
Raw
gaus_seidel.m
% TUGAS PRAKTIKUM KONVERSI ENERGI
% STUDI ALIRAN DAYA PADA SISTEM 3 BUS
% METODE GAUSS-SEIDEL
% -------------------------------------
% TEKNIK ELEKTRO
% UNIVERSITAS MUHAMMADIYAH TANGERANG
% -------------------------------------
% MOCHAMMAD DIMAS EDITIYA - 2220201067
% ARYA HAFIZH ALFARIZI - 2120201003
% TEGAR JUNI KADIYANTO - 2120201023
%Dasar 100 MVA
epsilon = 0.00001;
x = 1;
% Data impedansi
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
% Admitansi pada saluran
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
% Beban dalam per unit
S2 = -(256.6 + j*110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan:
V1 = 1.05 + j*0.0;
Vk1 = conj(V1);
% Estimasi tegangan awal
V2 = 1.0 + j*0.0;
E3 = 1.0 + j*0.0;
iter = 0;
%Format display will be written here
%--
%--
disp('_____________________________________________________________________________________________');
disp('| | | |');
disp('| | DAYA PADA | ALIRAN DAYA PADA SALURAN |');
disp('| Iterasi | SLACK BUS | |');
disp('| Ke: | (BUS 1) | |');
disp('| |---------------------------------------------------------------------------------');
disp('| | S1 | S12 | S21 |');
disp('_____________________________________________________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
Vk2 = conj(V2);
Ek3 = conj(E3);
V2 = 1/(y12 + y23)*((conj(S2))/(Vk2) + (y12)*(V1)+(y23)*(E3));
V3 = 1/(y13 + y23)*((conj(S3))/(Ek3) + (y13)*(V1)+(y23)*(V2));
% Arus pada saluran
I12 = y12*(V1 - V2);
I21 = -I12;
% konjugasi pada saluran
Ik12 = conj(I12);
Ik21 = conj(I21);
%Daya dalam saluran dalam format komplek
S1 = (Vk1 * (V1 * (y12 + y13) - (y12*V2+y13*V3))) * 100;
% Aliran daya dalam bentuk kompleks
S12 = V1 * Ik12 * 100;
S21 = V2 * Ik21 * 100;
x = abs(V3-E3);
E3 = V3;
fprintf(' %i', iter), disp([S1, S12, S21]);
end
% Lanjutan
% Display command will be run here
%-----
disp('_____________________________________________________________________________________________');
disp(' ');
disp('==========================');
disp(' LANJUTAN I ');
disp('==========================');
%Dasar 100 MVA
epsilon = 0.00001;
x = 1;
% Data impedansi
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
% Admitansi pada saluran
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
% Beban dalam per unit
S2 = -(256.6 + j * 110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan :
V1 = 1.05 + j* 0.0;
VK1 = conj (V1);
%Estimasi tegangan awal untuk:
V2 = 1.0 + j * 0.0;
E3 = 1.0 + j * 0.0;
iter = 0;
%Format display will be written here
%--
%--
disp('_____________________________________________________________________________________________');
disp('| | |');
disp('| | ALIRAN DAYA PADA SALURAN |');
disp('| Iterasi | |');
disp('| Ke: | |');
disp('| |---------------------------------------------------------------------------------');
disp('| | S13 | S31 | S23 |');
disp('_____________________________________________________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
Vk2 = conj(V2);
Ek3 = conj(E3);
V2 = 1/(y12 + y23)*((conj(S2))/(Vk2) + (y12)*(V1)+(y23)*(E3));
V3 = 1/(y13 + y23)*((conj(S3))/(Ek3) + (y13)*(V1)+(y23)*(V2));
% Arus pada saluran
I13 = y13*(V1 - V3);
I31 = -I13;
I23 = y23 * (V2 - V3);
% konjugasi pada saluran
Ik13 = conj(I13);
Ik31 = conj(I31);
Ik23 = conj(I23);
% Aliran daya dalam bentuk kompleks
S13 = V1 * Ik13 * 100;
S31 = V3 * Ik31 * 100;
S23 = V2 * Ik23 * 100;
x = abs(V3-E3);
E3 = V3;
fprintf(' %i', iter), disp([S13, S31, S23]);
end
% Lanjutan
% Display command will be run here
%-----
disp('_____________________________________________________________________________________________');
disp(' ');
disp('==========================');
disp(' LANJUTAN II ');
disp('==========================');
%Dasar 100 MVA
epsilon = 0.00001;
x = 1;
% Data impedansi
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
% Admitansi pada saluran
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
% Beban dalam per unit
S2 = -(256.6 + j * 110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan:
V1 = 1.05 + j*0.0;
Vk1 = conj(V1);
% Estimasi tegangan awal
V2 = 1.0 + j*0.0;
E3 = 1.0 + j*0.0;
iter = 0;
%Format display will be written here
%--
%--
disp('_____________________________________________________________________________________________');
disp('| | | |');
disp('| | DAYA PADA | ALIRAN DAYA PADA SALURAN |');
disp('| Iterasi | PADA SALURAN | |');
disp('| Ke: | | |');
disp('| |---------------------------------------------------------------------------------');
disp('| | S32 | SL12 | SL13 |');
disp('_____________________________________________________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
Vk2 = conj(V2);
Ek3 = conj(E3);
V2 = 1/(y12 + y23)*((conj(S2))/(Vk2) + (y12)*(V1)+(y23)*(E3));
V3 = 1/(y13 + y23)*((conj(S3))/(Ek3) + (y13)*(V1)+(y23)*(V2));
% Arus pada saluran
I12 = y12*(V1 - V2);
I21 = -I12;
I13 = y13*(V1 - V3);
I31 = -I13;
I23 = y23 * (V2 - V3);
I32 = -I23;
% konjugasi pada saluran
Ik12 = conj(I12);
Ik21 = conj(I21);
Ik13 = conj(I13);
Ik31 = conj(I31);
Ik32 = conj(I32);
% Aliran daya dalam bentuk kompleks
S12 = V1 * Ik12 * 100;
S21 = V2 * Ik21 * 100;
S13 = V1 * Ik13 * 100;
S31 = V3 * Ik31 * 100;
S32 = V3 * Ik32 * 100;
% Rugi-rugi daya dalam bentuk bilangan kompleks
SL12 = S12 + S21;
SL13 = S13 + S31;
x = abs(V3-E3);
E3 = V3;
fprintf(' %i ', iter), disp([S32, SL12, SL13]);
end
disp('------------------------------------------------------------');
disp(' ');
disp('==========================');
disp(' LANJUTAN III ');
disp('==========================');
%Dasar 100 MVA
epsilon = 0.00001;
x = 1;
% Data impedansi
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
% Admitansi pada saluran
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
% Beban dalam per unit
S2 = -(256.6 + j * 110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan:
V1 = 1.05 + j*0.0;
Vk1 = conj(V1);
% Estimasi tegangan awal
V2 = 1.0 + j*0.0;
E3 = 1.0 + j*0.0;
iter = 0;
%Format display will be written here
%--
%--
disp('_______________________________________________________');
disp('| | |');
disp('| | RUGI-RUGI DAYA PADA SALURAN |');
disp('| Iterasi | |');
disp('| Ke: | |');
disp('| |-------------------------------------------');
disp('| | SL23 |');
disp('_______________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
Vk2 = conj(V2);
Ek3 = conj(E3);
V2 = 1/(y12 + y23)*((conj(S2))/(Vk2) + (y12)*(V1)+(y23)*(E3));
V3 = 1/(y13 + y23)*((conj(S3))/(Ek3) + (y13)*(V1)+(y23)*(V2));
% Arus pada saluran
I23 = y23 * (V2 - V3);
I32 = -I23;
% konjugasi pada saluran
Ik23 = conj(I23);
Ik32 = conj(I32);
% Aliran daya dalam bentuk kompleks
S23 = V2 * Ik23 * 100;
S32 = V3 * Ik32 * 100;
% Rugi-rugi daya dalam bentuk bilangan kompleks
SL23 = S23 + S32;
x = abs(V3-E3);
E3 = V3;
fprintf(' %i ', iter), disp([SL23]);
end
Raw
newton_raphson.m
% TUGAS PRAKTIKUM KONVERSI ENERGI
% STUDI ALIRAN DAYA PADA SISTEM 3 BUS
% METODE NEWTON-RAPHSON
% -------------------------------------
% TEKNIK ELEKTRO
% UNIVERSITAS MUHAMMADIYAH TANGERANG
% -------------------------------------
% MOCHAMMAD DIMAS EDITIYA - 2220201067
% ARYA HAFIZH ALFARIZI - 2120201003
% TEGAR JUNI KADIYANTO - 2120201023
%Dasar 100 MVA
epsilon = 0.001;
x = 1;
% Data impedansi:
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
% Admitansi pada saluran:
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
% Admitansi pada Ybus:
Y11 = (y12 + y13);
Y12 = -y12;
Y13 = -y13;
Y21 = Y12;
Y22 = y12+y23;
Y23 = -y23;
Y31 = -y13;
Y32 = -y23;
Y33 = y13+y23;
% Pembentukan sudut:
theta11 = atan2(imag(Y11), real(Y11));
theta12 = atan2(imag(Y12), real(Y12));
theta13 = atan2(imag(Y13), real(Y13));
theta21 = atan2(imag(Y21), real(Y21));
theta22 = atan2(imag(Y22), real(Y22));
theta23 = atan2(imag(Y23), real(Y23));
theta31 = atan2(imag(Y31), real(Y31));
theta32 = atan2(imag(Y32), real(Y32));
theta33 = atan2(imag(Y33), real(Y33));
% Beban dalam per unit:
S2 = -(256.6 + j*110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan:
V1 = 1.05 + j*0.0;
% Estimasi tegangan awal:
V2 = 1.0 + j*0.0;
V3 = 1.0 + j*0.0;
% Harga delta1:
delta1 = 0;
% Estimasi data awal untuk:
delta2 = 0;
delta3 = 0;
iter = 0;
disp('_____________________________________________________________________________________________');
disp('| | | |');
disp('| | DAYA PADA | ALIRAN DAYA PADA SALURAN |');
disp('| Iterasi | SLACK BUS | |');
disp('| Ke: | (BUS 1) | |');
disp('| |---------------------------------------------------------------------------------');
disp('| | S1 | S12 | S21 |');
disp('_____________________________________________________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
% Harga cos dan sin dari sudut:
sda = cos(theta11);
sdb = cos(theta12-delta1+delta2);
sdc = cos(theta13-delta1+delta3);
sd1 = cos(theta21-delta2+delta1);
sd2 = cos(theta22);
sd3 = cos(theta23-delta2+delta3);
sd4 = cos(theta31-delta3+delta1);
sd5 = cos(theta32-delta3+delta2);
sd6 = cos(theta33);
sd7 = sin(theta21-delta2+delta1);
sd8 = sin(theta22);
sd9 = sin(theta23-delta2+delta3);
sd10 = sin(theta31-delta3+delta1);
sd11 = sin(theta32-delta3+delta2);
sd12 = sin(theta33);
sdd = sin(theta11);
sde = sin(theta12-delta1+delta2);
sdf = sin(theta13-delta1+delta3);
% Harga daya aktif dan daya reaktif pada bus 2 dan bus 3;
P2 = (abs(V2) * abs(V1) * abs(Y21) * sd1) + (abs(V2^2) * abs(Y22) * sd2) + (abs(V2) * abs(V3) * abs(Y23) * sd3);
P3 = (abs(V3) * abs(V1) * abs(Y31) * sd4) + (abs(V3) * abs(V2) * abs(Y32) * sd5) + (abs(V3^2) * abs(Y33) * sd6);
Q2 = -(abs(V2) * abs(V1) * abs(Y21) * sd7) - (abs(V2^2) * abs(Y22) * sd8) - (abs(V2) * abs(V3) * abs(Y23) * sd9);
Q3 = -(abs(V3) * abs(V1) * abs(Y31) * sd10) - (abs(V3) * abs(V2) * abs(Y32) * sd11) - (abs(V3^2) * abs(Y33) * sd12);
% Elemen matriks jacobian:
J(1,1) = (abs(V2) * abs(V1) * abs(Y21) * sd7) + (abs(V2) * abs(V3) * abs(Y23) * sd9);
J(1,2) = -(abs(V2) * abs(V3) * abs(Y23) * sd9);
J(1,3) = (abs(V1) * abs(Y21) * sd1) + (2 * abs(V2) * abs(Y22) * sd2) + (abs(V3) * abs(Y23) *sd3);
J(1,4) = abs(V2) * abs(Y23) * sd3;
J(2,1) = -(abs(V3) * abs(V2) * abs(Y32) * sd11);
J(2,2) = (abs(V3) * abs(V1) * abs(Y31) * sd10) + (abs(V3) * abs(V2) * abs(Y32) * sd11);
J(2,3) = abs(V3) * abs(Y32) * sd5;
J(2,4) = (abs(V1) * abs(Y31) * sd4) + (abs(V2) * abs(Y32) * sd5) + (2 * abs(V3) * abs(Y33) *sd6);
J(3,1) = (abs(V2) * abs(V1) * abs(Y21) * sd1) + (abs(V2) * abs(V3) * abs(Y23) * sd3);
J(3,2) = -(abs(V2) * abs(V3) * abs(Y23) * sd3);
J(3,3) = -(abs(V1) * abs(Y21) * sd7) - (2 * abs(V2) * abs(Y22) * sd8) - (abs(V3) * abs(Y23) * sd9);
J(3,4) = -(abs(V2) * abs(Y23) * sd9);
J(4,1) = -(abs(V3) * abs(V2) * abs(Y32) * sd5);
J(4,2) = (abs(V3) * abs(V1) * abs(Y31) * sd4) + (abs(V3) * abs(V2) * abs(Y32) * sd5);
J(4,3) = -(abs(V3) * abs(Y32) * sd11);
J(4,4) = -(abs(V1) * abs(Y31) * sd10) - (abs(V2) * abs(Y32) * sd11) - (2 * abs(V3) * abs(Y33) * sd12);
% Sisa daya
deltaP2 = real(S2) - P2;
deltaP3 = real(S3) - P3;
deltaQ2 = imag(S2) - Q2;
deltaQ3 = imag(S3) - Q3;
% Sisa daya (power residual) dalam bentuk matrik kolom
deltaS = [deltaP2; deltaP3; deltaQ2; deltaQ3;];
% Perhitungan harga tegangan dan sudut tegangan pada bus 2 dan 3;
deltaX = J\deltaS;
delta2 = delta2 + deltaX(1,1);
delta3 = delta3 + deltaX(2,1);
V2 = V2 + deltaX(3,1);
E3 = (V3 + deltaX(4,1));
% Tegangan dalam bentuk bilangan kompleks pada bus 2 dan 3;
Vkom2 = abs(V2) * (cos(delta2) + j*sin(delta2));
Vkom3 = abs(E3) * (cos(delta3) + j*sin(delta3));
% Arus pada saluran
I12 = y12 * (V1 - Vkom2);
I21 = -I12;
% Konjugat dari arus saluran
Ik12 = conj(I12);
Ik21 = conj(I21);
% Daya aktif dan daya reaktif pada bus 1:
P1 = (abs(V1^2) * abs(Y11) * sda) + (abs(V1) * abs(Vkom2) * abs(Y12) * sdb) + (abs(V1) * abs(Vkom3) * abs(Y13) * sdc);
Q1 = -(abs(V1^2) * abs(Y11) * sdd) - (abs(V1) * abs(Vkom2) * abs(Y12) * sde) - (abs(V1) * abs(Vkom3) * abs(Y13) * sdf);
% Daya kompleks pada bus 1:
S1 = (P1 - j*Q1) * 100;
% Aliran daya dalam bentuk bilangan kompleks pada saluran:
S12 = V1 * Ik12 * 100;
S21 = Vkom2 * Ik21 * 100;
x = abs(E3 - V3);
V3 = E3;
fprintf(' %i ', iter), disp([S1, S12, S21]);
end;
% Lanjutan
% Display command will be run here
%-----
disp('_____________________________________________________________________________________________');
disp(' ');
disp('==========================');
disp('! LANJUTAN I !');
disp('==========================');
% Dasar 100 MVA
epsilon = 0.001;
x = 1;
%Data impedansi pada saluran :
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
%Admintasi pada saluran :
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
%Admintasi Ybus :
Y11 = (y12 + y13);
Y12 = -y12;
Y13 = -y13;
Y21 = Y12;
Y22 = y12+y23;
Y23 = -y23;
Y31 = -y13;
Y32 = -y23;
Y33 = y13+y23;
%Pembentukan sudut :
theta11 = atan2(imag(Y11),real(Y11));
theta12 = atan2(imag(Y12),real(Y12));
theta13 = atan2(imag(Y13),real(Y13));
theta21 = atan2(imag(Y21),real(Y21));
theta22 = atan2(imag(Y22),real(Y22));
theta23 = atan2(imag(Y23),real(Y23));
theta31 = atan2(imag(Y31),real(Y31));
theta32 = atan2(imag(Y32),real(Y32));
theta33 = atan2(imag(Y33),real(Y33));
%Beban dalam per unit:
S2 = -(256.6 +j*110.2)/100;
S3 = -(138.6 +j*45.2)/100;
%Bus 1 sebagai slack bus dengan:
V1 = 1.05 + j*0.0;
%Estimasi tegangan awal untuk:
V2 = 1.0 + j*0.0;
V3 = 1.0 + j*0.0;
%Harga delta1:
delta1 = 0;
%Estimasi delta awal untuk:
delta2 = 0;
delta3 = 0;
iter = 0;
disp('_____________________________________________________________________________________________');
disp('| | |');
disp('| | |');
disp('| Iterasi | ALIRAN DAYA PADA SALURAN |');
disp('| Ke: | |');
disp('| |---------------------------------------------------------------------------------');
disp('| | S13 | S31 | S23 |');
disp('_____________________________________________________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
%Harga cos dan sin dari sudut:
sda = cos(theta11);
sdb = cos(theta12-delta1+delta2);
sdc = cos(theta13-delta1+delta3);
sd1 = cos(theta21-delta2+delta1);
sd2 = cos(theta22);
sd3 = cos(theta23-delta2+delta3);
sd4 = cos(theta31-delta3+delta1);
sd5 = cos(theta32-delta3+delta2);
sd6 = cos(theta33);
sd7 = sin(theta21-delta2+delta1);
sd8 = sin(theta22);
sd9 = sin(theta23-delta2+delta3);
sd10 = sin(theta31-delta3+delta1);
sd11 = sin(theta32-delta3+delta2);
sd12 = sin(theta33);
sdd = sin(theta11);
sde = sin(theta12-delta1+delta2);
sdf = (theta13-delta1+delta3);
%Harga daya aktif dan daya reaktif pada bus 2 dan 3:
P2 = (abs(V2)*abs(V1)*abs(Y21)*sd1)+(abs(V2^2)*abs(Y22)*sd2)+(abs(V2)*abs(V3)*abs(Y23)*sd3);
P3 = (abs(V3)*abs(V1)*abs(Y31)*sd4)+(abs(V3)*abs(V2)*abs(Y32)*sd5)+(abs(V3^2)*abs(Y33)*sd6);
Q2 = -(abs(V2)*abs(V1)*abs(Y21)*sd7)-(abs(V2^2)*abs(Y22)*sd8)-(abs(V2)*abs(V3)*abs(Y23)*sd9);
Q3 = -(abs(V3)*abs(V1)*abs(Y31)*sd10)-(abs(V3)*abs(V2)*abs(Y32)*sd11)-(abs(V3^2)*abs(Y33)*sd12);
%Elemen-elemen matriks Jacobian:
J(1,1) = (abs(V2)*abs(V1)*abs(Y21)*sd7)+(abs(V2)*abs(V3)*abs(Y23)*sd9);
J(1,2) = -(abs(V2)*abs(V3)*abs(Y23)*sd9);
J(1,3) = (abs(V1)*abs(Y21)*sd1)+(2*abs(V2)*abs(Y22)*sd2)+(abs(V3)*abs(Y23)*sd3);
J(1,4) = abs(V2)*abs(Y23)*sd3;
J(2,1) = -(abs(V3)*abs(V2)*abs(Y32)*sd11);
J(2,2) = (abs(V3)*abs(V1)*abs(Y31)*sd10)+(abs(V3)*abs(V2)*abs(Y32)*sd11);
J(2,3) = abs(V3)*abs(Y32)*sd5;
J(2,4) = (abs(V1)*abs(Y31)*sd4)+(abs(V2)*abs(Y32)*sd5)+(2*abs(V3)*abs(Y33)*sd6);
J(3,1) = (abs(V2)*abs(V1)*abs(Y21)*sd1)+(abs(V2)*abs(V3)*abs(Y23)*sd3);
J(3,2) = -(abs(V2)*abs(V3)*abs(Y23)*sd3);
J(3,3) = -(abs(V1)*abs(Y21)*sd7)-(2*abs(V2)*abs(Y22)*sd8)-(abs(V3)*abs(Y23)*sd9);
J(3,4) = -(abs(V2)*abs(Y23)*sd9);
J(4,1) = -(abs(V3)*abs(V2)*abs(Y32)*sd5);
J(4,2) = (abs(V3)*abs(V1)*abs(Y31)*sd4)+(abs(V3)*abs(V2)*abs(Y32)*sd5);
J(4,3) = -(abs(V3)*abs(Y32)*sd11);
J(4,4) = -(abs(V1)*abs(Y31)*sd10)-(abs(V2)*abs(Y32)*sd11)-(2*abs(V3)*abs(Y33)*sd12);
%Sisa daya (power residuals):
deltaP2 = real(S2) - P2;
deltaP3 = real(S3) - P3;
deltaQ2 = imag(S2) - Q2;
deltaQ3 = imag(S3) - Q3;
%Sisa daya (power residuals) dalam bentuk matriks kolom:
deltaS = [deltaP2; deltaP3; deltaQ2; deltaQ3];
%Perhitungan harga tegangan dan sudut tegangan pada bus 2 dan 3:
deltaX = J\deltaS;
delta2 = delta2 + deltaX(1,1);
delta3 = delta3 + deltaX(2,1);
V2 = V2 + deltaX(3,1);
E3 = (V3 + deltaX(4,1));
%Tegangan dalam bentuk bilangan komplek pada bus 2 dan 3:
Vkom2 = abs(V2)*(cos(delta2)+j*sin(delta2));
Vkom3 = abs(E3)*(cos(delta3)+j*sin(delta3));
%Arus pada saluran:
I13 = y13*(V1-Vkom3);
I31 = -I13;
I23 = y23*(Vkom2-Vkom3);
%Konjugat dari arus pada saluran:
Ik13 = conj(I13);
Ik31 = conj(I31);
Ik23 = conj(I23);
%Aliran daya dalam bentuk bilangan komplek pada saluran:
S13 = V1*Ik13*100;
S31 = Vkom3*Ik31*100;
S23 = Vkom2*Ik23*100;
x = abs(E3-V3);
V3 = E3;
fprintf(' %i ', iter), disp([S13, S31, S23]);
end;
% Lanjutan
% Display command will be run here
%-----
disp('_____________________________________________________________________________________________');
disp(' ');
disp('==========================');
disp('! LANJUTAN II !');
disp('==========================');
% Dasar 100 MVA
epsilon = 0.001;
x = 1;
%Data impedansi pada saluran:
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
%Admitansi pada saluran:
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
%Admitansi Ybus:
Y11 = (y12+y13);
Y12 = -y12;
Y13 = -y13;
Y21 = Y12;
Y22 = y12+y23;
Y23 = -y23;
Y31 = -y13;
Y32 = -y23;
Y33 = y13+y23;
% Pembentukan sudut :
theta11 = atan2(imag(Y11),real(Y11));
theta12 = atan2(imag(Y12),real(Y12));
theta13 = atan2(imag(Y13),real(Y13));
theta21 = atan2(imag(Y21),real(Y21));
theta22 = atan2(imag(Y22),real(Y22));
theta23 = atan2(imag(Y23),real(Y23));
theta31 = atan2(imag(Y31),real(Y31));
theta32 = atan2(imag(Y32),real(Y32));
theta33 = atan2(imag(Y33),real(Y33));
% Beban dalam per unit :
S2 = -(256.6 + j*110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan :
V1 = 1.05 + j*0.0;
% Estimasi tegangan awal untuk :
V2 = 1.0 + j*0.0;
V3 = 1.0 + j*0.0;
% Harga delta 1:
delta1 = 0;
% Estimasi delta awal:
delta2 = 0;
delta3 = 0;
iter = 0;
disp('_____________________________________________________________________________________________');
disp('| | | |');
disp('| | ALIRAN DAYA | ALIRAN DAYA PADA SALURAN |');
disp('| Iterasi | PADA SALURAN | |');
disp('| Ke: | | |');
disp('| |---------------------------------------------------------------------------------');
disp('| | S32 | SL12 | SL13 |');
disp('_____________________________________________________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
% Harga cos dan sin dari sudut :
sda = cos(theta11);
sdb = cos(theta12-delta1+delta2);
sdc = cos(theta13-delta1+delta3);
sd1 = cos(theta21-delta2+delta1);
sd2 = cos(theta22);
sd3 = cos(theta23-delta2+delta3);
sd4 = cos(theta31-delta3+delta1);
sd5 = cos(theta32-delta3+delta2);
sd6 = cos(theta33);
sd7 = sin(theta21-delta2+delta1);
sd8 = sin(theta22);
sd9 = sin(theta23-delta2+delta3);
sd10 = sin(theta31-delta3+delta1);
sd11 = sin(theta32-delta3+delta2);
sd12 = sin(theta33);
sdd = sin(theta11);
sde = sin(theta12-delta1+delta2);
sdf = sin(theta13-delta1+delta3);
%Harga daya aktif dan daya reaktif pada bus 2 dan 3:
P2 = (abs(V2)*abs(V1)*abs(Y21)*sd1)+(abs(V2^2)*abs(Y22)*sd2)+(abs(V2)*abs(V3)*abs(Y23)*sd3);
P3 = (abs(V3)*abs(V1)*abs(Y31)*sd4)+(abs(V3)*abs(V2)*abs(Y32)*sd5)+(abs(V3^2)*abs(Y33)*sd6);
Q2 = -(abs(V2)*abs(V1)*abs(Y21)*sd7)-(abs(V2^2)*abs(Y22)*sd8)-(abs(V2)*abs(V3)*abs(Y23)*sd9);
Q3 = -(abs(V3)*abs(V1)*abs(Y31)*sd10)-(abs(V3)*abs(V2)*abs(Y32)*sd11)-(abs(V3^2)*abs(Y33)*sd12);
%Elemen-elemen matriks Jacobian:
J(1,1) = (abs(V2)*abs(V1)*abs(Y21)*sd7)+(abs(V2)*abs(V3)*abs(Y23)*sd9);
J(1,2) = -(abs(V2)*abs(V3)*abs(Y23)*sd9);
J(1,3) = (abs(V1)*abs(Y21)*sd1)+(2*abs(V2)*abs(Y22)*sd2)+(abs(V3)*abs(Y23)*sd3);
J(1,4) = abs(V2)*abs(Y23)*sd3;
J(2,1) = -(abs(V3)*abs(V2)*abs(Y32)*sd11);
J(2,2) = (abs(V3)*abs(V1)*abs(Y31)*sd10)+(abs(V3)*abs(V2)*abs(Y32)*sd11);
J(2,3) = abs(V3)*abs(Y32)*sd5;
J(2,4) = (abs(V1)*abs(Y31)*sd4)+(abs(V2)*abs(Y32)*sd5)+(2*abs(V3)*abs(Y33)*sd6);
J(3,1) = (abs(V2)*abs(V1)*abs(Y21)*sd1)+(abs(V2)*abs(V3)*abs(Y23)*sd3);
J(3,2) = -(abs(V2)*abs(V3)*abs(Y23)*sd3);
J(3,3) = -(abs(V1)*abs(Y21)*sd7)-(2*abs(V2)*abs(Y22)*sd8)-(abs(V3)*abs(Y23)*sd9);
J(3,4) = -(abs(V2)*abs(Y23)*sd9);
J(4,1) = -(abs(V3)*abs(V2)*abs(Y32)*sd5);
J(4,2) = (abs(V3)*abs(V1)*abs(Y31)*sd4)+(abs(V3)*abs(V2)*abs(Y32)*sd5);
J(4,3) = -(abs(V3)*abs(Y32)*sd11);
J(4,4) = -(abs(V1)*abs(Y31)*sd10)-(abs(V2)*abs(Y32)*sd11)-(2*abs(V3)*abs(Y33)*sd12);
%Sisa daya (power residuals):
deltaP2 = real(S2) - P2;
deltaP3 = real(S3) - P3;
deltaQ2 = imag(S2) - Q2;
deltaQ3 = imag(S3) - Q3;
%Sisa daya (power residuals) dalam bentuk matriks kolom:
deltaS = [deltaP2; deltaP3; deltaQ2; deltaQ3];
%Perhitungan harga tegangan dan sudut tegangan pada bus 2 dan 3:
deltaX = J\deltaS;
delta2 = delta2 + deltaX(1,1);
delta3 = delta3 + deltaX(2,1);
V2 = V2 + deltaX(3,1);
E3 = (V3 + deltaX(4,1));
%Tegangan dalam bentuk bilangan komplek pada bus 2 dan 3:
Vkom2 = abs(V2)*(cos(delta2)+j*sin(delta2));
Vkom3 = abs(E3)*(cos(delta3)+j*sin(delta3));
%Arus pada saluran:
I12 = y12*(V1-Vkom2);
I21 = -I12;
I13 = y13*(V1-Vkom3);
I31 = -I13;
I23 = y23*(Vkom2-Vkom3);
I32 = -I23;
%Konjugat dari arus pada saluran:
Ik12 = conj(I12);
Ik21 = conj(I21);
Ik13 = conj(I13);
Ik31 = conj(I31);
Ik32 = conj(I32);
%Aliran daya dalam bentuk bilangan komplek pada saluran:
S12 = V1*Ik12*100;
S21 = Vkom2*Ik21*100;
S13 = V1*Ik13*100;
S31 = Vkom3*Ik31*100;
S32 = Vkom3*Ik32*100;
SL12 = S12 + S21;
SL13 = S13 + S31;
x = abs(E3-V3);
V3 = E3;
fprintf(' %i ', iter), disp([S32, SL12, SL13]);
end
% Lanjutan
% Display command will be run here
%-----
disp('_____________________________________________________________________________________________');
disp(' ');
disp('==========================');
disp('! LANJUTAN III !');
disp('==========================');
% Dasar 100 MVA
epsilon = 0.001;
x = 1;
%Data impedansi pada saluran:
z12 = 0.02 + j*0.04;
z13 = 0.01 + j*0.03;
z23 = 0.0125 + j*0.025;
%Admitansi pada saluran:
y12 = 1/z12;
y13 = 1/z13;
y23 = 1/z23;
%Admitansi Ybus:
Y11 = (y12+y13);
Y12 = -y12;
Y13 = -y13;
Y21 = Y12;
Y22 = y12+y23;
Y23 = -y23;
Y31 = -y13;
Y32 = -y23;
Y33 = y13+y23;
% Pembentukan sudut :
theta11 = atan2(imag(Y11),real(Y11));
theta12 = atan2(imag(Y12),real(Y12));
theta13 = atan2(imag(Y13),real(Y13));
theta21 = atan2(imag(Y21),real(Y21));
theta22 = atan2(imag(Y22),real(Y22));
theta23 = atan2(imag(Y23),real(Y23));
theta31 = atan2(imag(Y31),real(Y31));
theta32 = atan2(imag(Y32),real(Y32));
theta33 = atan2(imag(Y33),real(Y33));
% Beban dalam per unit :
S2 = -(256.6 + j*110.2)/100;
S3 = -(138.6 + j*45.2)/100;
% Bus 1 sebagai slack bus dengan :
V1 = 1.05 + j*0.0;
% Estimasi tegangan awal untuk :
V2 = 1.0 + j*0.0;
V3 = 1.0 + j*0.0;
% Harga delta 1:
delta1 = 0;
% Estimasi delta awal:
delta2 = 0;
delta3 = 0;
iter = 0;
disp('_______________________________________________________');
disp('| | |');
disp('| | RUGI-RUGI DAYA PADA SALURAN |');
disp('| Iterasi | |');
disp('| Ke: | |');
disp('| |-------------------------------------------');
disp('| | SL23 |');
disp('_______________________________________________________');
format short g
while x >= epsilon
iter = iter + 1;
% Harga cos dan sin dari sudut :
sda = cos(theta11);
sdb = cos(theta12-delta1+delta2);
sdc = cos(theta13-delta1+delta3);
sd1 = cos(theta21-delta2+delta1);
sd2 = cos(theta22);
sd3 = cos(theta23-delta2+delta3);
sd4 = cos(theta31-delta3+delta1);
sd5 = cos(theta32-delta3+delta2);
sd6 = cos(theta33);
sd7 = sin(theta21-delta2+delta1);
sd8 = sin(theta22);
sd9 = sin(theta23-delta2+delta3);
sd10 = sin(theta31-delta3+delta1);
sd11 = sin(theta32-delta3+delta2);
sd12 = sin(theta33);
sdd = sin(theta11);
sde = sin(theta12-delta1+delta2);
sdf = sin(theta13-delta1+delta3);
%Harga daya aktif dan daya reaktif pada bus 2 dan 3:
P2 = (abs(V2)*abs(V1)*abs(Y21)*sd1)+(abs(V2^2)*abs(Y22)*sd2)+(abs(V2)*abs(V3)*abs(Y23)*sd3);
P3 = (abs(V3)*abs(V1)*abs(Y31)*sd4)+(abs(V3)*abs(V2)*abs(Y32)*sd5)+(abs(V3^2)*abs(Y33)*sd6);
Q2 = -(abs(V2)*abs(V1)*abs(Y21)*sd7)-(abs(V2^2)*abs(Y22)*sd8)-(abs(V2)*abs(V3)*abs(Y23)*sd9);
Q3 = -(abs(V3)*abs(V1)*abs(Y31)*sd10)-(abs(V3)*abs(V2)*abs(Y32)*sd11)-(abs(V3^2)*abs(Y33)*sd12);
%Elemen-elemen matriks Jacobian:
J(1,1) = (abs(V2)*abs(V1)*abs(Y21)*sd7)+(abs(V2)*abs(V3)*abs(Y23)*sd9);
J(1,2) = -(abs(V2)*abs(V3)*abs(Y23)*sd9);
J(1,3) = (abs(V1)*abs(Y21)*sd1)+(2*abs(V2)*abs(Y22)*sd2)+(abs(V3)*abs(Y23)*sd3);
J(1,4) = abs(V2)*abs(Y23)*sd3;
J(2,1) = -(abs(V3)*abs(V2)*abs(Y32)*sd11);
J(2,2) = (abs(V3)*abs(V1)*abs(Y31)*sd10)+(abs(V3)*abs(V2)*abs(Y32)*sd11);
J(2,3) = abs(V3)*abs(Y32)*sd5;
J(2,4) = (abs(V1)*abs(Y31)*sd4)+(abs(V2)*abs(Y32)*sd5)+(2*abs(V3)*abs(Y33)*sd6);
J(3,1) = (abs(V2)*abs(V1)*abs(Y21)*sd1)+(abs(V2)*abs(V3)*abs(Y23)*sd3);
J(3,2) = -(abs(V2)*abs(V3)*abs(Y23)*sd3);
J(3,3) = -(abs(V1)*abs(Y21)*sd7)-(2*abs(V2)*abs(Y22)*sd8)-(abs(V3)*abs(Y23)*sd9);
J(3,4) = -(abs(V2)*abs(Y23)*sd9);
J(4,1) = -(abs(V3)*abs(V2)*abs(Y32)*sd5);
J(4,2) = (abs(V3)*abs(V1)*abs(Y31)*sd4)+(abs(V3)*abs(V2)*abs(Y32)*sd5);
J(4,3) = -(abs(V3)*abs(Y32)*sd11);
J(4,4) = -(abs(V1)*abs(Y31)*sd10)-(abs(V2)*abs(Y32)*sd11)-(2*abs(V3)*abs(Y33)*sd12);
%Sisa daya (power residuals):
deltaP2 = real(S2) - P2;
deltaP3 = real(S3) - P3;
deltaQ2 = imag(S2) - Q2;
deltaQ3 = imag(S3) - Q3;
%Sisa daya (power residuals) dalam bentuk matriks kolom:
deltaS = [deltaP2; deltaP3; deltaQ2; deltaQ3];
%Perhitungan harga tegangan dan sudut tegangan pada bus 2 dan 3:
deltaX = J\deltaS;
delta2 = delta2 + deltaX(1,1);
delta3 = delta3 + deltaX(2,1);
V2 = V2 + deltaX(3,1);
E3 = (V3 + deltaX(4,1));
%Tegangan dalam bentuk bilangan komplek pada bus 2 dan 3:
Vkom2 = abs(V2)*(cos(delta2)+j*sin(delta2));
Vkom3 = abs(E3)*(cos(delta3)+j*sin(delta3));
%Arus pada saluran:
I23 = y23*(Vkom2-Vkom3);
I32 = -I23;
%Konjugat dari arus pada saluran:
Ik23 = conj(I23);
Ik32 = conj(I32);
%Aliran daya dalam bentuk bilangan komplek pada saluran:
S23 = Vkom2*Ik23*100;
S32 = Vkom3*Ik32*100;
SL23 = S23 + S32;
x = abs(E3-V3);
V3 = E3;
fprintf(' %i ', iter), disp([SL23]);
end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment