TrabalhoC2

Prévia do material em texto

Trabalho Avaliativo C2 ASEP
Alunos: Matheus Scalzer, Pablo do Vale, William Welff
1) Utilizando a ferramenta computacional MatLab, achar a solução para o sistema elétrico
de potência abaixo, aplicando os métodos matemáticos Gauss-Seidel e Newton-Raphson. 
Deve-se:
a) considerar o erro para cada método.
b) aplicar a técnica flat start na primeira iteração.
c) considerar a barra 3 como flutuante, e:
d)calcular a matriz 
e)calcular o valor real da solução considerando 
f)indicar o método matemático que convergiu primeiro e justificar a resposta.
g)indicar o método matemático que ocupou mais processamento e justificar a resposta.
h)comparar as duas soluções encontradas (uma em cada método) e descrever a análise 
como conclusão do trabalho.
Obs.1: Todo o passo-passo utilizado durante a aplicação de cada método matemático, deve
estar comentado no código do cálculo numérico. 
Obs.2: As respostas encontradas deverão estar grifadas.
% Limpar as variáveis
clc; clear;
% Definição de números significativos (5)
format shortG
1
Método: Gauss-Seidel
% Definição da Y Barra
Y_barra = [ -10i-5i 5i 10i;
 5i -10i-5i 10i;
 10i 10i -10i-10i]
Y_barra = 3×3 complex
 0 - 15i 0 + 5i
 0 + 5i 0 - 15i
 0 + 10i 0 + 10i
% Definição Z Barra
Z_Barra = zeros(3);
for linha = 1:3
 for coluna = 1:3
 Z_Barra(linha,coluna) = Y_barra(linha,coluna)^(-1);
 end
end
Z_Barra
Z_Barra = 3×3 complex
 0 + 0.066667i 0 - 0.2i
 0 - 0.2i 0 + 0.066667i
 0 - 0.1i 0 - 0.1i
% Dados Fixos
V3_esp = 1.1;
Theta3_esp = 0;
V2_esp = 1.3;
P2_esp = 2;
P1 = -4.5;
Q1 = -0.5;
% Flat Start primeira interação
V1 = 1+0i;
Theta2 = 0;
[V2_real, V2_imaj] = pol2cart(Theta2, V2_esp);
V2 = V2_real + V2_imaj*1i;
V1_ant = V1;
V2_ant = V2;
2
conjugado_v2 = conj(V2_ant);
conjugado_v1 = conj(V1_ant);
% Interações
 % As interações só irão parar, quando o erro
 % calculado no teste de convergência for menor
 % que a referência dada no exercício (0.0004)
% Marcação tempo de início para comparação
time_ini = clock;
i = 0;
erro_ref = 0.004;
while true
 
 % Estimando Q2
 
 Q2_est = -imag(conjugado_v2 * (Y_barra(2,1) * V1_ant + ...
 Y_barra(2,2) * V2_ant + Y_barra(2,3) * V3_esp));
 
 % Calculando V2
 
 V2_aux = (1 / Y_barra(2,2)) * (((P2_esp - Q2_est*1i)/conjugado_v2)/ ...
 (Y_barra(2,1) * V1_ant) - (Y_barra(2,3) * V3_esp));
 
 % Ângulo de V2 novo
 
 Theta2_aux = angle(V2_aux);
 % Novo valor de V2 calculado, com o ângulo estimado e real
 % especificado
 
 [V2_real, V2_imaj] = pol2cart(Theta2_aux, V2_esp);
 V2_novo = V2_real + V2_imaj*1i;
 % Novo valor de V1
 
 V1_novo = (1 / Y_barra(1,1)) * (((P1 - Q1*1i)/conjugado_v1) - ...
 (Y_barra(1,2) * V2_novo + Y_barra(1,3) * V3_esp));
 % Teste de Convergência de V1 e V2
 
 d_V1 = abs(V1_novo - V1_ant); % Delta V1
 d_V2 = abs(V2_novo - V2_ant); % Delta V2
 
 fprintf(' -------- Iteração %d ---------',i)
 fprintf('V1 calculado: %.5f<%.3f° pu',V1_novo,round(rad2deg(angle(V1_novo)),3))
 fprintf('V2 calculado: %.5f<%.3f° pu',V2_novo,round(rad2deg(angle(V2_novo)),3))
3
 fprintf('Convergência V1: %.6f',d_V1)
 fprintf('Convergência V2: %.6f',d_V2)
 fprintf('\n')
 % Caso o valores do Delta V1 e V2 for inferior ao erro de referência
 % as interações são encerradas
 
 if (d_V2 < erro_ref) && (d_V1 < erro_ref)
 break;
 end
 % Atualização dos valores para nova interação
 
 i=i+1;
 V2_ant = V2_novo;
 V1_ant = V1_novo;
 conjugado_v2 = conj(V2_novo);
 conjugado_v1 = conj(V1_novo);
 
