TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
1. Introducción
En el presente trabajo se desarrolló el problema 7.4 del libro “Power Generation, Operation
and Control” [1] mediante el método metaheurístico de Optimización por Enjambre de
Partículas. Para esto se implementó un programa en el software matemático MATLAB. Como
primer paso, se desarrolló el ejemplo 7.E del libro [1], a fin de validar la coherencia de los
resultados obtenidos mediante el modelo implementado en Matlab de este método
metaheurístico. Cabe mencionar que en el libro de referencia se desarrolla el ejemplo 7.E
mediante el método matemático de programación dinámica.
2. Ejemplo 7.E
2.1 Enunciado del Ejemplo 7.E
Realizar el despacho económico para el sistema de la figura 1.
Figura 1 Sistema del ejemplo 7.E.
Características de la planta de vapor “S”:
Figura 2 Función de costos incrementales.
Características de la planta hidroeléctrica “H”:
Figura 3 Función de caudal vs Potencia generada.
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
Volumen inicial=Volumen final=10000 acre-ft
6000≤Volumen≤18000
Caudal natural=1000 acre-ft/h
Característica de la demanda
Tabla 1 Característica de la demanda
2.2 Formulación Matemática del problema
𝑗=6
𝑀𝑖𝑛 { ∑ 𝐹(𝑃𝑠(𝑗) ) × 4}
𝑗=1
s.a.
𝑃𝑠(𝑗) + 𝑃𝐻(𝑗) = 𝐷𝑗
200 ≤ 𝑃𝑠(𝑗) ≤ 1200
0 ≤ 𝑃𝐻(𝑗) ≤ 200
6000 ≤ 𝑉(𝑗) ≤ 18000
𝑉(0) = 10000
𝑉(6) = 10000
𝑉(𝑗) = 𝑉(𝑗−1) + 4(1000 − 𝑞(𝑗) )
𝑞(𝑗) = 260 + 10 × 𝑃ℎ(𝑗)
3. Modelamiento en Matlab
Para simplificar el código MATLAB implementamos dos funciones “ofun(x)” y “Resultados(x)”.
La primera función calcula el valor de la función objetivo para el vector “x” ingresado. La
segunda función grafica los resultados de generación y volúmenes para el vector “x” ingresado.
Las restricciones se modelaron con una penalización de R/. 10 000, con el fin de que estos
resultados tengan un valor de función objetivo muy alto respecto al resultado óptimo. Esto con
el fin de evitar que sean consideradas como mejores posiciones durante el desarrollo del
método heurístico PSO. El código implementado en Matlab se presenta en el anexo 1.
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
3.1 Análisis de los resultados obtenidos
En las figuras 4 y tabla 2 se presentan los resultados del ejemplo 7.E obtenidos mediante el
método matemático de Programación Dinámica. En la figura 4 se presenta el volumen del
embalse de la planta hidroeléctrica asociada a la operación óptima del sistema. Mientras en la
tabla 2 se resalta el valor óptimo de la función objetivo.
Figura 4 Comportamiento del volúmen del embalse de la planta hidroeléctrica asociado a la operación óptima del
sistema.
Tabla 2 Función Objetivo optima.
A continuación se presenta los resultados obtenidos mediante el método heurístico de
Enjambre de Partículas para el mismo problema. En la figura 5, se presenta la potencia
generada óptima por cada planta y para cada período. En la figura 6, se presenta el
comportamiento del volumen del embalse de la planta hidroeléctrica asociada a la operación
óptima del sistema. En la tabla 3, presentamos la función objetivo obtenido con el método
PSO. Nótese que los resultados obtenidos del método heurístico son coherentes con el
resultado obtenido del método matemático. Además, es importante hacer notar que los
resultados son diferentes en 0.14 % esto debido a que no necesariamente un método
heurístico llega a la solución óptima global del problema, sin embargo es importante precisar
que para este caso se obtuvo un resultado bastante aproximado al optimo real.
Programación Dinámica
PSO
Diferencia de resultados (%)
Función Objetivo Óptima
81738
81853
0.14%
Tabla 3 Funciones Objetivos optimas obtenidas con metodo matemático ymetodo heurístico.
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
Resumen de Resultados - Generación
1000
Planta de Vapor (MW )
Planta Hidroeléctrica (MW )
900
800
700
(MW)
600
500
400
300
200
100
0
1
2
3
4
Periodo "j"
5
6
Figura 5 Generación de la planta hidroeléctrica y Generación de la planta de vapor para cada período.
1.5
x 10
4
Resumen de Resultados - Volumen(j)
1.4
Volúmen (acre-ft)
1.3
1.2
1.1
1
0.9
0.8
0.7
0.6
1
2
3
4
Periodo "j"
5
6
7
Figura 6 Volúmen del embalse de la planta hidroeléctrica para cada período.
Luego, podemos decir que se obtuvo resultados coherentes de nuestro modelo implementado
en Matlab, basados en la verificación realizada con los resultados del ejemplo 7.E. Por lo tanto,
luego de haber verificado la validez de nuestro modelo, se procedió a solucionar el problema
7.4.
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
4. Problema 7.4
4.1 Enunciado del Problema 7.4
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
4.2 Formulación Matemática del problema
𝑗=6
𝑀𝑖𝑛 { ∑ 𝐹(𝑃𝑠(𝑗) ) × 4}
𝑗=1
s.a.
𝑃𝑠(𝑗) + 𝑃𝐻(𝑗) = 𝐷𝑗
200 ≤ 𝑃𝑠(𝑗) ≤ 1200
0 ≤ 𝑃𝐻(𝑗) ≤ 2000 (0.9 +
𝑉(𝑗) + 𝑉(𝑗−1)
)
200000
6000 ≤ 𝑉(𝑗) ≤ 18000
𝑉(0) = 10000
𝑉(6) = 10000
𝑉(𝑗) = 𝑉(𝑗−1) + 4 {1000 − (260 + 10𝑃𝐻(𝑗) (1.1 −
𝑉(𝑗) + 𝑉(𝑗−1)
))}
200000
4.3 Modelamiento en Matlab
Se actualizó los nuevos datos en el modelo implementado en 2.3. El código implementado en
Matlab se presenta en el anexo 2.
4.4 Análisis de los resultados obtenidos
A continuación se presenta los resultados obtenidos mediante el método heurístico de
Enjambre de Partículas. En la figura 7, se presenta la potencia generada óptima por cada
planta y para cada período. En la figura 8, se presenta el comportamiento del volumen del
embalse de la planta hidroeléctrica asociada a la operación óptima del sistema.
En la tabla 4, presentamos la función objetivo obtenido con el método PSO. Además se
presenta la generación óptima de las plantas para cada período “j”. Es importante hacer notar
que no necesariamente un método heurístico llega a la solución óptima global del problema,
sin embargo es importante precisar que en el punto 3.4 se obtuvo un resultado bastante
aproximado al óptimo real.
Período
1
2
3
4
Ps
600
1000
545.7
500
PH
0
0
354.3
0
Costo total de Operación (Función Objetivo Óptima)= R/. 89 401
5
200
200
Tabla 4. Generación Óptima de las planta por cada período “j” y Función Objetivo Optima.
6
300
0
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
Resumen de Resultados - Generación
1000
Planta de Vapor (MW )
Planta Hidroeléctrica (MW )
900
800
700
(MW)
600
500
400
300
200
100
0
1
2
3
4
Periodo "j"
5
6
Figura 7 Generación óptima de la planta de vapor y de la planta hidroeléctrica.
1.8
x 10
4
Resumen de Resultados - Volumen(j)
1.6
Volúmen (acre-ft)
1.4
1.2
1
0.8
0.6
0.4
1
2
3
4
Periodo "j"
5
6
7
Figura 8 Volúmen del embalse de la planta hidroeléctrica asociada a la operación óptima del sistema para cada
período “j”.
Suárez Limache José Alberto
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
Anexo 1
Código MATLAB – Ejemplo 7.E
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
v(i,j)=w*v(i,j)+c1*rand()*(pbest(i,j)x(i,j))+c2*rand()*(gbest(1,j)-x(i,j));
end
end
%Actualización de la posición de las "n" particulas.
for i=1:n
for j=1:m
x(i,j)=x(i,j)+v(i,j);
end
end
%Corrección de las posiciones de las "n" particulas en
caso violen algún límite.
for i=1:n
for j=1:m
if x(i,j)<LB(j)
x(i,j)=LB(j);
elseif x(i,j)>UB(j)
x(i,j)=UB(j);
end
end
end
%Calculamos la función objetivo para las nuevas
posiciones de las "n" partículas.
for i=1:n
f(i,1)=ofun(x(i,:));
end
%Se grabará las mejores posiciones por cada partícula.
for i=1:n
if f(i,1)<f0(i,1)
pbest(i,:)=x(i,:); %pbest: mejores posiciones.
f0(i,1)=f(i,1);
%f0: función objetivo de
las mejores posiciones.
end
end
[fmin,index]=min(f0);
%Mejor posición de las mejores
posiciones locales
%Actualizando la mejor posición global
if fmin<fmin0
gbest=pbest(index,:);
fmin0=fmin;
end
%Presentamos los resultados de la mejor partícula.
if ite==1
disp('Iteración
Mejor partícula Función
Objetivo');
end
disp(sprintf('%8g %8g %20g',ite,index,fmin0));
ite=ite+1;
end
%--------------------------------------------------------------------Funcion_Objetivo_Optima(run)=ofun(gbest);
Posicion_Optima(run,:)=gbest;
disp('----------------------------------------')
end
% pso main program-----------------------------------------------------end
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
disp('****************************************************************
**');
disp('Resumen de resultados');
[Funcion_Objetivo_Optima_Final,bestrun]=min(Funcion_Objetivo_Optima);
Funcion_Objetivo_Optima_Final
bestrun
Variables_optimas=Posicion_Optima(bestrun,:)
disp('****************************************************************
**');
toc
%Graficando Volumen del embalse
Resultados(Variables_optimas);
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
function f=ofun(x)
%Esta función calculará la función objetivo de la posición "x"
%Función Obejtivo "Min"
of=(700+(4.8*x(1))+(x(1)*x(1)/2000))*4+(700+(4.8*x(2))+(x(2)*x(2)/2000
))*4+(700+(4.8*x(3))+(x(3)*x(3)/2000))*4+(700+(4.8*x(4))+(x(4)*x(4)/20
00))*4+(700+(4.8*x(5))+(x(5)*x(5)/2000))*4+(700+(4.8*x(6))+(x(6)*x(6)/
2000))*4;
%Calculamos la generación de la planta hidroeléctrica.
Ph1=600-x(1);
Ph2=1000-x(2);
Ph3=900-x(3);
Ph4=500-x(4);
Ph5=400-x(5);
Ph6=300-x(6);
V0=10000;
%Calculamos los volúmenes.
if Ph1==0
V1=V0+4*(1000); %q1=0;
else
V1=V0+4*(1000-(260+10*Ph1)); %q1=260+10*Ph1;
end
if Ph2==0
V2=V1+4*(1000); %q2=0;
else
V2=V1+4*(1000-(260+10*Ph2)); %q2=260+10*Ph2;
end
if Ph3==0
V3=V2+4*(1000); %q3=0;
else
V3=V2+4*(1000-(260+10*Ph3)); %q3=260+10*Ph3;
end
if Ph4==0
V4=V3+4*(1000); %q4=0;
else
V4=V3+4*(1000-(260+10*Ph4)); %q4=260+10*Ph4;
end
if Ph5==0
V5=V4+4*(1000); %q5=0;
else
V5=V4+4*(1000-(260+10*Ph5)); %q5=260+10*Ph5;
end
if Ph6==0
V6=V5+4*(1000); %q6=0;
else
V6=V5+4*(1000-(260+10*Ph6)); %q6=260+10*Ph6;
end
%Restricciones
c0=[];
c0(1)=6000-V1;
c0(2)=V1-18000;
c0(3)=6000-V2;
c0(4)=V2-18000;
c0(5)=6000-V3;
c0(6)=V3-18000;
c0(7)=6000-V4;
c0(8)=V4-18000;
c0(9)=6000-V5;
%
%
%
%
%
%
%
%
%
<=0
<=0
<=0
<=0
<=0
<=0
<=0
<=0
<=0
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
c0(10)=V5-18000;% <=0
c0(11)=9999-V6; % <=0
c0(12)=V6-10001;% <=0
%Creamos la matriz "c" con "ceros" en caso cumpla las restricciones y
"unos" en caso viole las restricciones.
for i=1:length(c0)
if c0(i)>0
c(i)=1;
else
c(i)=0;
end
end
%Penalizamos las violaciones a las restricciones.
penalicacion=100000;
%Función Objetivo de la posición "x"
f=of+penalicacion*sum(c);
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
function f=Resultados(x)
%Graficando reusltados de generación
D=[600 1000 900 500 400 300]; %Actualizar la demanda
for i=1:length(x)
gen(i,1)=x(i);
end
for e=1:length(x)
gen(e,2)=D(e)-x(e);
end
bar(gen,'stack')
title('Resumen de Resultados - Generación')
xlabel('Periodo "j"')
ylabel('(MW)')
legend('Planta de Vapor (MW)','Planta Hidroeléctrica (MW)')
grid on
%Graficando Volúmenes
V(1)=10000;
%Calculamos los volúmenes.
for i=1:length(x)
if gen(i,2)==0
V(i+1)=V(i)+4*(1000); %q1=0;
else
V(i+1)=V(i)+4*(1000-(260+10*gen(i,2))); %q1=260+10*Ph1;
end
end
figure
plot(V)
title('Resumen de Resultados - Volumen(j)')
xlabel('Periodo "j"')
ylabel('Volúmen (acre-ft)')
grid on
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
Anexo 2
Código MATLAB – Problema 7.4
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
v(i,j)=w*v(i,j)+c1*rand()*(pbest(i,j)x(i,j))+c2*rand()*(gbest(1,j)-x(i,j));
end
end
%Actualización de la posición de las "n" particulas.
for i=1:n
for j=1:m
x(i,j)=x(i,j)+v(i,j);
end
end
%Corrección de las posiciones de las "n" particulas en
caso violen algún límite.
for i=1:n
for j=1:m
if x(i,j)<LB(j)
x(i,j)=LB(j);
elseif x(i,j)>UB(j)
x(i,j)=UB(j);
end
end
end
%Calculamos la función objetivo para las nuevas
posiciones de las "n" partículas.
for i=1:n
f(i,1)=ofun(x(i,:));
end
%Se grabará las mejores posiciones por cada partícula.
for i=1:n
if f(i,1)<f0(i,1)
pbest(i,:)=x(i,:); %pbest: mejores posiciones.
f0(i,1)=f(i,1);
%f0: función objetivo de
las mejores posiciones.
end
end
[fmin,index]=min(f0);
%Mejor posición de las mejores
posiciones locales
%Actualizando la mejor posición global
if fmin<fmin0
gbest=pbest(index,:);
fmin0=fmin;
end
%Presentamos los resultados de la mejor partícula.
if ite==1
disp('Iteración
Mejor partícula Función
Objetivo');
end
disp(sprintf('%8g %8g %20g',ite,index,fmin0));
ite=ite+1;
end
%--------------------------------------------------------------------Funcion_Objetivo_Optima(run)=ofun(gbest);
Posicion_Optima(run,:)=gbest;
disp('----------------------------------------')
end
% pso main program-----------------------------------------------------end
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
disp('****************************************************************
**');
disp('Resumen de resultados');
[Funcion_Objetivo_Optima_Final,bestrun]=min(Funcion_Objetivo_Optima);
Funcion_Objetivo_Optima_Final
bestrun
Variables_optimas=Posicion_Optima(bestrun,:)
disp('****************************************************************
**');
toc
%Graficando Volumen del embalse
Resultados(Variables_optimas);
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
function f=ofun(x)
%Esta función calculará la función objetivo de la posición "x"
%Función Obejtivo "Min"
of=(770+(5.28*x(1))+(x(1)*x(1)*0.00055))*4+(770+(5.28*x(2))+(x(2)*x(2)
*0.00055))*4+(770+(5.28*x(3))+(x(3)*x(3)*0.00055))*4+(770+(5.28*x(4))+
(x(4)*x(4)*0.00055))*4+(770+(5.28*x(5))+(x(5)*x(5)*0.00055))*4+(770+(5
.28*x(6))+(x(6)*x(6)*0.00055))*4;
%Calculamos la generación de la planta hidroeléctrica.
Ph1=600-x(1);
Ph2=1000-x(2);
Ph3=900-x(3);
Ph4=500-x(4);
Ph5=400-x(5);
Ph6=300-x(6);
V0=10000;
%Calculamos los volúmenes.
if Ph1==0
V1=V0+4*(1000); %q1=0;
else
V1=(5000*V0+14800000-220000*Ph1+Ph1*V0)/(5000-Ph1);
%q1=260+10*Ph1*(1.1-(V0+V1)/200000);
end
if Ph2==0
V2=V1+4*(1000); %q2=0;
else
V2=(5000*V1+14800000-220000*Ph2+Ph2*V1)/(5000-Ph2);
%q2=260+10*Ph2*(1.1-(V1+V2)/200000);
end
if Ph3==0
V3=V2+4*(1000); %q3=0;
else
V3=(5000*V2+14800000-220000*Ph3+Ph3*V2)/(5000-Ph3);
%q3=260+10*Ph3*(1.1-(V2+V3)/200000);
end
if Ph4==0
V4=V3+4*(1000); %q4=0;
else
V4=(5000*V3+14800000-220000*Ph4+Ph4*V3)/(5000-Ph4);
%q4=260+10*Ph4*(1.1-(V3+V4)/200000);
end
if Ph5==0
V5=V4+4*(1000); %q5=0;
else
V5=(5000*V4+14800000-220000*Ph5+Ph5*V4)/(5000-Ph5);
%q5=260+10*Ph5*(1.1-(V4+V5)/200000);
end
if Ph6==0
V6=V5+4*(1000); %q6=0;
else
V6=(5000*V5+14800000-220000*Ph6+Ph6*V5)/(5000-Ph6);
%q6=260+10*Ph6*(1.1-(V5+V6)/200000);
end
%Restricciones
c0=[];
c0(1)=6000-V1;
% <=0
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
c0(2)=V1-18000; % <=0
c0(3)=6000-V2; % <=0
c0(4)=V2-18000; % <=0
c0(5)=6000-V3; % <=0
c0(6)=V3-18000; % <=0
c0(7)=6000-V4; % <=0
c0(8)=V4-18000; % <=0
c0(9)=6000-V5; % <=0
c0(10)=V5-18000;% <=0
c0(11)=10000-V6; % <=0 No es necesario modelar una restricción de
igualdad
%Creamos la matriz "c" con "ceros" en caso cumpla las restricciones y
"unos" en caso viole las restricciones.
for i=1:length(c0)
if c0(i)>0
c(i)=1;
else
c(i)=0;
end
end
%Penalizamos las violaciones a las restricciones.
penalicacion=100000;
%Función Objetivo de la posición "x"
f=of+penalicacion*sum(c);
TÓPICOS EN SISTEMAS EXPERTOS E INTELIGENCIA ARTIFICIAL
function f=Resultados(x)
%Graficando reusltados de generación
D=[600 1000 900 500 400 300]; %Actualizar la demanda
for i=1:length(x)
gen(i,1)=x(i);
end
for e=1:length(x)
gen(e,2)=D(e)-x(e);
end
bar(gen,'stack')
title('Resumen de Resultados - Generación')
xlabel('Periodo "j"')
ylabel('(MW)')
legend('Planta de Vapor (MW)','Planta Hidroeléctrica (MW)')
grid on
%Graficando Volúmenes
V(1)=10000;
%Calculamos los volúmenes.
for i=1:length(x)
if gen(i,2)==0
V(i+1)=V(i)+4*(1000); %q1=0;
else
V(i+1)=V(i)+4*(1000-(260+10*gen(i,2))); %q1=260+10*Ph1;
end
end
figure
plot(V)
title('Resumen de Resultados - Volumen(j)')
xlabel('Periodo "j"')
ylabel('Volúmen (acre-ft)')
grid on