Ailili1997

Coursework Task 2

Nov 11th, 2019
138
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. %===========================================%
  2. %Power Electronics in Distribution Networks %
  3. %Coursework Task 2 %
  4. %===========================================%
  5. %Clear Data
  6. clear
  7. clc
  8. %System Parameters
  9. V1=1.01;
  10. ZL12=0.05+1i*0.01;
  11. ZL13=0.04;
  12. ZL34=0.01;
  13. Zload3=2.21;
  14. Zload4=2.08+1i*0.52;
  15. %SOP Parameters
  16. %PG=0.6;
  17. PG=0;
  18. PSOP=0;
  19. QSOP2=0;
  20. QSOP4=0;
  21. %Admittance Matrix
  22. YL12=1/ZL12;
  23. YL13=1/ZL13;
  24. YL34=1/ZL34;
  25. Yload3=1/Zload3;
  26. Yload4=1/Zload4;
  27.  
  28. Y=[YL12+YL13,-YL12,-YL13,0; -YL12,YL12,0,0;-YL13,0,YL13+YL34+Yload3,-YL34;0,0,-YL34,YL34+Yload4];
  29.  
  30. G12=real(Y(1,2));
  31. B12=imag(Y(1,2));
  32.  
  33. G13=real(Y(1,3));
  34. B13=imag(Y(1,3));
  35.  
  36. G14=real(Y(1,4));
  37. B14=imag(Y(1,4));
  38.  
  39. G22=real(Y(2,2));
  40. B22=imag(Y(2,2));
  41.  
  42. G23=real(Y(2,3));
  43. B23=imag(Y(2,3));
  44.  
  45. G24=real(Y(2,4));
  46. B24=imag(Y(2,4));
  47.  
  48. G33=real(Y(3,3));
  49. B33=imag(Y(3,3));
  50.  
  51. G34=real(Y(3,4));
  52. B34=imag(Y(3,4));
  53.  
  54. G44=real(Y(4,4));
  55. B44=imag(Y(4,4));
  56.  
  57. %Define a vector voltage_at_2 to store the calculated V2 at each loop
  58. voltage_at_2=[];
  59.  
  60. %Loop for plotting the figure when PG=0
  61. for n= 1:30000
  62. %resolution defined as 0.0001
  63. PSOP=-1.0001+0.0001*n;
  64. f=@(V2,V3,V4,theta2,theta3,theta4) [
  65. V1*V2*(G12*cos(theta2)+B12*sin(theta2))+V2*V2*G22+V2*V3*(G23*cos(theta2-theta3)+B23*sin(theta2-theta3))+V4*V2*(G24*cos(theta2-theta4)+B24*sin(theta2-theta4))-PG+PSOP;
  66. V1*V3*(G13*cos(theta3)+B13*sin(theta3))+V2*V3*(G23*cos(theta3-theta2)+B23*sin(theta3-theta2))+V3*V3*G33+V4*V3*(G34*cos(theta3-theta4)+B34*sin(theta3-theta4));
  67. V1*V4*(G14*cos(theta4)+B14*sin(theta4))+V2*V4*(G24*cos(theta4-theta2)+B24*sin(theta4-theta2))+V4*V3*(G34*cos(theta4-theta3)+B34*sin(theta4-theta3))+V4*V4*G44-PSOP;
  68. V1*V2*(G12*sin(theta2)-B12*cos(theta2))-V2*V2*B22+V2*V3*(G23*sin(theta2-theta3)-B23*cos(theta2-theta3))+V2*V4*(G24*sin(theta2-theta4)-B24*cos(theta2-theta4))+QSOP2;
  69. V1*V3*(G13*sin(theta3)-B13*cos(theta3))+V2*V3*(G23*sin(theta3-theta2)-B23*cos(theta3-theta2))-V3*V3*B33+V3*V4*(G34*sin(theta3-theta4)-B34*cos(theta3-theta4));
  70. V1*V4*(G14*sin(theta4)-B14*cos(theta4))+V2*V4*(G24*sin(theta4-theta2)-B24*cos(theta4-theta2))+V3*V4*(G34*sin(theta4-theta3)-B34*cos(theta4-theta3))-V4*V4*B44+QSOP4;];
  71.  
  72. fp=@(x) f(x(1),x(2),x(3),x(4),x(5),x(6));
  73. disp(x);
  74. theta2=x(4);
  75. theta3=x(5);
  76. theta4=x(6);
  77. %Define vectors to to store the calculated V2, V3, V4 at each loop
  78. voltage_at_2(n)=x(1);
  79. voltage_at_3(n)=x(2);
  80. voltage_at_4(n)=x(3);
  81. end
  82.  
  83.  
  84. %Loop for plotting the figure when PG=0.6
  85. for n= 1:30000
  86. PG=0.6;
  87. PSOP=-1.0001+0.0001*n;
  88. f=@(V2,V3,V4,theta2,theta3,theta4) [
  89. V1*V2*(G12*cos(theta2)+B12*sin(theta2))+V2*V2*G22+V2*V3*(G23*cos(theta2-theta3)+B23*sin(theta2-theta3))+V4*V2*(G24*cos(theta2-theta4)+B24*sin(theta2-theta4))-PG+PSOP;
  90. V1*V3*(G13*cos(theta3)+B13*sin(theta3))+V2*V3*(G23*cos(theta3-theta2)+B23*sin(theta3-theta2))+V3*V3*G33+V4*V3*(G34*cos(theta3-theta4)+B34*sin(theta3-theta4));
  91. V1*V4*(G14*cos(theta4)+B14*sin(theta4))+V2*V4*(G24*cos(theta4-theta2)+B24*sin(theta4-theta2))+V4*V3*(G34*cos(theta4-theta3)+B34*sin(theta4-theta3))+V4*V4*G44-PSOP;
  92. V1*V2*(G12*sin(theta2)-B12*cos(theta2))-V2*V2*B22+V2*V3*(G23*sin(theta2-theta3)-B23*cos(theta2-theta3))+V2*V4*(G24*sin(theta2-theta4)-B24*cos(theta2-theta4))+QSOP2;
  93. V1*V3*(G13*sin(theta3)-B13*cos(theta3))+V2*V3*(G23*sin(theta3-theta2)-B23*cos(theta3-theta2))-V3*V3*B33+V3*V4*(G34*sin(theta3-theta4)-B34*cos(theta3-theta4));
  94. V1*V4*(G14*sin(theta4)-B14*cos(theta4))+V2*V4*(G24*sin(theta4-theta2)-B24*cos(theta4-theta2))+V3*V4*(G34*sin(theta4-theta3)-B34*cos(theta4-theta3))-V4*V4*B44+QSOP4;];
  95.  
  96. fp=@(x) f(x(1),x(2),x(3),x(4),x(5),x(6));
  97. [x, fval, info] = fsolve (fp, [1.01;1.01;1.01;0;0;0]);
  98. disp(x);
  99. theta2=x(4);
  100. theta3=x(5);
  101. theta4=x(6);
  102. voltage_at_22(n)=x(1);
  103. voltage_at_33(n)=x(2);
  104. voltage_at_44(n)=x(3);
  105. end
  106.  
  107.  
  108. %Loop for plotting the upper and lower boundary
  109. for n=1:30000
  110. ub(n)=1.02;
  111. lb(n)=0.98;
  112. end
  113.  
  114. figure(1)
  115. index=-1.0000:0.0001:1.9999;
  116. plot(index,voltage_at_2)
  117. hold on
  118. plot(index,ub)
  119. hold on
  120. plot(index,lb)
  121. hold on
  122. plot(index,voltage_at_3)
  123. hold on
  124. plot(index,voltage_at_4)
  125. hold on
  126. plot(index,voltage_at_22)
  127.  
  128. hold off
RAW Paste Data