end
 -------- Iteração 0 ---------
V1 calculado: 1.13251<-16.095° pu
V2 calculado: 1.29752<-3.542° pu
Convergência V1: 0.352617
Convergência V2: 0.080361
 -------- Iteração 1 ---------
V1 calculado: 1.06866<-13.269° pu
V2 calculado: 1.29919<-2.026° pu
Convergência V1: 0.098304
Convergência V2: 0.034402
 -------- Iteração 2 ---------
V1 calculado: 1.07398<-14.517° pu
V2 calculado: 1.29873<-2.531° pu
Convergência V1: 0.026615
Convergência V2: 0.011445
 -------- Iteração 3 ---------
V1 calculado: 1.06940<-14.308° pu
V2 calculado: 1.29882<-2.446° pu
Convergência V1: 0.007049
Convergência V2: 0.001922
 -------- Iteração 4 ---------
V1 calculado: 1.06982<-14.402° pu
V2 calculado: 1.29878<-2.485° pu
Convergência V1: 0.002021
Convergência V2: 0.000883
Theta2_novo = angle(V2_novo);
Theta1_novo = angle(V1_novo);
% Subsistema 2
% Solução
4
P3 = abs(V3_esp * ((V1_novo * Y_barra(1,1) * sin(Theta1_novo - Theta1_novo)) + ...
 (V2_novo * Y_barra(1,2) * sin(Theta1_novo - Theta2_novo)) + ...
 (V3_esp * Y_barra(1,3) * sin(Theta1_novo - Theta3_esp))));
 
perdas = P3 + P2_esp + P1;
fprintf(['P3: %.5f pu\n'''...
 'Perdas: %.5f pu'],P3,perdas)
P3: 4.48507 pu
'Perdas: 1.98507 pu
% Convertendo os valores para a base real
Sbase = 100; % MVA
Vbase = 230; % kV
% Tensão
V3 = abs(V3_esp) * Vbase;
V2 = abs(V2_novo) * Vbase;
V1 = abs(V1_novo) * Vbase;
% Pontencia
P1 = P1 * Sbase;
P2 = P2_esp * Sbase;
P3 = P3 * Sbase;
perdas = perdas * Sbase;
% Final do processamento
time_fim = clock;
tempo_exec = time_fim - time_ini;
% Resultado fianl
fprintf(['V1: %.2f kV\n' ...
 'V2: %.2f kV\n' ...
 'V3: %.2f kV\n' ...
 'P1: %.4f MVA\n' ...
 'P2: %.4f MVA\n' ...
 'P3: %.4f MVA\n' ...
 'Perdas: %.4f MVA'],V1,V1,V1,P1,P2,P3,perdas)
V1: 254.04 kV
V2: 254.04 kV
V3: 254.04 kV
P1: -450.0000 MVA
P2: 200.0000 MVA
P3: 448.5074 MVA
Perdas: 198.5074 MVA
5
% Processamento
fprintf(['Tempo de processamento: %.1f ms\n' ...
 'Quantidade de iterações: %d'],tempo_exec(6)*1000,i)
Tempo de processamento: 53.0 ms
Quantidade de iterações: 4
% Limpando as variáveis
clc; clear;
% Definição de números significativos (5)
format shortG
Método: Newton-Raphson
% Definição Y Barra
Y_barra = [ -10i-5i 5i 10i;
 5i -10i-5i 10i;
 10i 10i -10i-10i]
Y_barra = 3×3 complex
 0 - 15i 0 + 5i
 0 + 5i 0 - 15i
 0 + 10i 0 + 10i
% Definição Z Barra
Z_barra = zeros(3);
for linha = 1:3
 for coluna = 1:3
 Z_barra(linha,coluna) = Y_barra(linha,coluna)^(-1);
 end
end
Z_barra
Z_barra = 3×3 complex
 0 + 0.066667i 0 - 0.2i
 0 - 0.2i 0 + 0.066667i
 0 - 0.1i 0 - 0.1i
% Ddefinição G Barra
G_barra = real(Y_barra)
G_barra = 3×3
 0 0 0
 0 0 0
 0 0 0
6
% Definição B Barra
B_barra = imag(Y_barra)
B_barra = 3×3
 -15 5 10
 5 -15 10
 10 10 -20
% Dados Fixos
V3_esp = 1.1;
Theta3_esp = 0;
V2_esp = 1.3;
P2_esp = 2;
P1_esp = -4.5;
Q1_esp = -0.5;
% Flat Start primeira interação
V1 = 1;
Theta1 = 0;
Theta2 = 0;
% Interações
% Marcaçãod o tempo de início para comparação de
% processamento de cada método
time_ini = clock;
i = 0;
erro_ref = 0.0007;
while true
 
% Calculando P1, P2 e Q1
P1_aux = V1 * (V1 * (G_barra(1,1) * cos(Theta1 - Theta1) + ...
B_barra(1,1) * sin(Theta1 - Theta1)) + V2_esp * (G_barra(1,2)* ...
cos(Theta1 - Theta2) + B_barra(1,2) * sin(Theta1 - Theta2)) + ...
V3_esp * (G_barra(1,3) * cos(Theta1 - Theta3_esp) + B_barra(1,3) * ...
sin(Theta1 - Theta3_esp)));P2_aux = V2_esp * (V1 * (G_barra(2,1) * ...
cos(Theta2 - Theta1) + B_barra(2,1) * sin(Theta2 - Theta1)) + ...
V2_esp * (G_barra(2,2) * cos(Theta2 - Theta2) + B_barra(2,2) * ...
sin(Theta2 - Theta2)) + V3_esp * (G_barra(2,3) * cos(Theta2 - ...
Theta3_esp) + B_barra(2,3) * sin(Theta2 - Theta3_esp)));
7
Q1_aux = V1 * (V1 * (G_barra(1,1)* sin(Theta1 - Theta1) - B_barra(1,1) * ...
cos(Theta1 - Theta1)) + V2_esp * (G_barra(1,2) * sin(Theta1 - Theta2) - ...
B_barra(1,2) * cos(Theta1 - Theta2)) + V3_esp * (G_barra(1,3) * ...
sin(Theta1 - Theta3_esp) - B_barra(1,3) * cos(Theta1 - Theta3_esp)));
% Calculando os deltas
d_P1 = P1_esp + P1_aux;
d_P2 = P2_esp + P2_aux;
d_Q1 = Q1_esp + Q1_aux;
% Matriz de deltas
m_deltas = [d_P1;
 d_P2;
 d_Q1];
fprintf(' -------- Iteração %d ---------',i)
fprintf('P1: Delta: %.6f',d_P1)
fprintf('P2: Delta: %.6f',d_P2)
fprintf('Q1: Delta: %.6f',d_Q1)
fprintf('\n')
% Teste de convergência
if (abs(d_P1) < erro_ref) && (abs(d_P2) < erro_ref) && (abs(d_Q1) < erro_ref)
break;
end
% Construindo a matriz Jacobiana
H11 = V1 * (V2_esp * (-G_barra(1,2) * sin(Theta1 - Theta2) + B_barra(1,2) * ...
cos(Theta1 - Theta2)) + V3_esp * (-G_barra(1,3) * sin(Theta1 - Theta3_esp) + ...
B_barra(1,3) * cos(Theta1 - Theta3_esp)));
H12 = V1 * V2_esp * (G_barra(1,2) * sin(Theta1 - Theta2) - B_barra(1,2) * ...
cos(Theta1 - Theta2));
H21 = V2_esp * V1 * (G_barra(2,1) * sin(Theta2 - Theta1) - B_barra(2,1) * ...
cos(Theta2 - Theta1));H22 = V2_esp * (V1 * (-G_barra(2,1) * sin(Theta2 - Theta1) ...
+ B_barra(2,1) * cos(Theta2 - Theta1)) + V3_esp * (-G_barra(2,3) * ...
sin(Theta2 - Theta3_esp) + B_barra(2,3) * cos(Theta2 - Theta3_esp)));
N11 = 2 * V1 * G_barra(1,1) + (V2_esp * (G_barra(1,2) * cos(Theta1 - Theta2) + ...
B_barra(1,2) * sin(Theta1 - Theta2)) + V3_esp * (G_barra(1,3) * ...
cos(Theta1 - Theta3_esp) + B_barra(1,3) * sin(Theta1 - Theta3_esp)));
N21 = V2_esp * (G_barra(2,1) * cos(Theta2 - Theta1) + B_barra(2,1) * ...
sin(Theta2 - Theta1));
M11 = V1 * ((V2_esp * (G_barra(1,2) * cos(Theta1 - Theta2) + B_barra(1,2) * ...
sin(Theta1 - Theta2)) + V3_esp * (G_barra(1,3) * cos(Theta1 - Theta3_esp) + ...
B_barra(1,3) * sin(Theta1 - Theta3_esp))));
8
M12 = - V1 * V2_esp * (G_barra(1,2) * cos(Theta1 - Theta2) + B_barra(1,2) * ...
sin(Theta1 - Theta2));
L11 = -2 * V1 * B_barra(1,1) + (V2_esp * (G_barra(1,2) * sin(Theta1 - Theta2) ...
- B_barra(1,2) * cos(Theta1 - Theta2)) + V3_esp * (G_barra(1,3) * ...
sin(Theta1 - Theta3_esp) - B_barra(1,3) * cos(Theta1 - Theta3_esp)));
% Matriz Jacobiana
m_jacob = - [H11 H12 N11;
H21 H22 N21;
M11 M12 L11];
% Calculando o inverso da matriz Jacobiana e os novos valores de Theta1,
% Theta2e V1 para a próxima iteração
m_jacob_inversa = inv(m_jacob);
% delta de incognitas
d_incognitas = m_jacob_inversa * m_deltas;
Theta1 = Theta1 + d_incognitas(1);
Theta2 = Theta2 + d_incognitas(2);
V1 = V1 + d_incognitas(3);
i=i+1;
end
 -------- Iteração 0 ---------
P1: Delta: -4.500000
P2: Delta: 2.000000
Q1: Delta: -3.000000
 -------- Iteração 1 ---------
P1: Delta: 1.018502
P2: Delta: -0.392780
Q1: Delta: 1.578266
 -------- Iteração 2 ---------
P1: Delta: 0.036569
P2: Delta: -0.013617
Q1: Delta: 0.092422
 -------- Iteração 3 ---------
P1: Delta: 0.000072
P2: Delta: -0.000026
Q1: Delta: 0.000383
% Subsistema 2
% Solução
P3 = V3_esp * (V1 * (G_barra(3,1) * cos(Theta3_esp - Theta1) + B_barra(3,1) * ...
sin(Theta3_esp - Theta1)) + V2_esp * (G_barra(3,2) * cos(Theta3_esp - Theta2) + ...
B_barra(3,2) * sin(Theta3_esp - Theta2)) + V3_esp * (G_barra(3,3) * ...
cos(Theta3_esp - Theta3_esp) + B_barra(3,3) * sin(Theta3_esp - Theta3_esp)));
9
Q2 = V2_esp * (V1 * (G_barra(2,1) * sin(Theta2 - Theta1) - B_barra(2,1) * ...
cos(Theta2 - Theta1)) + V2_esp * (G_barra(2,2) * sin(Theta2 - Theta2) - ...
B_barra(2,2) * cos(Theta2 - Theta2)) + V3_esp * (G_barra(2,3) * ...
sin(Theta2 - Theta3_esp) - B_barra(2,3) * cos(Theta2 - Theta3_esp)));
Q3 = V3_esp * (V1 * (G_barra(3,1) * sin(Theta3_esp - Theta1) - B_barra(3,1) * ...
cos(Theta3_esp - Theta1)) + V2_esp * (G_barra(3,2) * sin(Theta3_esp - Theta2) - ...
B_barra(3,2) * cos(Theta3_esp - Theta2)) + V3_esp * (G_barra(3,3) * ...
sin(Theta3_esp - Theta3_esp) - B_barra(3,3) * cos(Theta3_esp - Theta3_esp)));
fprintf(['P3: %.3f pu\n' ...
 'Q2: %.3f pu\n' ...
 'Q3: %.3f pu'],P3, Q2, Q3)
P3: -2.500 pu
Q2: 3.675 pu
Q3: -2.632 pu
% Convertendo os para a base real
Sbase = 100; % MVA
Vbase = 230; % kV
V3 = V3_esp * Vbase;
P3 = P3 * Sbase;
Q2 = Q2 * Sbase;
Q3 = Q3 * Sbase;
time_fim = clock;
tempo_exec = time_fim - time_ini;
% Resusltado Final
fprintf(['V3: %.2f kV\n' ...
 'P2: %.3f MVA\n' ...
 'Q2: %.3f MVAr\n' ...
 'Q3: %.3f MVAr\n'],V3,P3,Q2,Q3)
V3: 253.00 kV
P2: -250.005 MVA
Q2: 367.495 MVAr
Q3: -263.188 MVAr
% Processamento
fprintf(['Tempo de processamento: %.1f ms\n' ...
 'Quantidade de iterações: %d'],tempo_exec(6)*1000,i)
Tempo de processamento: 33.0 ms
Quantidade de iterações: 3
Análise Avaliativa
10
Portanto, considerando os dois indicadores, foi comprovado que o método Newton-Raphson é o melhor.
Apesar da complexidade para codificação do algoritmo e a pequena diferença de tempo de execução, a
diferença de desempenho tende a ser bem mais significativa para sistemas muito maiores.
11