Universidad Nacional Autónoma de México
Universidad Nacional Autónoma de México
TESIS
QUE PARA OPTAR POR EL GRADO DE:
DOCTOR EN INGENIERÍA CIVIL
PRESENTA:
OMAR ANTONIO DE LA CRUZ COURTOIS
TUTORES PRINCIPALES
DR. RAMÓN DOMÍNGUEZ MORA
INSTITUTO DE INGENIERÍA
DERECHOS RESERVADOS ©
PROHIBIDA SU REPRODUCCIÓN TOTAL O PARCIAL
Todo el material contenido en esta tesis esta protegido por la Ley Federal
del Derecho de Autor (LFDA) de los Estados Unidos Mexicanos (México).
Índice de figuras
1.1. Capacidad hidroeléctrica hasta 1987 [Garcia, 1992] . . . . . . . . . . . . . 36
1.2. Comparativo entre hidroléctricas y termoeléctricas hasta 1987 [Garcia, 1992] 36
1.3. Región factible en la programación lineal . . . . . . . . . . . . . . . . . . . 38
1.4. Representación de los estados de una presa . . . . . . . . . . . . . . . . . 41
1.5. Ejemplo de una polı́tica óptima de dos presas . . . . . . . . . . . . . . . . . 49
2.6. Esquema de una variable aleatoria . . . . . . . . . . . . . . . . . . . . . . . 57
2.7. Función de distribución continua y discreta . . . . . . . . . . . . . . . . . . 60
2.8. Área bajo la curva y la función de densidad . . . . . . . . . . . . . . . . . . 61
2.9. Partición de un intervalo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
2.10.Suma de Riemann-Stieltjes . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
3.11.Regiones Hidrológicas de México I-VI [CFE, 2009] . . . . . . . . . . . . . . 93
3.12.Regiones Hidrológicas de México VII-XIII [CFE, 2009] . . . . . . . . . . . . 94
3.13.Mapa de Regiones Hidrológicas de México [CONAGUA, 2009] . . . . . . . 95
3.14.Superficie por estados de la Región Hidrológica 30 [CONAGUA, 1994] . . . 95
3.15.Superficie por estados de la Región Hidrológica 30 [INEGI, 2007] . . . . . . 96
3.16.Región Hidrológica 30 [INEGI, 2009] . . . . . . . . . . . . . . . . . . . . . . 96
3.17.Sistema hidroeléctrico de la cuenca del rı́o Grijalva [CFE, 1985] . . . . . . 99
3.18.Perfil del sistema hidroeléctrico Grijalva [Arganis, 2004] . . . . . . . . . . . 100
3.19.Sistema hidroeléctrico de la cuenca del rı́o Grijalva [INEGI, 2009] . . . . . . 101
3.20.Datos de la cortina de la Presa Malpaso . . . . . . . . . . . . . . . . . . . . 103
3.21.Datos del embalse de la Presa Malpaso . . . . . . . . . . . . . . . . . . . . 103
3.22.Datos de la casa de máquinas de la Presa Malpaso . . . . . . . . . . . . . 103
3.23.Datos de la potencia y generación de la Presa Malpaso . . . . . . . . . . . 104
3.24.Datos de la cortina de la Presa La Angostura . . . . . . . . . . . . . . . . . 105
3.25.Datos del embalse de la Presa La Angostura . . . . . . . . . . . . . . . . . 105
3.26.Datos de la casa de máquinas de la Presa La Angostura . . . . . . . . . . 105
3.27.Datos de la potencia y generación de la Presa La Angostura . . . . . . . . 106
3.28.Datos de la cortina de la Presa Chicoasén . . . . . . . . . . . . . . . . . . . 109
3.29.Datos del embalse de la Presa Chicoasén . . . . . . . . . . . . . . . . . . . 109
3.30.Datos de la casa de máquinas de la Presa Chicoasén . . . . . . . . . . . . 109
3.31.Datos de la potencia y generación de la Presa Chicoasén . . . . . . . . . . 109
3.32.Datos de la cortina de la Presa Peñitas . . . . . . . . . . . . . . . . . . . . 111
3.33.Datos del embalse de la Presa Peñitas . . . . . . . . . . . . . . . . . . . . . 111
3.34.Datos de la casa de máquinas de la Presa Peñitas . . . . . . . . . . . . . . 112
3.35.Datos de la potencia y generación de la Presa Peñitas . . . . . . . . . . . . 112
3.36.Hidrologı́a de todas las presas [Arellano, 1997] . . . . . . . . . . . . . . . . 112
3.37.Datos del embalse de todas las presas [Arellano, 1997] . . . . . . . . . . . 113
3.38.Datos de la obra de toma de todas las presas [Arellano, 1997] . . . . . . . 114
3.39.Datos de la cortina de todas las presas [Arellano, 1997] . . . . . . . . . . . 114
3.40.Datos del vertedor de todas las presas [Arellano, 1997] . . . . . . . . . . . 115
3.41.Impacto de los desastres naturales en el PIB estatal de Chiapas [INEGI, 2007]116
3.42.Eventos extremos hidrometeorológicos de 1995 a 1999 . . . . . . . . . . . 117
ÍNDICE DE FIGURAS 4
5.109.Ingresos totales para los meses de Junio y Julio. Años 1952-1990 . . . . . 169
5.110.Ingresos totales para los meses de Junio y Julio. Años 1991-2017 . . . . . 169
5.111.Tabla de frecuencias de la Etapa 5 . . . . . . . . . . . . . . . . . . . . . . . 170
5.112.Histograma de la Etapa 5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 170
5.113.Tabla de probabilidades de la Etapa 5 . . . . . . . . . . . . . . . . . . . . . 171
5.114.Ingresos totales para los meses de Enero a Mayo. Años 1952-1971 . . . . 171
5.115.Ingresos totales para los meses de Enero a Mayo. Años 1972-1990 . . . . 172
5.116.Ingresos totales para los meses de Enero a Mayo. Años 1991-2003 . . . . 172
5.117.Ingresos totales para los meses de Enero a Mayo. Años 2004-2017 . . . . 173
5.118.Tabla de frecuencias de la Etapa 6 . . . . . . . . . . . . . . . . . . . . . . . 173
5.119.Histograma de la Etapa 6 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 174
5.120.Tabla de probabilidades de la Etapa 6 . . . . . . . . . . . . . . . . . . . . . 174
5.121.Ingresos totales para los meses de noviembre y diciembre. Años 1952-1990 175
5.122.Ingresos totales para los meses de noviembre y diciembre. Años 1991-2017 175
5.123.Tabla de frecuencias de la Etapa 1 . . . . . . . . . . . . . . . . . . . . . . . 176
5.124.Histograma de la Etapa 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 176
5.125.Tabla de probabilidades de la Etapa 1 . . . . . . . . . . . . . . . . . . . . . 177
5.126.Ingresos totales para el mes de octubre. Años 1952-2017 . . . . . . . . . . 177
5.127.Tabla de frecuencias de la Etapa 2 . . . . . . . . . . . . . . . . . . . . . . . 178
5.128.Histograma de la Etapa 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 178
5.129.Tabla de probabilidades de la Etapa 2 . . . . . . . . . . . . . . . . . . . . . 179
5.130.Ingresos totales para el mes de septiembre. Años 1952-2017 . . . . . . . . 179
5.131.Tabla de frecuencias de la Etapa 3 . . . . . . . . . . . . . . . . . . . . . . . 180
5.132.Histograma de la Etapa 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 180
5.133.Tabla de probabilidades de la Etapa 3 . . . . . . . . . . . . . . . . . . . . . 181
5.134.Ingresos totales para el mes de agosto. Años 1952-2017 . . . . . . . . . . 181
5.135.Tabla de frecuencias de la Etapa 4 . . . . . . . . . . . . . . . . . . . . . . . 182
5.136.Histograma de la Etapa 4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 182
5.137.Tabla de probabilidades de la Etapa 4 . . . . . . . . . . . . . . . . . . . . . 183
5.138.Ingresos totales para los meses de Junio y Julio. Años 1952-1990 . . . . . 183
5.139.Ingresos totales para los meses de Junio y Julio. Años 1991-2017 . . . . . 184
5.140.Tabla de frecuencias de la Etapa 5 . . . . . . . . . . . . . . . . . . . . . . . 184
5.141.Histograma de la Etapa 5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 185
5.142.Tabla de probabilidades de la Etapa 5 . . . . . . . . . . . . . . . . . . . . . 185
5.143.Ingresos totales para los meses de Enero a Mayo. Años 1952-1971 . . . . 186
5.144.Ingresos totales para los meses de Enero a Mayo. Años 1972-1990 . . . . 186
5.145.Ingresos totales para los meses de Enero a Mayo. Años 1991-2003 . . . . 187
5.146.Ingresos totales para los meses de Enero a Mayo. Años 2004-2017 . . . . 187
5.147.Tabla de frecuencias de la Etapa 6 . . . . . . . . . . . . . . . . . . . . . . . 188
5.148.Histograma de la Etapa 6 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 188
5.149.Tabla de probabilidades de la Etapa 6 . . . . . . . . . . . . . . . . . . . . . 189
5.150.Escurrimientos históricos del Grijalva ([CFE, 2009]) . . . . . . . . . . . . . 189
5.151.Trayectoria del huracán Roxanne ([WU, 2018]) . . . . . . . . . . . . . . . . 190
5.152.Trayectoria del huracán Opal ([WU, 2018]) . . . . . . . . . . . . . . . . . . . 191
5.153.Trayectoria de la tormenta Tropical Frances [WU, 2018] . . . . . . . . . . . 191
ÍNDICE DE FIGURAS 7
Abreviaciones
Notación
Índice
1. Introducción 19
1.1. Descripción del problema . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
1.2. Justificación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
1.3. Hipótesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
1.4. Antecedentes históricos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
1.5. Terminologı́a previa . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
1.5.1. Aspectos hidrológicos . . . . . . . . . . . . . . . . . . . . . . . . . . 27
1.5.2. Almacenamiento en vasos . . . . . . . . . . . . . . . . . . . . . . . 29
1.5.3. Potencia y energı́a . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
1.6. Métodos de optimización . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
1.6.1. Programación lineal y no lineal . . . . . . . . . . . . . . . . . . . . . 37
1.6.2. Método de Lloyd y Morán . . . . . . . . . . . . . . . . . . . . . . . . 39
1.6.3. Valor esperado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42
1.6.4. Optimización multiobjetivo . . . . . . . . . . . . . . . . . . . . . . . . 43
1.6.5. Principio del máximo de Pontryagin . . . . . . . . . . . . . . . . . . 44
1.6.6. Programación dinámica . . . . . . . . . . . . . . . . . . . . . . . . . 45
1.6.7. Método de Berezowsky-Domı́nguez-Fuentes . . . . . . . . . . . . . 46
1.6.8. Programación dinámica estocástica . . . . . . . . . . . . . . . . . . 47
1.6.9. Método de Domı́nguez . . . . . . . . . . . . . . . . . . . . . . . . . 47
1.7. Objetivos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
1.7.1. Objetivo general . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
1.7.2. Objetivos especı́ficos . . . . . . . . . . . . . . . . . . . . . . . . . . 49
2. Modelación estocástica 51
2.1. Nociones de probabilidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
2.1.1. La σ-álgebra . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
2.1.2. Espacios de probabilidad . . . . . . . . . . . . . . . . . . . . . . . . 52
2.1.3. Probabilidad condicional y eventos independientes . . . . . . . . . . 54
2.1.4. Variables aleatorias . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
2.1.5. Función de distribución . . . . . . . . . . . . . . . . . . . . . . . . . 59
2.1.6. Función de densidad y probabilidad . . . . . . . . . . . . . . . . . . 60
2.1.7. Integral de Riemann-Stieltjes . . . . . . . . . . . . . . . . . . . . . . 62
2.1.8. Esperanza matemática y varianza . . . . . . . . . . . . . . . . . . . 64
2.2. Procesos estocásticos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
2.3. Cadenas de Markov y probabilidades de transición . . . . . . . . . . . . . . 69
2.4. Programación dinámica . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
2.5. El algoritmo de la programación dinámica . . . . . . . . . . . . . . . . . . . 72
2.6. Control estocástico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
2.7. Modelo de control markoviano . . . . . . . . . . . . . . . . . . . . . . . . . 74
2.7.1. Espacios métricos, separables y completos . . . . . . . . . . . . . . 74
2.7.2. Polı́ticas de control . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
2.7.3. El modelo del control descontado con restricciones . . . . . . . . . 79
ÍNDICE 15
4. Metodologı́a 128
4.1. Cálculo de probabilidades . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132
4.1.1. Registros históricos . . . . . . . . . . . . . . . . . . . . . . . . . . . 132
4.1.2. Histogramas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132
4.2. Cálculo de extracciones máximas y restricciones . . . . . . . . . . . . . . . 133
4.2.1. Extracciones máximas . . . . . . . . . . . . . . . . . . . . . . . . . . 133
4.2.2. Restricciones . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133
4.3. Modelo de Markov . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
4.3.1. Espacio de estados . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
4.3.2. Espacio de acciones . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
4.3.3. Espacio de acciones admisibles . . . . . . . . . . . . . . . . . . . . 134
4.3.4. Kernel de transición . . . . . . . . . . . . . . . . . . . . . . . . . . . 134
4.3.5. Funciones de ganancia y costo . . . . . . . . . . . . . . . . . . . . . 135
4.4. Cálculo de beneficios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136
ÍNDICE 16
5. Resultados 141
5.1. Registros históricos de escurrimientos . . . . . . . . . . . . . . . . . . . . . 141
5.1.1. Registros históricos de La Angostura . . . . . . . . . . . . . . . . . 141
5.1.2. Registros históricos de Malpaso . . . . . . . . . . . . . . . . . . . . 147
5.2. Probabilidades de ingreso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 158
5.2.1. Etapa 1. La Angostura . . . . . . . . . . . . . . . . . . . . . . . . . . 158
5.2.2. Etapa 2. La Angostura . . . . . . . . . . . . . . . . . . . . . . . . . . 162
5.2.3. Etapa 3. La Angostura . . . . . . . . . . . . . . . . . . . . . . . . . . 164
5.2.4. Etapa 4. La Angostura . . . . . . . . . . . . . . . . . . . . . . . . . . 167
5.2.5. Etapa 5. La Angostura . . . . . . . . . . . . . . . . . . . . . . . . . . 169
5.2.6. Etapa 6. La Angostura . . . . . . . . . . . . . . . . . . . . . . . . . . 171
5.2.7. Etapa 1. Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 175
5.2.8. Etapa 2. Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 177
5.2.9. Etapa 3. Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 179
5.2.10. Etapa 4. Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 181
5.2.11. Etapa 5. Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 183
5.2.12. Etapa 6. Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . . . . 186
5.2.13. Eventos extremos . . . . . . . . . . . . . . . . . . . . . . . . . . . . 189
5.3. Extracción mensual máxima . . . . . . . . . . . . . . . . . . . . . . . . . . . 198
5.4. Curvas guı́a . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 200
5.4.1. Curva guı́a de La Angostura . . . . . . . . . . . . . . . . . . . . . . 201
5.4.2. Curva guı́a de Malpaso . . . . . . . . . . . . . . . . . . . . . . . . . 202
5.4.3. Curva guı́a de La Angostura, quincena 1-6 . . . . . . . . . . . . . . 203
5.4.4. Curva guı́a de La Angostura, quincenas 7-12 . . . . . . . . . . . . . 204
5.4.5. Curva guı́a de La Angostura, quincenas 13-24 . . . . . . . . . . . . 205
5.4.6. Curva guı́a de Malpaso, quincenas 1-12 . . . . . . . . . . . . . . . . 206
5.4.7. Curva guı́a de Malpaso, quicenas 13-20 . . . . . . . . . . . . . . . . 207
5.4.8. Curva guı́a de Malpaso, quincenas 21-24 . . . . . . . . . . . . . . . 209
5.5. Polı́ticas de control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 210
5.5.1. Espacio de estados . . . . . . . . . . . . . . . . . . . . . . . . . . . 210
5.5.2. Espacio de acciones . . . . . . . . . . . . . . . . . . . . . . . . . . . 212
5.5.3. Espacio de acciones admisibles . . . . . . . . . . . . . . . . . . . . 213
5.5.4. Kernel de transición . . . . . . . . . . . . . . . . . . . . . . . . . . . 214
5.5.5. Función de ganancia y función de costo . . . . . . . . . . . . . . . . 216
5.5.6. Funciones de restricción . . . . . . . . . . . . . . . . . . . . . . . . . 219
ÍNDICE 17
RESUMEN
En el presente trabajo se aplicó la teorı́a de los modelos de control markovianos a
fin de obtener polı́ticas de operación óptima para un sistema de presas que funcionan
en serie. Para ello se tomó en cuenta los registros hidrológicos históricos, lo que permitió
maximizar la energı́a electrica que se genera, ası́ como minimizar los derrames y déficits.
ABSTRACT
In the present work the theory of Markovian control models was applied in order to
obtain optimal operation policies for a system of dams that works in series. For this, histo-
rical hydrological records were taken into account, which allowed to maximize the electric
energy that is generated, as well as to minimize spills and deficits.
1 INTRODUCCIÓN 19
1. Introducción
1.1. Descripción del problema
A medida que crece la población de una región, las demandas de recursos hidráuli-
cos aumenta, y con ello las demandas de energı́a eléctrica, es por ello que se necesita
revisar continuamente la gestión y manejo de la generación eléctrica de una zona dada.
En nuestro paı́s un porcentaje de esta generación se debe a los embalses hidroeléctri-
cos, con lo cual se debe estudiar la forma en la que se genere energı́a hidroeléctrica en
función del nivel de almacenamiento del vaso. Con esto en mente se necesita resolver
dos principales problemas: satisfacer las demandas de energı́a de una región y evitar
que surjan problemas de inundaciones aguas abajo del embalse.
Actualmente, poder resolver los problemas anteriores conlleva a una diversidad de es-
tudios, desde ambientales, hasta sociales, y por supuesto, los técnicos. Es por ello que
el manejo actual de un embalse difiere del manejo para el cual fue proyectado en el pa-
sado, con lo cual se deben actualizar las formas de gestionar estos recursos. Los niveles
de un embalse hidroeléctrico serán los que determinarán en gran medida la distribución
de la energı́a y la posibilidad de desbordes. Estos niveles son los que se conocen como
las polı́ticas de operación en el embalse. Los operadores de las presas hidroeléctricas
son los encargados de mantener estos niveles de tal manera que cumplan los requisitos
mencionados anteriormente, y esto se realiza a través de polı́ticas previas que la CFE y
actualmente el Instituto de Ingenierı́a de la UNAM proponen.
b) Reducir la probabilidad de los derrames por el vertedor que pueden causar inunda-
ciones aguas abajo.
Otro de los factores que presentan los sistemas hidroeléctricos es el hecho que los
volúmenes de ingreso son aleatorios y en ocasiones presentan correlaciones, depen-
diendo la época del año. Los objetivos presentados anteriormente tienen una estrecha
interdependencia y para resolver todos se deben tratar al mismo tiempo. En otras pala-
bras, para maximizar la energı́a se necesita mayor nivel en el embalse, lo que provoca
mayor probabilidad en derrames por el vertedor, y por el contrario, garantizar una seguri-
dad en los niveles para evitar inundaciones, conlleva a que las demandas de energı́a no
se satisfagan por tener los niveles muy bajos.
Desde las últimas dos décadas la Comisión Federal de Electricidad (CFE) ha ido evo-
lucionando en la determinación de sus polı́ticas óptimas. Antes se realizaba una polı́tica
de operación quincenal de predespacho, empı́rica en su mayorı́a, restringida por la Co-
misión Nacional del Agua (CONAGUA), donde se establecı́an únicamente los niveles
máximos por cada época del año para evitar inundaciones. Estos niveles limitantes oca-
sionaron pérdidas con algunos centenares de millones de pesos anuales, por el hecho de
ser demasiado conservadoras. La CFE reevaluó los hechos y optó por métodos analı́ti-
cos, motivando y orientando a los operadores a su uso y aplicación. Los modelos técnicos
se han ido desarrollando en conjunto con la experiencia de los operadores, formando un
lenguaje en conjunto matemático-empı́rico para tener una terminologı́a afı́n tanto a los
investigadores como a los operadores existentes.
Es aquı́ donde surgen las técnicas de programación dinámica que más adelante se
explicarán. En este trabajo, el nivel de complejidad del problema se adecuó para modelar
matemáticamente la situación actual. Se realizaron simulaciones numéricas para saber
qué ocurrirı́a en diversos escenarios, para calibrar los mismos modelos matemáticos pro-
puestos. Los registros históricos hidrológicos jugaron un papel fundamental para tener un
completo panorama y preveer lo más exacto posible situaciones extremas de manera que
se pueda actuar en la forma más eficiente. La aplicación se hizo en presas del Sistema
Grijalva.
pretende contribuir a una revisión de las polı́ticas actuales mediante nuevas técnicas de
optimización a fin de reducir tiempos de cálculo que se traducen en costos indirectos y
proponer nuevas matrices de polı́ticas óptimas.
1.2. Justificación
Las dificultades en mantener un nivel deseado en el embalse de una presa se pue-
den clasificar en dos tipos, a corto y a largo plazo. A corto plazo se presentan cuando
existen demandas inusuales o excesivas en intervalos cortos de tiempo. A largo plazo se
presentan en determinar qué cantidad de agua se debe almacenar para los siguientes
dı́as, o semanas o meses, cuando se presentan épocas de sequı́as o de lluvias escasas.
El hecho de contar con energı́a disponible al momento se debe balancear con pre-
decir cuánta energı́a se requerirá en un futuro cercano o lejano, tomando en cuenta los
futuros escurrimientos o incorporaciones de agua al embalse. Además la generación de
energı́a y la regulación de los almacenamientos están restringidos por los distintos ele-
mentos operativos de un sistema hidroeléctrico. Si a esto se le agrega que un sistema
puede operar en forma múltiple, el número de variables a considerar aumenta, y ahora
se debe tener en cuenta los demás embalses que se ven afectados por una decisión to-
mada aguas arriba. Finalmente los picos ocasionados por demandas de energı́a también
conllevan a tener reservas de energı́a que se traduce en reservas de agua para satisfacer
los requerimientos solicitados.
Debe entonces existir un equilibrio entre maximizar la generación energética pero sin
poner en riesgo a los poblados que viven en las márgenes del sistema hidroelétrico, sobre
todo cuando el ingreso del agua dependerá de factores climatológicos que no dependen
del hombre, y que su posible “control”se encontrará limitado por el tiempo.
En este sentido, disponer de una polı́tica que determine el nivel de un embalse, estará
ligado a la disposición de la información de las precipitaciones cercanas a la zona, ya que
los niveles de aguas se verán afectados, no sólo por la cantidad de agua que llueva en
el lugar, sino por el agua que se junte en la cuenca propia a través de escurrimientos.
Esta información dependerá de las mediciones meteorológicas directas o indirectas que
se tengan, nuevamente del lugar de estudio.
cible del clima con respecto al tiempo, y el modelo que resuelva el problema deberá ser
capaz entonces de cubrir tanto los requerimientos energéticos como los requerimientos
de seguridad operacional.
1.3. Hipótesis
La modelación estocástica permite determinar las polı́ticas óptimas en los niveles de
un sistema hidroeléctrico que funciona en cascada, disminuyendo la incertidumbre que
se genera cuando se suponen variables continuas.
en concreto polı́ticas de operación para embalses del estado de Chiapas. Arnold y Sim-
mons (1988), Palmer y Holmes (1988), Savic y Simonovic (1991), Srinivasan y Engel
(1994), Ford y Killen (1995), Liang (1996) y Shresta (1996) desarrollaron un modelo tri-
dimensional, a partir de la componente hidrológica, la componente de optimización y la
componente de decisión; utilizando la teorı́a de decisiones para resultados racionales
y reales. Kitanidis y Foufoua-Georgiu (1987) trabajaron en embalses múltiples y propu-
sieron una metodologı́a a partir de análisis de errores y programación dinámica gradiente.
Domı́nguez et. al (1993) comenzaron a estudiar el Sistema Hidroeléctrico del rı́o Gri-
jalva, dividiéndolo en dos partes: polı́ticas de regulación de avenidas para el uso del
vertedor y polı́ticas de los niveles de embalse como un volumen útil de generación de
energı́a eléctrica; dichas polı́ticas se determinaron en forma mensual maximizando la
generación eléctrica y minimizando los derrames y posibles déficits. Johnson, Stedinger,
Schoemaker, Li, y Tejada-Guibert (1993) propusieron nuevos algoritmos de simulación
numérica para programación dinámica en estados continuos a través de interpolación
lineal. Avilés (1994) mediante la ingenierı́a eléctrica realizó trabajos de optimización en
lı́nea. Eschenbach, Schoemaker y Caffey (1995) crearon algoritmos paralelos para pro-
gramación dinámica estocástica en estados continuos y variables de control. Jacobs,
Freeman, Grygier, Morton, Schultz, Staschus y Stedinger (1995) realizaron un esquema
de programación para la generación hidroeléctrica bajo condiciones de incertidumbre.
1 INTRODUCCIÓN 24
Belaineh (1999) toma en cuenta el uso de agua superficial y subterránea para generar
un modelo de optimización de polı́ticas de operación lineales, proveyendo a los usuarios
la simulación de los caudales formados en el embalse. ReVelle (1999) obtiene polı́ticas
de operación en un embalse y embalses en paralelo por medio de modelos determinı́sti-
cos de programación lineal; propone modelos de simulación desde el disenõ hasta la
operación del sistema hidroeléctrico. Huang Wen-Chen (1999) trabajó con sistemas de
apoyo de decisiones (DSS por sus siglas en inglés) para la operación de embalses a
largo plazo y lo aplicaron en Taiwán, commplementando los trabajos de los años ochenta
anteriores.
Sánchez y Wagner (2003) determinaron reglas de operación óptima para dos em-
balses utilizando un algoritmo genético en el Instituto Mexicano de Tecnologı́a del Agua.
Sánchez y Wagner (2004) modelaron numéricamente la operación óptima de un hidrosis-
tema de aguas superficiales. Ailing Li (2004) trabajó en un estudio a gran escala mediante
un método de sistema de descomposición para determinar polı́ticas óptimas de un siste-
ma hidroeléctrico. Empleó en sistemas de presas en series en China, obteniendo buenos
resultados en China.
Markus (2006) utilizó la programación estocástica para econtrar polı́ticas óptimas mi-
nimizando los costos de producción. En su modelo incluyó tanto la carga de la energı́a
como las entradas diversas de agua y utilizó herramientas de descomposición jerarqui-
zado, llamada método de aproximación media de la muestra. Domı́nguez et. al (2006)
generaron nuevas polı́ticas de operación que toman los eventos hidrológicos del 2005.
1 INTRODUCCIÓN 26
Ngo (2006) empleó métodos numéricos para optimizar variables de decisión para simu-
lación de polı́ticas de embalses. Las aplicaciones fueron hechas en Vietnam empleando
software especializado y modelos de flujo de rı́os.
Deepti y Maria (2010) realizaron una encuesta sobre diferentes modelos de optimi-
zación y simulación de embalses. Además incorporaron algoritmos basados tanto en la
simulación como en la modelación en conjunto, tales como programación dinámica, pro-
gramación no lineal, algoritmos genéticos, lógica difusa y computación evolutiva.
Vijendra y Yadab (2018) desarrollaron los algoritmos llamados TLBO (Teaching Lear-
ned Based Optimization) y JA (Jaya Algorithm) para trabajar múltiples embalses con-
trolando las variables poblacionales, demandas y climatológicas. Propusieron además
algoritmos computacionales para incorporar los cálculos en presas hidroeléctricas que
operan en cascada.
d) Tunel de desfogue. Conducto que transporta el agua desde la obra de toma hasta
las turbinas.
a) Parteaguas. Lı́nea imaginaria formada por los puntos de mayor nivel topográfico y
que la delimita con otras cuencas.
Para indicar el grado de respuesta de una cuenca se cuenta con la pendiente del cauce
principal. Como la pendiente varı́a a lo largo del cauce, generalmente dependiendo del
autor, se define la pendiente media, como el desnivel entre los extremos de la corriente
dividido entre su longitud media en planta. Existen diversas expresiones y metodologı́as
para el cálculo de dichas pendientes [Naghettini, 2017].
Este escurrimiento, fluye sobre el terreno donde se sigue infiltrando, e incluso se eva-
pora en pequeñas cantidades, y se produce mientras el agua no llegue a cauces bien
definidos. Cuando llega al cauce principal se convierte en escurrimiento en corrientes. El
flujo sobre el terreno y el escurrimiento en corrientes, es lo que forma el escurrimiento
superficial, a diferencia del escurrimiento subsuperficial, que se ocasiona cuando hay in-
filtraciones, y fluye paralelo al escurrimiento anterior. Si por el contrario, se infiltra hasta
niveles inferiores al freático, se denomina escurrimiento subterráneo.
Los picos pueden variar desde algunos litros por segundo hasta miles de metros cúbi-
cos por segundo, y el tiempo base puede durar desde minutos hasta algunos dı́as. El área
bajo el hidrograma representa el volumen total escurrido y se expresa
Z t1
V = Qdt. (1.1)
t0
Para fines de diseño conviene separar el gasto base del gasto directo, para lo cual existen
diversos métodos, desde trazar lı́neas horizontales aproximadas, hasta algunos métodos
analı́tidos y expresiones de tipo exponencial que involucran el área de la cuenca y tiem-
pos de vaciado de la misma.
1 INTRODUCCIÓN 29
Para obtener fı́sicamente los escurrimientos, existen los métodos de aforo. Aforar una
corriente significa determinar a través de mediciones el gasto que pasará por una sección
dada. En México se utilizan las secciones de control, las relaciones sección pendiente y
la relación sección velocidad, para determinar el aforo. Una sección de control es una
sección donde existe una relación entre el tirante y el gasto, entre las más comunes son
las de tirante crı́tico y los vertedores. Una relación sección-pendiente se emplea para
determinar el gasto máximo que se presentó durante una avenida reciente en un cauce.
Una relación sección-velocidad mide las velocidades en distintos puntos de una sec-
ción transversal para luego calcular el gasto mediante la ecuación de continuidad.
cide con la cresta o punto más alto del vertedor. El vertedor funciona para desalojar los
volúmenes excedentes de agua que puedan poner en peligro la seguridad de la obra, y
por ello la seguridad de los poblados aguas abajo de la central. Existen presas donde los
niveles del NAMO son variables a lo largo del año ya que las descargas de vertedores
son controladas por compuertas. El volumen de agua que se almacena entre el NAMO y
el NAMINO se conoce como capacidad útil, y con este se satisfacen las demandas de
agua solicitadas.
El nivel más alto que debe alcanzar el agua en el vaso bajo cualquier situación es el
NAME; el volumen de agua que queda entre el NAMO y el NAME, se le llama superal-
macenamiento y funciona para controlar las avenidas que se presentan cuando el nivel
en el vaso está cercano al NAMO. Cuando existe la marea producida por el viento en la
presa, ası́ como el oleaje, este es controlado por el bordo libre, el cual es el espacio que
se localiza entre la corona (máxima elevación de la cortina) y el NAME. De igual forma,
en el caso de asentamientos de la presa, esta diferencia funciona como un sistema de
seguridad.
Para las estimaciones del volumen útil existen diversos métodos como la curva masa
y el algoritmo de pico secuente [Naghettini, 2017]. El hecho de utilizar alguno de ellos
depende del tipo de demanda, generalmente utilizado el primero para demandas cons-
tantes y el segundo para demandas variables. Las curvas masas representan en forma
gráfica los volúmenes acumulados en función del tiempo. La pendiente de la curva masa
representa el gasto que pasa por el lugar. Es claro ver que si la pendiente de la curva es
mayor que la pendiente de la curva de escurrimiento, el gasto demandado es mayor que
el aportado por el rı́o.
caso de que se repitan exactamente las aportaciones que se utilizan como datos. Otra
forma de resolver el problema anterior, es mediante un análisis numérico, el cual es más
conveniente cuando la demanda no es contante. Existe el llamado algoritmo del pico se-
cuente para este fin, es decir, cuando las demandas son variables.
A través de los métodos anteriores es posible tener una serie de opciones prelimina-
res de volumen útil. Aunque, ya cuando la presa se encuentra funcionando, los valores
medios se conservan, y la ocurrencia de varios años secos durante la vida útil puede
ocasionar un déficit que haga que la obra deje de ser rentable, o incluso, la ocurrencia de
muchos años húmedos puede ocasionar derrames que pudieran aprovecharse aumen-
tando el volumen útil.
Ya que muchos factores no se toman en cuenta para los cálculos anteriores, tales
como evaporaciones o infiltraciones, siempre es necesario corroborar el valor del volu-
men útil mediante simulaciones del funcionamiento del vaso. Para ello, se cuenta con la
ecuación de continuidad, la cual se expresa de la siguiente manera:
X − D = ∆V, (1.2)
donde Ecp son las entradas por cuenca propia, Et son las entradas por escurrimientos
desde otras cuencas y Ell son las entradas por precipitaciones directas sobre el vaso.
D = Sd + Se + Si + Sderr , (1.4)
Las entradas por cuenca propia son los volúmenes de escurrimiento superficial gene-
rados en la cuenca no controlada que descarga de manera directa a la presa, la cual está
delimitada por el lugar de la cortina y las presas que se encuentran situadas aguas arriba.
Estas entradas se obtienen a partir de los datos recopilados en las estaciones hidrométri-
cas de la zona, y en algunas ocasiones, la misma central hidroeléctrica cuenta con su
propia estación hidrométrica para conocer de manera más precisa los escurrimientos. No
obstante, en caso de no contar con una estación en el lugar exacto, se puede extrapolar
información de estaciones cercanas.
1 INTRODUCCIÓN 32
donde
Los factores que sirven para corregir están en función del área de la cuenca de apor-
tación a la estación i, de la posición de la cuenca y de sus caracterı́sticas. Existen dos
formas de calcular estos factores de corrección.
Vll
F1 = , si n = 1, (1.6)
Vi
donde Vll es el volumen de lluvia que cae en la cuenca propia durante el intervalo de
tiempo analizado y Vi es el volumen de lluvia que cae en la cuenca asociada a la estación
hidrométrica durante el intervalo de tiempo asociado. En algunas ocasiones se puede
utilizar también la expresión
Acp
F1 = , (1.7)
Ai
donde Acp es el área de la cuenca propia y Ai es el área de la cuenca correspondiente a
la estación hidrométrica.
Las entradas por transferencia desde otras cuencas provienen de descargas contro-
ladas o descargas libres de otras presas que se localicen aguas arriba de la presa en
cuestión o en otras cuencas.
Las entradas por lluvia directa sobre el vaso se calculan a partir de aparatos que
registran la cantidad de lluvia directa, y se considera volumen por unidad de área, y
esto se grafica como una altura de precipitación. El área que tenga la superficie libre
del vaso se multiplica por esta altura de precipitación para obtener el volumen de lluvia,
1 INTRODUCCIÓN 33
durante el intervalo de tiempo estudiado. Donde el área está dada de acuerdo con la
curva elevación-área. Es decir
Ell = hpA, (1.10)
donde hp es la altura de precipitación y A es el área promedio del vaso en el intervalo de
tiempo que se estudia.
Luego de considerarse las entradas al vaso, se cuantifican cada una de las salidas. El
volumen extraı́do para satisfacer las demandas, está constituido por la ley de demandas
bajo análisis, y ésta depende del tipo de aprovechamiento que se trate, el cual puede ser
para agua potable, para riego, para generación de energı́a eléctrica, entre otros. También
depende de la relación beneficio/costo de la central. Las leyes de demanda están dadas
por la CFE, a partir de los datos de la población.
Análogo a las evaporaciones, son las infiltraciones, las cuales en general no suelen
ser consideradas por ser cantidades menores de afectación. Sin embargo, también existe
literatura al respecto para su cálculo y análisis.
Finalmente, se considera el volumen derramado por los desbordes que existan en las
obras de excedencia, generalmente a través de los vertedores. Cada central hidroeléctri-
ca también tiene polı́ticas de operación de compuertas, que dependerá de los niveles
caracterı́sticos, de simulaciones previas y del NAMO.
P = γQH, (1.12)
La potencia teórica es la que se presenta antes de la llegada del agua a las turbinas
hidráulicas. La potencia que se transmite a partir de las turbinas se conoce como la
potencia real, y está dada por la siguiente expresión
Pr = ηγQH, (1.14)
donde
P otencia real
η= . (1.15)
P otencia teórica
Generalmente se manejan los resultados en kilowatts, y sabiendo que γ = 1000 kg/m3
para el agua, se tiene que
P = 9.81ηQH. (1.16)
Nuevamente, ya que la potencia es instantánea, lo importante es el tiempo que dicha
potencia se puede sostener, en otras palabras, la energı́a que la central sea capaz de
proporcionar durante un perı́odo de tiempo establecido. Esto es lo que se conoce como
generación. De esta manera, se define lo que es el factor de generación como la relación
entre el volumen V utiliado por la central hidroeléctrica, la carga media Hm y la energı́a
E producida durante el tiempo T . En términos matemáticos
9.81η Hm
fg = , (1.17)
3600
o bien
H
fg = . (1.18)
450
Por otro lado, la energı́a eléctrica producida por una turbina se obtiene mediante el
principio de la utilización del gasto Q y la altura H, los cuales intervienen para producir
energı́a mecánica en el rodete de una turbina. Esta energı́a se transmite al generador,
que a su vez produce electricidad de bajo voltaje, dependiendo de las condiciones de
cada central hidroeléctrica. Este bajo voltaje tiene el propósito de evitar aislamientos
caros entre los alambres de la bobina del generador. Una vez que sale la energı́a del
generador es necesario subir el voltaje para hacer la transmisión a través de los cables.
La potencia se calcula como
P = fp IV, (1.19)
donde fp es el factor de potencia, I la intensidad y V el voltaje. Cuando se sube el voltaje
se puede reducir la intensidad, con lo cual se pueden colocar cables de menor diámetro,
economizando los costos.
1 INTRODUCCIÓN 35
Una vez que la energı́a se encuentra en los cables, este voltaje se eleva mediante
transformadores llamados subestaciones elevadoras. A lo largo de su trayectoria, se
cuenta con varias subestaciones para ir evitando pérdidas de carga. Llegada la carga a
la ciudad, se emplean subestaciones reductoras, las cuales se encargan de distribuir-
las mediante los transformadores en las calles de la ciudad para consumo doméstico.
Además, la potencia que se necesita se grafica mediante una curva llamada curva
de demanda, la cual midel la potencia demandada en función del tiempo. El área bajo
esta curva es la energı́a solicitada en un tiempo determinado. Si en particular, se grafica
el funcionamiento de una central hidroeléctrica, se obtiene una curva de operación. A
fin de garantizar el servicio siempre se necesita una potencia mayor que la solicitada, en
cuyo caso se habla de una reserva. Una central hidroeléctrica cuenta con tres tipos de
reserva: reserva rodante, reserva para mantenimiento y reserva para reparación.
La forma de la curva de operación indica si la central trabaja mucho o poco tiempo con
su potencia máxima. En este sentido, una planta es de pico si trabaja a mayor capacidad
durante las horas de máxima demanda, aun dejando de funcionar fuera de dichas horas.
Si por el contrario, se trabaja con una potencia que no tome muchas variaciones, la
potencia se dice de base. La distinción entre ambos tipos de potencia lo establece el
factor de planta que se define como
P otencia media Gtot
fp = = , (1.20)
P otencia max T0 × Pmax
donde Gtot es la generación total en el perı́odo T0 , en otras palabras, es el área bajo la
curva de operación, y Pm es la potencia media, que se calcula como el cociente de Gtot
y T0 . En términos generales, si fp ≤ 0.40 la planta es de pico, mientras que si fp > 0.40 la
planta se considera de base.
Hasta la década de los años 80, la situación en el paı́s estaba en desarrollo energético
con respecto a las centrales hidroeléctricas. Diversos autores opinan que un paı́s se
considera industrializado si cuenta con más centrales termoeléctricas que hidroeléctricas.
La Figura 1 presenta la capacidad hidroeléctrica del paı́s hasta el año 1987.
1 INTRODUCCIÓN 36
La planeación adecuada del desarrollo eléctrico del paı́s debe conllevar a una correcta
gestión de sus recursos. Se deben estudiar los factores de crecimiento para escoger
las mejores alternativas, sin importar diversas pecualiaridades que llevan a tomar malas
decisiones. Un paı́s debe analizar la capacidad económica, el avance tecnológico, el
potencial energético, el nivel cultural, entre otros, para un avance correcto. Los siguientes
datos que se presentan en la Figura 2 se tomaron de los informes de operación de la
Comisión Federal de Electricidad (CFE).
Figura 1.2: Comparativo entre hidroléctricas y termoeléctricas hasta 1987 [Garcia, 1992]
Una vez que se ha escogido el lugar para construir la cortina, de acuerdo con los estu-
dios previos, las dimensiones de la presa dependerán de los volúmenes que aporta el rı́o
y de las demandas de energı́a que se requiere. La altura y el tipo de cortina dependerá
de la geologı́a y topografı́a de la región. El comportamiento de los cauces se conocerá
mejor entre mayor se conozcan los registros históricos de las estaciones hidrométricas.
Estos escurrimientos se podrán conocer a partir de métodos probabilı́sticos, a fin de po-
der simular el funcionarmiento y obtener polı́ticas de demandas, y criterios de operación.
Con las simulaciones se podrá conocer el volumen de almacenamiento, los niveles de
operación, capacidad de las obras de excedencia, capacidad de la obra de toma, la re-
gulación del vaso, la generación eléctrica y la potencia a instalar.
Sabiendo los datos anteriores se conoce la capacidad del embalse y los niveles de
operación, que son dos de los parámetros más importantes. Para ello, se necesita tener
las aportaciones del poryecto junto a las caracterı́sticas topográficas del vaso y las ex-
tracciones. Cuando los almacenamientos no son grandes, las curvas masas describen
gráficamente las aportaciones contra las demandas. Para centrales hidroeléctricas, las
simulaciones computacionales han permitido realizar una mayor cantidad de cálculos en
1 INTRODUCCIÓN 37
Los datos de entrada al embalse son de dos tipos: el registro histórico de los escurri-
mientos, y los escurrimientos generados por técnicas estocásticas basados en los históri-
cos. Entre más datos tenga un registro histórico, más confiable serán las simulaciones
hechas. El manejo correcto de estos datos, es lo que llevará a una estimación adecuada;
para ello se pueden identificar perı́odos de lluvia o perı́odos de estiaje. Se deben pues
analizar todas las tendencias posibles para diferentes épocas del año, para determinar
los posibles volúmenes que entrarán al vaso. A partir de las técnicas estocásticas se pue-
den simular escurrimientos futuros del perı́odo que se desee. Con ello, se podrá estimar
el comportamiento del vaso; dichos registros simulados es lo que se conoce como regis-
tros sintéticos.
El volumen que debe tener el vaso se determina por los escurrimientos del cauce,
su uso principal, algunos usos secundarios y las limitaciones de la altura de la cortina.
Para su uso principal, se debe plantear la energı́a de base a generar ası́ como la energı́a
de picos. Como se desea mayor potencia, esto significa que se debe tener una mayor
carga para solventar los máximos energéticos. Pero el hecho de mantener un gasto esta-
ble, también es otro factor que hará que la altura se pueda disminuir para mantener una
energı́a hasta cierto punto constante.
min(cX);
Sujeto a :
AX = b
X≥0
1 INTRODUCCIÓN 38
optf (x)
Sujeto a:
hj (x) = 0, j = 1, . . . , m
gj (x) ≥ 0, j = m + 1, . . . , p
x∈ / E n.
donde opt puede representar una maximización o minimización; f (x) es una función con-
tinua, que representa al objetivo; h1 (x), . . . , hm (x) son funciones continuas de restricción;
y X = [x1 x2 · · · xn ] es un vector columna en el espacio euclideano de n dimensiones E n .
Los problemas de programación no lineal se pueden clasificar de la siguiente manera
[Arganis, 2004]:
a) Restringidos;
b) No restringidos;
c) Continuos;
d) Discretos;
e) Diferenciables;
es una cadena de Markov de primer orden y a los valores posibles de S se les llama
estados de proceso. Los procesos estocásticos se pueden entender como modelos fı́si-
cos probabilı́sticos que evolucionan con el tiempo [Salas, 2013]. Además, si el número
de estados de una cadena de Markov es finito, se puede definir la probabilidad de pasar
de un estado cualquiera i en el tiempo t, a un estado cualquiera j en el instante t + 1.
Si las probabilidades se arreglan en una matriz T cuyos elementos aij son las probabili-
dades de pasar del estado i al estado j, en un intervalo de tiempo (es decir de T = t a
T = t + 1) a la matriz T se le llama la matriz de transición. Una matriz de transición se
llama estocástica cuando cumple con las propiedades siguientes:
a) 0 ≤ aij ≤ 1 ∀i, j = 1, . . . , m.
P
b) j aij = 1.
k
Un vector de estado asociado a una etapa cualquiera k, denotado por P , es aquel
cuyos elementos definen las probabilidades de que el proceso se encuentre en cada uno
de los estados en la etapa k. Una cadena de Markov queda totalmente determinada si
se conoce la matriz de transición y el vector de estado para la etapa inicial, con ellos
puede conocerse el vector de estado en etapas posteriores. En términos matemáticos,
los vectores de estado se obtienen multiplicando las matrices de transición como se
observa a continuación:
1 0
P = P T
2 1 0
P = P T = P T2
3 2 0
P = P T = P T3
..
.
n+1 n 0
P = P T = P Tn
b) Los valores de los elementos de cualquier columna j son idénticos entre sı́ e iguales
a la probabilidad a largo plazo asociada al estado j.
1 INTRODUCCIÓN 41
Además, si se conoce P
la matriz de transición en un instante dado, se consideran las
ecuaciones P E = P ET y i Pi = 1, y se resuelve el sistema se puede obtener el vector
de equilibrio.
y si suponemos que los ingresos se distribuyen de acuerdo con una función de probabi-
lidad F (x), se obtiene que:
En el método del valor esperado se pretende minimizar el costo esperado de las fu-
turas operaciones de la distribución del almacenamiento del embalse para la generación
de energı́a. El modelo matemático es muy similar al modelo de un inventario dado por
Bellman [Bellman, 1957], con la diferencia que en un modelo de inventario las salidas de
mercancı́a son las variables aleatorias desconocidas, mientras que en un sistema hidro-
eléctrico, las entradas de agua, son las variables aleatorias que juegan ese papel. Si vi ,
xi y si denotan la cantidad de agua que llega, el volumen del rı́o en el i-ésimo intervalo
de tiempo y el volumen destinado a generar energı́a, entonces
Z ∞ Z ∞
Ei = ··· ki f (xN , . . . , xi+1 , xi |xi−1 , . . . , x0 )dxN · · · dxi (1.30)
−∞ −∞
1 INTRODUCCIÓN 43
es un mı́nimo para cada par (vi , xi−1 ) y por lo tanto Ei (vi , xi−1 ) se resuelve el proble-
ma de optimización. Además, las distribuciones que representan los niveles en el cauce
del rı́o se pueden considerar tomadas de datos estadı́sticos históricos o bien como las
siguientes funciones de probabilidad
mn ∆x ∆
Pi = prob{m x − ≤ xi ≤m x + |xi−1 } =n x. (1.32)
2 2
Los embalses de gran tamaño se emplean para diversos usos, tales como generación
hidroeléctrica, abastecimiento de agua, control de avenidas, recreación, entre otros. Des-
de hace algunas décadas Goicochea clasificó las técnicas multiobjetivos para embalses
en tres categorı́as [Goicochea, Duckeistein, 1982]:
a) Técnicas para generación del conjunto de soluciones no dominadas.
b) Métodos continuos con articulación o preferencias previas.
c) Métodos discretos con una articulación o preferencias progresivas.
Según Ko y Fontane [Ko, Fontane, 1992], el esquema básico del problema de optimi-
zación para embalses se plantea de la siguiente forma. La función objetivo está dada por
máx[F1 , . . . , Fm ] sujeto a las siguientes restricciones:
Vt+1 = Vt + It + SQt − Dt − Et (Vt , Vt+1 ) (1.33)
1 INTRODUCCIÓN 44
V (T ) = c1 X1 (T ) + c2 X2 (T ) + · · · + cn Xn (T ) (1.36)
sujeto a las condiciones iniciales X1 (0), X2 (0), . . . , Xn (0) [SJSU, 2001], donde los coefi-
cientes ci se conocen y el tiempo T también está dado.
El problema es completado cuando se dan las restricciones sobre las variables, intro-
duciendo las variables de control ui , dadas por
dX1
= f1 (X1 , . . . , Xn , u1 , . . . , um )
dt
dX2
= f2 (X1 , . . . , Xn , u1 , . . . , um )
dt
..
.
dXn
= fn (X1 , . . . , Xn , u1 , . . . , um )
dt
La programación dinámica tiene con ventajas que las funciones objetivos y las restric-
ciones no necesariamente deben ser lineales; la solución puede ser una función o toda
una familia de funciones, y se pueden considerar un gran número de etapas sin que la
dimensión de las ecuaciones crezca en gran medida.
Por otro lado, si los estados se representan mediante un vector Fi (xi1 , xi2 , . . . , xim ), el
valor óptimo de la función se puede calcular considerando la primer componente como
variable y el resto como constantes Fi (x0i1 , x0i2 , . . . , x0im ) y optimizando la primer entrada.
Posteriormente se considera la segunda componente como variable y las demás como
constantes, teniendo ya calculada la primer componente, Fi (x1i1 , xi2 , . . . , x0im ) y ası́ suce-
sivamente, hasta la última componente. El procedimiento deja de aplicarse cuando ya no
hay variaciones significativas entre dos aproximaciones consecutivas, es decir, cuando
Fi (x(k) ) ∼
= Fi (x(k−1) ).
1 INTRODUCCIÓN 46
donde
Asimismo, Domı́nguez et. al. (2000) reorganizó la ecuación anterior de la siguiente ma-
nera:
N S1 X
X N S2
∗
BnK1 ,K2 = φn,K1 ,K2 (i1 , i2 ) + qn,K1 (ii , j1 )qn,K2 (i2 , j2 )Bn+1 (j1 , j2 ) (1.47)
j1 =1 j2 =1
donde
N S1
X N S2
X
φn,K1 ,K2 (i1 , i2 ) = qn,K1 (i1 , j1 )bn,K1 (i1 , j1 ) + qn,K2 (i2 , j2 )bn,K1 ,K2 (i1 , j1 , i2 , j2 )
j1 =1 j2 =1
E = f ac ∗ vt ∗ CT ∗ ef ic (1.49)
1.7. Objetivos
1.7.1. Objetivo general
Realizar la modelación de la generación de energı́a de presas hidroeléctricas que for-
man un sistema en cascada empleando técnicas de control óptimo para sistemas mar-
kovianos a horizonte infinito a fin de determinar las polı́ticas óptimas de operación en el
embalse.
2. Modelación estocástica
2.1. Nociones de probabilidad
2.1.1. La σ-álgebra
En esta primera sección se comenzará con algunas definiciones básicas y la termino-
logı́a necesaria de las herramientas probabilı́sticas que se emplearán a fin de abordar el
problema final. Se empieza con la definición de lo que es una σ-álgebra.
i) Ω ∈ F.
ii) Si A ∈ F entonces Ac ∈ F.
∞
iii) Si A1 , A2 , . . . ∈ F, entonces ∪ An ∈ F.
n=1
A la pareja (Ω, F) se le llama espacio medible y a los elementos de F se les llama con-
juntos medibles o eventos, donde los conjuntos A que pertenecen a la σ-álgebra F son
los eventos, a los que se les calcularán las probabilidades.
Algunas de las propiedades más conocidas de las σ-álgebras son las siguientes:
i) ∅ ∈ F.
∞
ii) Si A1 , A2 , . . . ∈ F, entonces ∩ An ∈ F.
n=1
La σ-álgebra de Borel se puede extender en forma natural no sólo a los reales, sino a
R2 , R3 , y en general a Rn , de donde surge el concepto de σ-álgebra producto.
1. P (A) ≥ 0.
2. P (Ω) = 1.
∞ ∞
P
3. P ( ∪ Ak ) = P (Ak ) cuando A1 , A2 , . . . son disjuntos a pares.
k=1 k=1
i) P (Ω) = 1.
Luego, toda función P definida sobre una σ-álgebra que tome valores en el intervalo
[0, 1] se llama una medida de probabilidad.
Algunas de las propiedades principales de una medida de probabilidad son las si-
guientes:
1. P (∅) = 0.
n n
P
2. Si A1 , A2 , . . . , An ∈ F son disjuntos a pares, entonces P ( ∪ Ak ) = P (Ak ).
k=1 k=1
3. P (Ac ) = 1 − P (A).
6. 0 ≤ P (A) ≤ 1.
probabilidad condicional del evento A dado que ocurrió el evento B, la cual se denota
P (A|B) se define como [Loeve, 1977]
P (A ∩ B)
P (A|B) = . (2.3)
P (B)
De igual manera, el Teorema de Bayes, bajo las mismas suposiciones y además que
P (B) > 0 y dado un m ≥ 1 cumple que
P (B|Am )P (Am )
P (Am |B) = P∞ . (2.5)
n=1 P (B|An )P (An )
Dentro de las propiedades de los eventos a los cuales se les calcularán las probabili-
dades, existe una clase de eventos bastante útil, y son los eventos independientes. Dos
eventos A y B son independientes si P (A ∩ B) = P (A)P (B). Como una consecuencia
directa de esta definición se tiene que
Con esto se puede ver toda una relación entre los eventos que siempre ocurren y los
eventos que no ocurren, concluyendo que siempre se tiene independencia en cuales-
quiera de estos casos. Generalmente se descartan estos eventos llamados triviales y
el enfoque se centra en aquellos cuya probabilidad quede estrictamente en el intervalo
abierto (0, 1).
Debe notarse que al pedir que X −1 (B) ∈ F, se está pidiendo que sean eventos o
que sean elementos de la respectiva σ-álgebra, de tal modo que sean factibles para po-
der calcular sus probabilidades. Se recuerda de la sección anterior, que las medidas de
probabilidad P actúan únicamente en los eventos (elementos de F). Es aquı́, donde se
relacionan entonces los elementos no conocidos de Ω con el ya conocido campo de los
números reales R.
2 MODELACIÓN ESTOCÁSTICA 57
Algunas notaciones que se utilizan para los conceptos de variable aleatoria son X −1 (B),
X −1 B, (X ∈ B), {ω ∈ Ω : X(ω) ∈ B}. En particular, si B = [0, ∞), el conjunto anterior se
puede escribir por
{ω ∈ Ω : X(ω) ∈ [0, ∞)} = X −1 [0, ∞) = X ∈ [0, ∞) = (X ≥ 0). (2.10)
Un resultado equivalente a la definición de una variable aleatoria es la siguiente pro-
posición.
Proposición. Una función X : Ω → R es una variable aleatoria si y sólo si, para cada
x ∈ R se cumple que (X ≤ x) ∈ F.
Este resultado permite entonces determinar qué función es variable aleatoria, única-
mente describiendo los conjuntos (X ≤ x), que son las imágenes inversas, y esto para
cada x ∈ R, y entonces, ya no se necesita comprobar la propiedad que describı́a a las
variables aleatorias por medio de los conjuntos borelianos, los cuales pueden ser más
abstractos. En pocas palabras, sólamente basta comprobar que los intervalos (−∞, 0]
pertenezcan a la σ-álgebra dada.
Proposición. Las siguientes funciones cumplen con las caracterı́sticas de variables alea-
torias.
Lo anterior permite tratar con más de una variable aleatoria y poder combinarlas a fin
de formar nuevas variables aleatorias. Muchos de los fenómenos que se estudian invo-
lucran a más de un parámetro para calcular, y cada uno de ellos con una naturaleza no
determinista. Es importante saber pues, que muchas de las transformaciones algebraicas
conocidas no afectan la estructura de las variables aleatorias. Además, dado que muchos
experimentos y fenómenos son de naturaleza infinita, es importante conocer conjuntos
de variables aleatorias infinitas. La siguiente proposición resume este hecho.
b) Sea X1 , X2 , . . . una sucesión infinita de variables aleatorias tales que para cada ω ∈
Ω los números sup{X1 , X2 , . . .} e ı́nf{X1 , X2 , . . .} son finitos. Entonces las funciones
sup{Xn : n ≥ 0} e ı́nf{Xn : n ≥ 0} son variables aleatorias.
c) Sea X1 , X2 , . . . una sucesión infinita de variables aleatorias tales que para cada
ω ∈ Ω los números lı́m sup{X1 , X2 , . . .} y lı́m inf{X1 , X2 , . . .} son finitos. Entonces
las funciones lı́m supn→∞ Xn y lı́m inf n→∞ Xn son variables aleatorias.
En el inciso c) se establecen los conceptos de lı́mites, recordando que no son más que
intersecciones y uniones de conjuntos, de manera más formal
∞ ∞
lı́m sup An = ∩ ∪ Ak (2.11)
n→∞ n=1 k=n
∞ ∞
lı́m sup An = ∪ ∩ Ak (2.12)
n→∞ n=1 k=n
Además, si el lı́mite superior y el lı́mite inferior coinciden, se dice que existe el lı́mite, con
lo que ya tiene sentido el inciso d).
b) lı́mx→−∞ F (x) = 0.
Las funciones de distribución pueden ser continuas o discretas dependiendo los valo-
res que tomen en el eje real X. La figura 2.7 muestra cada uno de los dos casos, y aquı́
se observan en forma más clara las propiedades anteriores, como el hecho de que son
crecientes y toman el valor de 0 en menos infinito y el valor de 1 en más infinito.
Proposición. Sea X una variable aleatoria con función de distribución F (x). Para cua-
lesquiera números reales a < b se cumple que:
a) P (X < a) = F (a−).
b) P (X = a) = F (a) − F (a−).
2 MODELACIÓN ESTOCÁSTICA 60
d) P (a ≤ X ≤ b) = F (b) − F (a−).
donde la notación F (a−) indica que se toma el lı́mite lateral por la izquierda.
Una función de probabilidad o también llamada función de masa cumple con las si-
guientes propiedades:
P
b) x f (x) = 1.
Esto nos indica que el área bajo la curva de la función de densidad es precisamente la
probabilidad de todo el conjunto Ω, es decir, P (Ω) = 1. Geométricamente, el cálculo de
probabilidades de intervalos se reduce ahora al cálculo de áreas bajo curvas para las
funciones de densidad continuas, mientras que para las funciones discretas, las probabi-
lidades son las imágenes precisamente de la función.
En primer lugar se define lo que es una partición del intervalo [a, b] como un conjunto
finito de puntos P = {a = x0 , x1 , . . . , xn = b} donde xi−1 < xi para todo i = 1, 2, . . . , n.
Dentro de esta partición se toman los puntos ci ∈ [xi−1 , xi ]. Nótese que una partición
induce subintervalos dentro del intervalo original como se observa en la Figura 2.9.
y
g(xi ) = ı́nf{g(x) : xi−1 < x ≤ xi }.
Se define la suma superior e inferior de Riemann-Stieltjes de la siguiente manera:
n
X
Sn = h(xi )[F (xi ) − F (xi−1 )], (2.14)
i=1
y
n
X
Sn = h(xi )[F (xi ) − F (xi−1 )]. (2.15)
i=1
Y entonces se define Z b Z b
g(x)dF (x) = lı́m hN (x)dF (x) (2.18)
a N →∞ a
cuando este lı́mite existe. Además, se puede extender la definición de esta integral de la
siguiente forma
Z ∞ Z b
g(x)dF (x) = lı́m g(x)dF (x), (2.19)
−∞ a,b→∞ a
Teorema. Sea X con función de distribución FX (x) y sea g : R → R una función Borel
medible tal que g(X) tiene esperanza finita. Entonces
Z ∞
E[g(X)] = g(x)dFX (x). (2.22)
−∞
Se observa pues que la esperanza matemática comparte propiedades con las inte-
grales conocidas de Riemann, más aún, con la integral de Riemann-Stieltjes. Es por ello
que la esperanza tiene las siguientes propiedades:
b) E(cX) = cE(X).
d) Si X ≥ 0, entonces E(X) ≥ 0.
Por otra parte, existe otra medida que indica cuánta dispersión existe en los datos
de una variable aleatoria. Esto se conoce como la varianza. La varianza de una variable
aleatoria X, denotada por V ar(X), se define como V ar(X) = E(X − E(X))2 . Depen-
diendo el caso en que X es discreta o continua, la varianza se puede expresar como:
b) V ar(c) = 0.
c) V ar(cX) = c2 V ar(X).
d) V ar(X + c) = V ar(X).
a) Si X ≥ 0 entonces E(X|G) ≥ 0.
d) E(E(X|G)) = E(X).
conjunto A.
Procesos de Markov. Un proceso de Markov es un modelo donde una vez que se conoce
el estado actual del sistema, los estados pasados no influyen en los estados futuros
del sistema. Matemáticamente, una propiedad se dice markoviana si dados los estados
x0 , x1 , . . . , xn−1 , xn , xn+1 , se tiene que
Es interesante notar que los sistemas dinámicos deterministas se pueden expresar me-
diante modelos markovianos, pues la evolución del sistema se determina por una posi-
ción inicial y su ley de transición.
E(Xn+1 : X0 = x0 , . . . , Xn = xn ) = xn . (2.27)
Por otro lado supongamos que i y j son dos estados. A la probabilidad condicional
P (Xn+1 = j|Xn = i) (2.29)
se le llama probabilidad de transición del estado i en el tiempo n al estado j en el tiempo
n + 1, y se denota simplemente como pij (n, n + 1). Si ocurriese que pij (n, n + 1) no de-
pendiese de n, se dice que la cadena es estacionaria en el tiempo.
a) pij ≥ 0.
P
b) j pij = 1.
Estas dos propiedades forman lo que se llama una matriz estocástica. El segundo
inciso dice que a partir de cualquier estado i con probabilidad uno, la cadena pasa nece-
sariamente a algún elemento del espacio de estados al siguiente momento.
Finalmente, una distribución inicial para una cadena de Markov con espacio de estado
discreto, es una distribución de probabilidad sobre el conjunto {0, 1, 2, . . .}, es decir, es
una colección de números p0 , p1 , p2 , . . . no negativos y que suman uno, donde pi repre-
senta la probabilidad de que la cadena inicie en el estado i [Basu, 2003].
Para un estado inicial x0 , una polı́tica óptima π ∗ es aquella que minimiza el costo, esto
es
Jπ∗ (x0 ) = mı́n Jπ (x0 ), (2.39)
π∈Π
donde Π es el conjunto de polı́ticas admisibles. Pese a que se parte del estado inicial
x0 , la programación dinámica nos dice que es posible encontrar una polı́tica π ∗ que sea
óptima para cualquier estado inicial del sistema. El costo óptimo depende del estado
inicial x0 y se denota por J ∗ (x0 ), luego
En ocasiones J ∗ se ve como una función que asigna a cada estado inicial x0 un costo
óptimo J ∗ (x0 ) y se llama la función de costo óptimo.
Principio de optimalidad. Sea π ∗ = {µ∗0 , µ∗1 , . . . , µ∗N −1 } una polı́tica óptima para el proble-
ma básico y supóngase que el estado xi ocurre en el tiempo i con probabilidad positiva,
cuando se emplea la polı́tica π ∗ . Considérese además que se está en el estado xi al
tiempo i y se desea minimizar el costo del tiempo i al tiempo N dado por
( N −1
)
X
E gN (xN ) + gk (xk , µk (xk ), wk ) . (2.41)
k=i
Luego, la polı́tica truncada {µ∗i , µ∗i+1 , . . . , µ∗N −1 } es óptima para este último propósito.
En otras palabras, una polı́tica óptima que minimiza un costo en general, también es
óptima en cada una de sus partes que forman el total. Esta es la parte medular de la
programación dinámica. Este principio permite seccionar un problema en cierto número
de pasos, garantizando que la optimalidad no se pierda en cada uno de ellos. Matemáti-
camente, esto se enuncia a través de la siguiente proposición.
Proposición. Para cada estado inicial x0 , el costo óptimo J ∗ (x0 ) del problema básico es
igual a J0 (x0 ), donde la función J0 está dada por el último paso del siguiente algoritmo en
retroceso desde el tiempo N − 1 al tiempo 0:
Jk (xk ) = mı́n Ewk {gk (xk , uk , wk ) + Jk+1 (fk (xk , uk , wk ))} , k = 0, 1, , N − 1, (2.43)
uk ∈Uk (xk )
2 MODELACIÓN ESTOCÁSTICA 73
Cuando el espacio de estados se vuelve más abstracto, surgen los llamados pro-
blemas de control markoviano a tiempo discreto sobre espacios Borel medibles, éstos
serán descritos más adelante. Para optimizar estos procesos se puede emplear la llama-
da ganancia esperada promedio y la ganancia descontada promedio, donde es posible
considerar el problema con o sin restricciones sobre un número finito de costos espera-
dos promedio y costos descontados promedios, respectivamente.
Existen diversas técnicas para analizar los problemas con restricciones. Dentro de
estas se encuentra el método directo, la programación lineal y el método de desvaneci-
miento.
El inciso i) dice que las distancia siempre son positivas o bien cero; el inciso ii) dice
que el único caso en que la distancia es cero es cuando se mide de un punto a él mismo;
el inciso iii) consiste en que la distancia se puede medir tomando como punto inicial cual-
quier punto; y finalmente, el inciso iv) dice que la distancia entre dos puntos siempre es
menor si se mide directamente, a que se realice una pausa en un punto intermedio. Estas
propiedades surgen de manera natural y existen muchos ejemplos de métricas como la
distancia que medimos normalmente entre dos objetos.
Por otro lado, se recuerda que una función es biyectiva si es tanto inyectiva como
suprayectiva. Si existe una biyección entre un conjunto arbitrario y el conjunto de los na-
turales N, se dice que el conjunto tiene la cardinalidad de los naturales, o bien que el
conjunto es numerable. En este sentido, partiendo de lo precedente, un espacio métrico
se dice separable si posee un subconjunto denso y numerable.
Los espacios separables junto con los espacios completos forman lo que se conoce
como los espacios polacos. Estos obtienen su nombre de sus autores, que principal-
mente fueron la mayorı́a polacos, entre ellos Tarski, Sierpinski y Kuratowski. Los espa-
cios polacos son importantes ya que son la base de toda la teorı́a de conjuntos y por sus
amplias aplicaciones a la teorı́a de la probabilidad, que es en donde estará basado prin-
cipalmente todo este trabajo. Ası́ pues, en concreto, un espacio métrico se dice polaco si
es completo y separable [Marsden, 1974].
Ya que se tienen los estados posibles X, las acciones admisibles A(x) y las histo-
rias Ht , considérese a todas las posibles decisiones que se puedan tomar dentro de las
acciones admisibles, y denótese a este conjunto como F. Explicando un poco más, es-
te conjunto F, consta de todas las funciones f : X → A, tal que f (x) ∈ A(x); dicho
de otra manera, son funciones llamadas deterministas tomadas por un controlador, son
funciones seleccionadas intencionalmente dentro del posible conjunto de decisiones que
existen, son las acciones reales que modificarán la evolución del sistema, dado que se
conoce el estado del sistema [Borkar, 1994].
Con todo lo anterior ya se puede definir qué son las polı́ticas de control, y lo más
2 MODELACIÓN ESTOCÁSTICA 77
Ahora bien, dado que las probabilidades están sujetas a la ocurrencia de las historias
admisibles, se necesita el concepto de esperanza condicional, el cual depende a su vez
de las polı́ticas de control y las distribuciones iniciales. En este sentido el operador es-
peranza con respecto a Pνπ se denotará por Eνπ . Además si ν es la distribución inicial que
se concentra en el estado inicial x ∈ X, se denotarán Pνπ y Eνπ como Pxπ y Exπ , respectiva-
mente. Más aún, si π = ϕ es una polı́tica estacionaria, luego se denotarán Pνπ y Eνπ como
Pxϕ y Exϕ , respectivamente.
Por otro lado, sea ϕ ∈ ΦS un kernel estocástico sobre A dado X, c y r las funciones
de costo y ganancia, y Q la ley de transición. Luego se define, para cada x ∈ X,
Z Z
cϕ (x) := c(x, a)ϕ(da|x), rϕ (x) := r(x, a)ϕ(da|x), (2.44)
A A
y Z
Qϕ (·|x) := Q(·|x, a)ϕ(da|x). (2.45)
A
Estas ecuaciones anteriores se pueden interpretar como el cálculo de medidas con res-
pecto al conjunto de acciones A dado que se tiene el estado X. Se recuerda que las
funciones de costo y ganancia dependen tanto del estado x como de la acción a reali-
zada. La primera y segunda ecuación representan el costo y la ganancia en función del
estado inicial x cuando se tiene la polı́tica ϕ; la tercera ecuación, representa la ley de
transición dado que se tiene el estado x y se presenta la polı́tica ϕ.
convierten en
cf (x) = c(x, f (x)), rf (x) = r(x, f (x)) y Qf (B|x) = Q(B|x, f (x)). (2.46)
Si π = ϕ es una polı́tica estacionaria, la probabilidad de transición en el n-ésimo paso
será denotada por Qnϕ , esto significa
Qnϕ (B|x) := Pxϕ (xn ∈ B), n ∈ N0 , B ∈ B(X), x ∈ X, (2.47)
con Q1ϕ (·|x) := Qϕ (·|x) y Q0ϕ (·|x) = δx la medida de Dirac concentrada en el estado inicial
x. Además, se puede escribir recursivamente ([Borkar, 1994])
Z
Qnϕ (B|x) = Qϕ (B|y)Qn−1
ϕ (dy|x)
ZX
= Qn−1
ϕ (B|y)Qϕ (dy|x) n ≥ 1. (2.48)
X
Con esto, el cálculo de probabilidades de transición Pxϕ quedará ligado a los kérneles
estocásticos, y especı́ficamente a la ley de transición Q, y además se tendrá una relación
para conocer la transición en cualquier paso del proceso, a través de una fórmula recur-
siva, es decir, si se conoce la ley de transición en el n-ésimo paso, se podrá conocer la
ley en el paso n + 1.
Con toda esta terminologı́a, sea ν una distribución inicial de probabilidades sobre los
estados X. Si π, que se puede denotar como {ϕt }, es una polı́tica de Markov aleatori-
zada, es decir, {ϕt } ∈ ΠRM , entonces el proceso estocástico de estados {xt } también
es un proceso de Markov. Y no sólo ello, sino que depende del tiempo, es decir, es no
homogéneo y sus kerneles de transición son {Qϕt (·|·)}∞t=0 . Con esto se está diciendo que
el cálculo de probabilidades no depende de los estados pasados, sino sólo del estado
actual, lo que define que un proceso sea Markoviano; en términos matemáticos
Pνϕt (xt+1 ∈ B|x0 , x1 , . . . , xt ) = Pνϕt (xt+1 ∈ B|xt ) = Qϕt (B|xt ). (2.49)
En particular, si π = {ϕ} es una polı́tica de Markov estacionaria aleatorizada, también
vale la relación anterior, es decir
Pνϕ (xt+1 ∈ B|x0 , x1 , . . . , xt ) = Pνϕ (xt+1 ∈ B|xt ) = Qϕ (B|xt ). (2.50)
Para ello se supondrán las siguientes hipótesis que se retomarán en el cálculo de las
polı́ticas óptimas de control. Para cada estado x ∈ X
2 MODELACIÓN ESTOCÁSTICA 80
b) Tanto las funciones de costo c(x, a) como la función de ganancia r(x, a) son funcio-
nes continuas cuando dependen de a ∈ A(x).
d) Tantos los costos como las ganancias estarán acotadas por una función llamada de
peso W (x) y una constante no negativa M , esto es
|c(x, a)| ≤ M W (x), |r(x, a)| ≤ M W (x), para cada (x, a) ∈ K. (2.51)
El inciso a) permite que las acciones se tomen dentro de conjuntos acotados, inclu-
yendo los extremo posibles de un intervalo dado. El inciso b) y el inciso c) permitirá que
al tratar con conjuntos compactos, las funciones continuas siempre toman sus valores
máximos y mı́nimos, con lo cual se podrán tener los valores óptimos. Finalmente, el inci-
so d) permite que las soluciones no se excedan de las restricciones, al poder mantener
a los costos y a las ganancias por debajo de funciones acotadas.
Por otro lado sea N un número natural. El costo por etapas y la ganancia por etapas,
hasta la etapa N se definen como [Chang, 2006]
N
X N
X
r(xt , at ), y c(xt , at ). (2.52)
t=0 t=0
Como los estados y las acciones son variables aleatorias, los costos y ganancias no
se pueden expresar de manera determinista, pero sı́ se pueden calcular en promedio.
2 MODELACIÓN ESTOCÁSTICA 81
Ası́ pues, la ganancia promedio descontada y el costo promedio descontado dada una
polı́tica π ∈ Π y el estado inicial x ∈ X se definen respectivamente por [Costa, 2012]
"∞ # "∞ #
X X
Vα (x, π, r) := Exπ αt r(xt , at ) , y Vα (x, π, c) := Exπ αt c(xt , at ) . (2.54)
t=0 t=0
Ahora bien, dentro de las condiciones de un problema fı́sico, existen posibles res-
tricciones en cuanto a costos y a ganancias. Este tipo de restricciones pueden ser de
tipo fijo o bien, restricciones variables. Supóngase pues, el caso general, donde la res-
tricción sea una función. Dicha función se pedirá que sea una función medible, y res-
tringirá los posibles estados. Esta función se denotará por θα : X → R, y se supondrá
además, que se encuentra acotada por otra función W , llamada función de peso, es de-
cir, |θα (x)| ≤ kθα kW (x), para cada x ∈ X [Dutta, 1991].
Ya que las restricciones actuarán sobre variables aleatorias, nuevamente los cálculos
se realizarán en promedio. Ası́ pues, se define la restricción total esperada α-descontada
cuando el controlador usa la polı́tica π ∈ Π ([Chang, 2006]), dado el estado inicial x ∈ X
como "∞ #
X
θα (x, π) = (1 − α)Exπ αt θα (xt ) . (2.57)
t=0
Esta restricción total también satisface la ecuación de Poisson, es decir, si 0 < α < 1
y θα : X → R es una función de restricción, entonces para cada polı́tica estacionaria
2 MODELACIÓN ESTOCÁSTICA 82
Esta ecuación permitirá pues, que las restricciones también estén en función de los
kerneles de transición Q, que se tienen desde un principio.
es decir, que no exceda al valor más pequeño de los costos α-descontados ni al valor
más grande, de entre todas las posibles polı́ticas. Con esto se garantiza, que la función
de restricción, tome valores factibles y reales dentro del problema de optimización.
Para este valor constante, la función de restricción total esperada α-descontada será
θα (y, π) = θα para cada y ∈ X, y el conjunto Fθxα será simplemente Fθxα = {π ∈ Π :
Vα (x, π, c) ≤ θα }.
Asimismo, las funciones θα,mı́n (·) y θα,máx (·) satisfacen las ecuaciones de optimalidad
de Bellman, para cada y ∈ X
Z
θα,mı́n (y) = mı́n c(y, a) + α θα,mı́n (z)Q(dy|y, a) , (2.62)
a∈A(y) X
2 MODELACIÓN ESTOCÁSTICA 83
y Z
θα,máx (y) = máx c(y, a) + α θα,máx (z)Q(dy|y, a) . (2.63)
a∈A(y) X
En resumen, con toda esta terminologı́a, se dice que una polı́tica π ∗ ∈ Π es óptima
para el problema descontado con restricciones con estado inicial x ∈ X, si π ∗ ∈ Fθxα y
además
Vα (x, π ∗ , r) = sup Vα (x, π, r) para cada x ∈ X. (2.64)
π∈Fθxα
En este caso, Vα∗ (x, r) = Vα (x, π ∗ , r) se llama la ganancia α-descontada óptima para el
problema descontado con restricciones con estado inicial x.
De esta manera, las polı́ticas óptimas para el problema descontado con restricciones,
serán aquellas que maximicen el valor de la función α-descontada de ganancia r, es
decir, de entre todas las polı́ticas que no excedan a la restricción θα , aquella que cumpla
con que Vα (x, π, r) sea el mayor valor, es decir, la mayor ganancia, será la polı́tica óptima,
y se denotará por π ∗ . Y al valor Vα (x, π ∗ , r) el cual se recuerda que es
"∞ #
∗
X
Vα (x, π ∗ , r) = Exπ αt r(xt , at ) (2.65)
t=0
que a partir de ahora se denotará como Vα∗ (x, r) será la ganancia óptima.
Sea λ ≤ 0, es decir, un parámetro que puede tomar el valor de cero o ser negativo.
Considérese además las respectivas funciones de costo y ganancias, dado un estado x
y una acción a, c(x, a) y r(x, a), respectivamente. Además, sea θα una función de res-
tricción y defı́nase la siguiente función, que se llamará función de ganancia generalizada
[Feinberg, Kasyanov, Zadoianchuk, 2012]
Concepto Notación
PN
Costo por etapas c(xt , at )
Pt=0
N
Ganancia por etapas r(xt , at )
Pt=0
N
Costo descontado por etapas αt c(xt , at )
Pt=0
N t
Ganancia descontada por etapas t=0 α r(xt , at ) P
Costo promedio descontado Vα (x, π, c) = Exπ [ P∞ t
t=0 α c(xt , at )]
Ganancia promedio descontada Vα (x, π, r) = Exπ [ ∞ t
t=0 α r(xt , at )]
Función de restricción θα (x)
Restricción total esperada α-descontada θα (x, π) = (1 − α)Exπ [ ∞ t
P
t=0 α θα (xt )]
Restricción mı́nima θα,mı́n (x) = ı́nf π∈Π Vα (x, π, c)
Restricción máxima θα,máx (x) = supπ∈Π Vα (x, π, c)
Polı́tica óptima π∗
Ganancia óptima Vα (x, π ∗ , r) = supπ∈Fθx Vα (x, π, r)
α
Ganancia óptima Vα∗ (x, r) = Vα (x, π ∗ , r)
Esto significa que debe existir una polı́tica óptima para el problema α-descontado, y
que la ecuación de Bellman es factible a ser resuelta, dentro del conjunto de acciones
admisibles en A(x). De esta forma, se relaciona la ganancia generalizada rαλ con el Ker-
nel estocástico Q [Guo, Quanxin, 2006].
Ası́ pues, regresamos a las funciones de descuento de un principio Vα (·, π ∗ , rαλ ) para
poder calcular las polı́ticas óptimas. Esto sólo pone de manifiesto, la importancia de con-
siderar tales funciones a través de los respectivos factores α descontados.
Este resultado enlaza las condiciones necesarias para que se pueda resolver el pro-
blema descontado sin restricciones tanto para las polı́ticas aleatorizadas como para las
polı́ticas deterministas.
Este resultado permite calcular los máximos de las funciones α-descontadas de ga-
nancias a través de las derivadas al igual que en el cálculo. Si una función es diferen-
ciable, geométricamente el punto donde la pendiente es cero, es un punto candidato a
ser un punto óptimo. La función α-descontada de ganancia depende de el parámetro λ
dado por el multiplicador de Lagrange, que a su vez depende de la función generalizada
rαλ . Al derivarlo parcialmente con respecto a este λ, es lo mismo que tomar la función α-
descontada de costo con respecto a la polı́tica παλ y calcular su diferencia con la función
de restricción total θα , dentro del conjunto de polı́ticas óptimas. La proposición anterior
2 MODELACIÓN ESTOCÁSTICA 87
Teorema. [Guo, Quanxin, 2006] Bajo todas las hipótesis anteriores, sea α ∈ (0, 1) un
factor de descuento y x ∈ X un estado inicial. Entonces
Vα (x, π̂, c) = θα (x, π̂), y Vα (x, π̂, rαΛ ) = Vα∗ (x, rαΛ ). (2.76)
Por lo tanto π̂ es una polı́tica α-óptima para el problema descontado con restriccio-
nes y Vα∗ (x, rαΛ ) es el valor α-óptimo para el problema descontado con restricciones.
Más aún
Vα∗ (x, rαΛ ) = ı́nf Vα∗ (x, rαΛ ). (2.77)
Λ≤0
b) Si λ∗ < 0 es un punto crı́tico de λ → Vα∗ (x, rαλ ), esto es, si la derivada de la pro-
∗ ∗
posición anterior es igual a cero en λ = λ∗ , luego cada polı́tica παλ ∈ Πλα es
α-óptima para el problema descontado con restricciones. Además se tiene que
∗ ∗ ∗
Vα (x, παλ , c) = θα (x, παλ ), siendo Vα∗ (x, rαλ ) el valor α-óptimo para el problema des-
∗
contado con restricciones, el cual coincide con Vα (x, παλ , r). También se cumple que
∗ ∗
Vα∗ (x, rαλ ) = ı́nf Vα∗ (x, rαλ ). (2.78)
λ≤0
c) Si πα0 ∈ Π0α satisface Vα (x, πα0 , c) ≤ θα (x, πα0 ); esto es, πα0 ∈ Fθxα , luego πα0 es una
polı́tica α-óptima para el problema descontado con restricciones. En este caso
Vα∗ (x, rα0 ) = Vα∗ (x, r) es el valor α-óptimo para el problema descontado con res-
tricciones y coincide con Vα (x, πα0 , r). Además se tiene
El inciso a) del teorema anterior establece que si existe una polı́tica tal que el costo
que genera es igual a su restricción total y su ganancia generalizada, para cierto paráme-
tro Λ es la óptima, entonces, esta misma polı́tica que se ha encontrado también es ópti-
ma para el problema descontado con restricciones, y las ganancias óptimas coinciden
y de hecho, será el mı́nimo de todas las ganancias generalizadas posibles, variando al
parámetro λ.
Estas funciones están asociadas a los respectivos costos y ganancias que se pre-
sentan en los distintos estados, en forma acumulativa, dado que se realizan las acciones
respectivas. Toda la recopilación de respectiva información genera siempre un costo, pe-
ro a su vez produce una ganancia. Para el caso estocástico, los costos y las ganancias
serán en promedio, ya que no se dispone de información, más que lo concerniente a las
variables aleatorias correspondientes.
de manera respectiva.
2 MODELACIÓN ESTOCÁSTICA 89
Del análisis matemático se conoce que una función puede no tener lı́mite, pero siem-
pre existirá el lı́mite inferior respectivo, es por ello que las ecuaciones anteriores se pre-
sentan bajo estos lı́mites en el infinito. El factor n1 permitirá el cálculo promedio de las
respectivas funciones. Nótese además que para la ganancia se considera un lı́mite in-
ferior, que serı́a uno de los escenarios más desfavorables, el caso donde se obtenga la
mı́nima ganancia. Mientras que para el costo, se considera el lı́mite superior, en otras
palabras, el costo más elevado que se pueda presentar. Ası́, se garantiza que los resul-
tados siempre se encuentren acotados por las cuestiones menos apropiadas.
La siguiente condición será necesaria para evitar que las soluciones que se presen-
ten a las ecuaciones no tengan soluciones, y más aún, en términos matemáticos, sean
soluciones estables.
En este sentido, las funciones J(x, ϕ, r) y J(x, ϕ, c) son constantes, pues no dependen
del estado inicial x, con lo cual verifican que
n−1
1 ϕX
J(x, ϕ, r) = lı́m Ex r(xk , ak ) = µϕ (rϕ ). (2.85)
n→∞ n
k=0
Además, se denotará como g(ϕ, r) a la medida µϕ (rϕ ). Es decir, que la ganancia pro-
medio esperada es igual a la medida única establecida por la condición de ergodicidad
geométrica. Análogamente, para el costo se tiene que
n−1
1 ϕX
J(x, ϕ, c) = lı́m Ex c(xk , ak ) = µϕ (cϕ )., (2.86)
n→∞ n
k=0
y Z
θmáx = máx cϕ (y)µϕ (dy), (2.88)
ϕ∈Φ X
donde se puede notar que dadas las condiciones particulares, dichos parámetros son
finitos. Para evitar situaciones triviales, se supondrá una restricción a la constante θ de
tal forma que
θmı́n < θ < θmáx . (2.89)
Más adelante estas restricciones estaŕan representadas por coeficientes de pena-
lización en caso de cumplirse, a partir de los diversos criterios multiobjetivos que se
mencionaron en la sección de métodos de optimización.
sujeto a
π ∈ Π y J(x, π, c) ≤ θ para cada x ∈ X. (2.91)
Véase que la función a maximizar, es la ganancia promedio, donde las restricciones
son como siempre la polı́tica de control, el estado inicial x y la constante θ para el costo c.
Con este nuevo problema planteado, se deben definir las nuevas condiciones que
llevarán a resolverlo. En primer lugar, dada una polı́tica π ∈ Π, esta se dice ser admisi-
ble para el problema esperado con restricciones si esta satisface que J(x, π, c) ≤ θ. Al
conjunto de polı́ticas admisibles para el problema esperado con restricciones como
Más aún, una polı́tica admisible π ∗ se llama polı́tica óptima para el problema esperado
con restricciones si J(x, π, r) ≤ J(x, π ∗ , r) para todo estado x ∈ X y para cada polı́tica
admisible π ∈ Π. Con esto, se puede definir la función de valor óptimo para el problema
esperado con restricciones como
Por otro lado, sea θ : X → R una función acotada llamada la tasa de restricción. La
restricción esperada promedio cuando el controlador tiene una polı́tica π ∈ Π, dado el
estado inicial x ∈ X se define como
n−1
1 X
θ(x, π) = lı́m inf Exπ θ(xk ). (2.94)
n→∞ n
k=0
Nuevamente, como cada uno de los estados se representa mediante variables aleato-
rias, se calcula esta función a través de sus esperanzas matemáticas. Además, el lı́mite
se toma inferiormente, al no poderse garantizar que existe el lı́mite de la función. La fun-
ción de restricción esperada promedio se define con propiedades análogas al caso de
las medidas invariantes que se presentan en la propiedad ergódica geométrica, es decir,
aproximarán a las integrales estocásticas con las leyes de transición. En este sentido, si
ϕ ∈ Φ es una polı́tica estacionaria aleatorizada, luego θ(x, ϕ) = µϕ (θ), donde µϕ es la
medida de probabilidad invariante asociada a ϕ descrita anteriormente.
Bajo esta nueva notación, una polı́tica π ∈ Π se dice ser admisible para el problema
esperado con restricciones generalizado si esta satisface que J(x, π, c) ≤ θ(x, π) para
todo estado x ∈ X. Además, denotamos el conjunto de polı́ticas admisibles para el pro-
blema esperado con restricciones generalizado como
Siguiendo esta parte, una polı́tica admisible π ∗ se llama óptima para el problema espe-
rado con restricciones generalizado si J(x, π, c) ≤ J(x, π ∗ , r) para toda x ∈ X y para cada
polı́tica admisible π ∈ Π. De aquı́, se define la función de valor óptimo para el problema
esperado con restricciones generalizado como [Guo, Quanxin, 2006]
Para fines de estudio, una región hidrológica es la agrupación de varias cuencas hi-
drológicas con niveles similares de escurrimiento superficial. En las figuras 3.11 y 3.12
se muestran cada una de las regiones hidrológicas administrativas con sus subdivisiones
respectivas.
En total, la región hidrológica número 30 abarca una extensión de 102,641 km2 . Esta
región es la más húmeda del paı́s, y en ella se localizan los rı́os más caudalosos de
México; ejemplos de ellos son los rı́os Usumacinta y Grijalva, los cuales desembocan en
el Golfo de México.
mente, haciendo un total de 75.22 %, del territorio total estatal. El estado de campeche,
cuenta con los rı́os Usumacinta y Laguna de Términos, los cuales representan el 2.58 %
y 30.46 %, respectivamente, haciendo un total del 33.04 %. Finalmente, los estados de
Oaxaca y Veracruz cuentan con el Rı́o Grijalva-Tuxtla Gutiérrez, abarcando el 1.02 % y
3 DATOS DEL LUGAR DE ESTUDIO 96
Posteriormente, aguas abajo de la presa Peñitas recibe las aportaciones de los rı́os
Platanar y Camoapa, que dan origen al rı́o Mezcalapa; después se bifurca en los rı́os Sa-
maria por su margen izquierda y el rı́o Carrizal por su margen derecha; éste último cruza
la Ciudad de Villahermosa, capital del estado de Tabasco, donde recibe las aportaciones
del rı́o Pichucalco y La Sierra, que nacen en las montañas del bajo Grijalva. Después de
Villahermosa, continúa el rı́o Grijalva hasta confluir con el rı́o Usumacinta para finalmente
desembocar al Golfo de México [CFE, 1988].
La cuenca del rı́o Grijalva-Usumacinta, es una cuenca que divide fronteras; nace en
la república de Guatemala, cruza los estados de Chiapas y Tabasco, y una menor parte
recorre el estado de Campeche, para finalmente atravesar los estados de Oaxaca y Ve-
racruz.
La cuenca respectiva tiene una superficie aproximada de 131,157 km2 , de los cuales
aproximadamente 52,600 km2 corresponden a la cuenca del rı́o Grijalva y 78,757 km2
a la cuenca del rı́o Usumacinta. Además, comprende cuatro porciones geográficas bien
definidas que se conocen con los nombres de Alto Grijalva, Medio Grijalva, Bajo Grijalva
(Sierra) y Bajo Grijalva (Planicie) [CFE, 1976].
En el Bajo Grijalva (Sierra), se ubica la Sierra del Norte de Chiapas. Esta se compone
de una serie de serranı́as separadas por largos valles que bordean a los Altos y las Mon-
tañas del Oriente. De acuerdo a la Secretarı́a de Energı́a (SE), gracias a la localización
3 DATOS DEL LUGAR DE ESTUDIO 98
de estas montañas es que la humedad se puede interceptar, de tal modo que los vientos
del Golfo, junto con la respectiva orografı́a ocasiona que el clima sea húmedo todo el año
con lluvias en todo su perı́odo [SE, 2009].
El Bajo Grijalva (Planicie) cuenta con diversas planicies ubicadas en la Llanura Coste-
ra del Golfo, la cual está ocupada en su mayorı́a por el estado de Tabasco y está formada
por grandes cantidades de aluvión acarreado por los rı́os más caudalosos del paı́s. Di-
chos rı́os son el Usumacinta, el Grijalva, el Papaloapan y el Coatzacoalcos. Estos rı́os
atraviesan esta planicie, para finalmente desembocar en la parte sur del Golfo de México
[Rubio, Triana, 2006].
3.2.2. Hidrologı́a
Por otro lado, en la parte superior de la cuenca, se ubica una de las zonas de ma-
yor precipitación del paı́s, superando los 4,000 mm anuales. En esta zona se registran
las precipitaciones mayores, y más cuando existen las combinaciones climatológicas de
sistemas tropicales, con la entrada de frentes frı́os o bien, cuando se combinan con co-
rrientes de aire frı́o; lo anterior repercute en severas inundaciones aguas abajo de la
cuenca, principalmente el estado de Tabasco.
En la planicie del Bajo Grijalva, la precipitación oscila entre los 1,700 mm y los 2,300
mm. La influencia de sistemas atmosféricos es similar a la parte alta del Bajo Grijalva,
con la diferencia que la precipitación disminuye, ya que no existen las combinaciones
mencionadas anteriormente, de acuerdo a los estudios realizados [Arellano, 1997].
Un dato importante, es el hecho de que la invación de masas de aire frı́o del norte y
húmedos tropicales del Atlántico, junto con el océano Pacı́fico provocan la mayorı́a de las
precipitaciones anuales de la región. En la temporada de verano, las lluvias se vuelven
intensas. Entre otoño e invierno, los nortes soplan, con lluvias de larga duración y se
presentan los torrenciales. Esto provoca que los rı́os y lagunas alcancen sus máximos
niveles en los meses de septiembre, octubre y noviembre, con lo que la planicie se vuelve
un espejo de agua. En esta época, las inundaciones son frecuentes provocando daños a
diversas actividades agrı́colas, al igual que las poblaciones que tienen sus asentamientos
en la llanura costera del norte [CFE, 2009].
Figura 3.17: Sistema hidroeléctrico de la cuenca del rı́o Grijalva [CFE, 1985]
El Plan Integral del Grijalva aprobó la construcción de diversas obras en el cauce del
respectivo rı́o Grijalva, dentro de las cuales destacan las siguientes [INE, 2005]:
2. Obras de defensa contra posibles inundaciones, las cuales incluyen bordos de pro-
tección, encauzamiento de corrientes y rectificación de cauces en la planicie costera
del Estado de Tabasco y dentro del Estado de Chiapas.
Figura 3.19: Sistema hidroeléctrico de la cuenca del rı́o Grijalva [INEGI, 2009]
Los objetivos principales de la construcción [CFE, 1976] de la presa fueron los si-
guientes :
Controlar las avenidas máximas registradas en el rı́o Grijalva, a fin de reducir gas-
tos, para no ocasionar inundaciones severas en los poblados aguas abajo, particu-
larmente en La Chontalpa.
Una vez que se definió el lugar de su construcción, en 1958 comenzaron los trabajos
con la construcción del camino de acceso. Al mismo tiempo, se iniciaron los estudios
geológicos de la Boquilla, ası́ como los levantamientos topográficos detallados, y la veri-
ficación por tierra de los levantamientos fotogramétricos del vaso. Terminando el año de
1959, se habı́an definido los diseños generales, y se habı́an concluı́do los estudios de
materiales aprovechables para la construcción de la cortina, ası́ como las obras auxilia-
res. A mediados de 1960, ya se contaba con los planos de la obra de desvı́o, la cortina,
3 DATOS DEL LUGAR DE ESTUDIO 102
Cortina. Para la selección del tipo de cortina se estudiaron las siguientes alternativas.
Una posible cortina de tipo arco-bóveda de concreto, una sección de gravedad de con-
creto o una cortina de enrocamiento con corazón impermeable de arcilla. Los estudios
finales, detallaron que se debı́a implementar ésta última. La altura máxima de la cortina
es de 139 m; la elevación de la corona es de 192 msnm; el ancho de la corona es de 10 m;
la longitud de la corona es de 480 m; y el volumen total de la cortina es de 5.08x106 m3 .
Los respectivos niveles del embalse son para el NAME 188.00 msnm, el NAMO 182.50
msnm y NAMINO 144.00 msnm. Para la construcción de dichos elementos, se cerraron
tres depresiones naturales mediante la construcción de diques. Uno de ellos cerca de
la cortina y los otros dos en el parteaguas que divide las cuencas del rı́o La Venta y el
Uxpanapa. Los diques uno y tres se diseñaron con secciones similares a la cortina. El
dique número dos se localiza en una falla geológica regional importante; éste se diseñó
con una sección homogénea de arcilla y enrocamientos de protección en los taludes y
con un drenaje eficiente en la cimentación del talud aguas abajo [SE, 2009].
cada uno constituido por una turbina tipo Francis, con capacidad de 180 MW cada una.
En cuanto a la potencia y a la generación, se cuenta con una capacidad instalada de
1080 MW, una generación media anual de 4929 GWh y un factor de planta de 51.95 %
[CFE, 2009].
En resumen las figuras 3.20, 3.21, 3.22 y 3.23 muestran las caracterı́sticas de las presa
La Angostura.
Obra de desvı́o. Para esta obra, el cauce del rı́o fue desviado a través de dos túneles re-
vestidos de concreto de 13 m de diámetro interior. Uno de ellos se localiza a la izquierda
del margen del rı́o y el otro a la derecha; además se cuenta con dos ataguı́as, localizadas
aguas arriba, de 6 0m de altura y 30m, la que se encuentra aguas abajo. Las ataguı́as se
construyeron de arena, grava y arcilla [CONAGUA, 2009].
Obra de excedencias. La obra de excedencias cuenta con dos vertedores que se ubican
en la margen izquierda del cauce. Esta obra, cuenta con dos canales abiertos dotados de
tres compuertas radiales cada uno. La longitud de los canales es de aproximadamente
800 m. El vertedor está diseñado para una descarga máxima de 6,000 m3 /s. La elevación
de la cresta del vertedor es de 519.60 msnm. La estructura terminal, que se localiza en
la salida de cada canal, se constituye por una cubeta de lanzamiento, tipo salto de esquı́
[CONAGUA, 1994].
3 DATOS DEL LUGAR DE ESTUDIO 105
En resumen, las figuras 3.24, 3.25, 3.26 y 3.27 muestran la información anterior.
La cortina tiene una altura de 146.7m a partir del punto más bajo de la cimentación.
La elevación de la corona es de 543 msnm, el ancho de la corona es de 10 m y la longitud
de la corona es de 323.5 m. La corona tiene un bordo libre además de 3.50 m. La corona
fue construida a base de enrocamientos que se obtuvieron al momento de realizar las ex-
cavaciones de las obras a cielo abierto, ası́ como los canales de vertedores, además de
las excavaciones de obras subterráneas. La cortina también está compuesta por grava y
arena de los aluviones del rı́o, aguas abajo de la obra. El núcleo central de la cortina, es
impermeable y está formado por arcilla compactada. La cortina tiene un volumen total de
4.2 millones de m3 [Comisión del Grijalva, 1964].
Dentro del Plan Integral del Rı́o Grijalva, la Presa Chicoasén fue la tercera presa
en construirse. Dadas sus condiciones hidrológicas, y sus caracterı́sticas topográficas,
económicas y geológicas, permiten que en este sitio, la central hidroeléctrica sea la que
genera mayor energı́a eléctrica anual en el paı́s.
Obra de desvı́o. Para realizar la obra de desvı́o, el proceso constructivo no presentó di-
ficultades, ya que las condiciones geológicas permitieron la construcción de un túnel
auxiliar de 343 m de longitud de sección sin revestimiento, para poder desviar el cauce.
Esto ya que el caudal del rı́o se volvı́a reducido, y con ello se pudo construir la ataguı́a
aguas arriba, por la desviación previa. El túnel auxiliar tiene un ancho de 7 m y una altura
de 6 m, y se sitúa en la margen izquierda del cauce [CFE, 2009].
La presa hidroeléctrica Chicoasén tiene un gasto medio anual de 377 m3 /s. El embal-
se tiene una capacidad total de 1443 millones de m3 . El labio superior de la compuerta
se localiza a 394 msnm y la cresta del vertedor a 373 msnm con una longitud de 76 m,
3 DATOS DEL LUGAR DE ESTUDIO 110
La casa de máquinas tiene una longitud de 199 m, un ancho de 20.50 m y una altura
de 43 m. Para su acceso se cuenta con un túnel de 800 m de longitud aproximadamen-
te y una sección portal de 8.45x9.40 m. La obra civil se construyó en dos etapas para
poder alojar hasta ocho grupos de turbina generador; en la primera se instalaron 5 tur-
bogeneradoras tipo Francis de 300 MW cada una, y en la segunda se instalaron 3 más,
de igual capacidad. Actualmente, Chicoasén tiene una capacidad efectiva instalada de
generación de 2400 MW, pero está diseñada para 3324.87 GW [INE, 2005].
La cortina tiene una altura máxima de 53 msnm, una elevación de corona de 98 msnm,
un ancho de corona de 8m, una longitud de corona de 560 m y un volumen total de
3.24x106 m3 [Comisión del Grijalva, 1964]. El nivel del NAME es de 93.50 msnm y el NA-
MO es de 87.40 msnm.
Obra de desvı́o. Ya que el brazo izquierdo se bloqueó, el cauce del rı́o fluı́a en dirección
del brazo derecho. Con esto, se pudo avanzar con la cortina por la margen izquierda. En
dirección paralela al avance de la cortina, se excavó un canal de desvı́o a cielo abierto
de 35 m de ancho por la margen izquierda, con esto se permitió concluir con la cortina
[CFE, 1976].
de 18.7 m3 , siendo 9.350 m3 por cada uno de los vertedores [CFE, 1988].
Peñitas se localiza a una distancia de 72 km sobre los ejes de cortina aguas abajo de
Malpaso. Es la presa más pequeña en longitud, cuya área de cuenca está en los 1358
km2 . El labio superior de la compuerta está a 91.13 msnm. El aprovechamiento tiene una
capacidad total aproximada de 1630 millones de m3 , una capacidad útil de generación de
130 millones de m2 y una capacidad para control de avenidas de 961 millones de m3 . La
cresta del vertedor está a 76.50 msnm. Su aprovechamiento tiene una capacidad efectiva
instalada de generación de 420 MW, y su generación neta de diseño es de 1219.04 GWh
[SRH, 1971].
La figura 3.36 muestra los datos hidrológicos de todas las presas. Se observa que
en orden de tamaño de las cuencas, Peñitas es la mayor, seguida por Malpaso, luego
La Angostura y finalmente Chicoasén. Dadas además, las condiciones geográficas y la
topografı́a de la región, las adecuaciones técnicas permiten que Malpaso tenga el mayor
gasto medio anual, seguido de Chicoasén, luego La Angostura y finalmente Peñitas. No
obstante, este orden no se presenta para los gastos máximos, ya que existen avenidas
extraordinarias y dependen también de la ubicación de las cuencas, y de quién se en-
cuentre aguas abajo o aguas arriba; es por ello, que Malpaso presenta el mayor gasto
máximo registrado, seguido de Chicoasén, luego Peñitas y finalmente La Angostura. Las
precipitaciones medias anuales son iguales para todas las presas y los niveles de eva-
poración son mayores en La Angostura, seguido de Malpaso, luego Peñitas y finalmente
Chicoasén [SE, 2009]. Asimismo la figura 3.37 muestra los datos del embalse de todas
las presas.
Figura 3.37: Datos del embalse de todas las presas [Arellano, 1997]
Figura 3.38: Datos de la obra de toma de todas las presas [Arellano, 1997]
Se observa en la figura 3.38 los datos de la planta hidroeléctrica de todas las presas.
Malpaso cuenta actualmente con el mayor número de turbinas en función. Los estudios
hechos de manera previa, consideran que algunas de las presas puedan permitir el fun-
cionamiento de turbinas adicionales para una mayor generación de energı́a. En la figura
3.39 se observan los datos de la cortina de cada una de las presas.
Chicoasén es la presa con la mayor caı́da de agua, permitiendo una gran generación
hidroeléctrica, pese a que el embalse no sea tan grande como las otras presas hidro-
eléctricas. Los máximos volúmenes turbinables permitirán establecer las restricciones
fı́sicas en cuanto a la cantidad máxima de agua permitida, y por lo tanto, afectará a la
generación de energı́a.
3 DATOS DEL LUGAR DE ESTUDIO 115
Figura 3.40: Datos del vertedor de todas las presas [Arellano, 1997]
La presa Chicoasén es la que tiene la mayor altura, mientras que La Angostura tiene
la mayor elevación de la corona dado su situación topográfica. Chicoasén es además, el
embalse con el mayor bordo libre, ası́ como la que mayor cantidad de material se utilizó
para la cortina impermeable, la zona de transición y el enrocamiento, ası́ como el material
para el filtro. Finalmente, Chicoasén es también quien presenta las mayores dimensiones
de los vertedores como obra de excedencia, mientras que Malpaso no tiene vertedores
para este fin.
Del año 2000 al 2013, el impacto promedio anual de los desastres de origen hidrome-
teorológico en Chiapas fue de 2,905 millones de pesos, mientras que en el perı́odo de
1928 a 1999, fue de 572 millones de pesos. En términos del PIB estatal, en las décadas
de los años 80 y 90, los impactos de los desastres naturales representaban el 0.20 % del
PIB, mientras que a partir del año 2000 hasta el 2013, el porcentaje promedio se elevó a
1.12 % del PIB estatal [CENAPRED, 2014].
Para el año 2013, el impacto de los desastres naturales representó el 0.37 % del PIB
nacional. Desde 1997, exceptuando el 2010, el año 2013 registró la cifra más elevada
de impactos provocados. El impacto socioeconómico de los desastres toma relevancia
porque el monto llega hasta los 61,009 millones de pesos a nivel nacional, lo que sig-
nifica 203 veces el presupuesto del Fondo para la prevención de desastres naturales y
fideicomiso preventivo. Desde hace 15 años, 9 de cada 10 desastres naturales son de
origen hidrometeorológico, con ello, el 92 % de los impactos socioeconómicos fueron de
este tipo. Un ejemplo en concreto, fue el huracán Stan, el cual alcanzó el 8.87 % del PIB
estatal del 2005 [OCDE, 2013].
Se recuerda que una polı́tica de operación, también llamada una extracción, es una
regla que indica la cantidad de volumen de agua que debe extraerse en un futuro próxi-
mo, al tomar en cuenta el estado presente del sistema. Dicho estado, es definido en
términos de variables fı́sicamente medibles, como es el volumen almacenado en la presa
en diferentes épocas del año.
A lo largo de los años se ha tratado de determinar estas polı́ticas por diversos modelos
matemáticos, y distintas técnicas de investigación, como la investigación de operaciones,
la programación lineal, y actualmente la programación dinámica estocástica. Es ya co-
nocido el hecho, que el principio de optimización de Bellman establece que no importa
3 DATOS DEL LUGAR DE ESTUDIO 119
cual sea el estado inicial o la etapa inicial de un proceso secuencial de decisiones. Con
esto, las decisiones que restan deben corresponder a una polı́tica óptima con respecto
al estado que resulta de la primer decisión tomada [Bellman, 1957].
Para este año, se decidió tomar en cuenta evitar los derrames y posibles déficits, ası́
como tomar en cuenta la generación de energı́a a largo plazo. Es por ello, que la función
objetivo se modeló en forma mensual, a fin de determinar las polı́ticas de operación del
sistema de presas [Domı́nguez, Mendoza y Contrera, 1998]. Se propuso una función
objetivo de la forma
B = b − Cdef ∗ DEF − Cderr ∗ DERR, (3.1)
donde:
b = 100 + α[ EEdem
ge
− 1] si la energı́a generada (Egen ) es menor que la energı́a deman-
dada (Edem ); en caso contrario, se hace b = 100 + β[ E−genEdem
− 1], donde experimen-
talmente se determinaron los valores de α = 20 y β = 150.
Cdef y Cderr son coeficientes constantes que se consideran para penalizar el bene-
ficio B debido a la ocurrencia de un déficit o de un derrame.
manejo de volúmenes de ingreso de las presas y de extracción por las turbinas signifi-
cativamente mayores que su volumen útil. En este sentido, La Angostura funciona como
almacenamiento de Chicoasén y Malpaso, como almacenamiento de Peñitas. Ası́ pues,
la polı́tica de operación de las presas menores, están sujetas a tratar de turbinar volúme-
nes de ingreso sin modificar el nivel del agua en el vaso y se toman en cuenta el manejo
de las compuertas del vertedor [Domı́nguez y Mendoza, 2000].
Otra de las consideraciones previas, fue definir una polı́tica de operación inicial para
La Angostura. Esta polı́tica fue determinada, tomando en cuenta el volumen de almace-
namiento disponible en cada presa, y la ubicación en el sistema Grijalva. Luego, los ingre-
sos de la presa Malpaso se determinaron adicionando las aportaciones por cuenca propia
y las descargas obtenidas al simular la operación de La Angostura. Se realizó un análisis
por separado y luego en cascada, combinándose ambas presas [Domı́nguez, 1989].
Los datos iniciales para la determinación de las primeras polı́ticas óptimas son las
caracterı́sticas previamente mencionadas de las presas, las cuales se han detallado para
cada uno de los elementos que conforman el aprovechamiento del Grijalva. Tales ca-
racterı́sticas son la elevación inicial del vaso, la potencia instalada de cada una de las
turbinas, el gasto de diseño, el volumen de almacenamiento inicial, el máximo volumen
mensual turbinable, la carga de diseño, las curvas de elevaciones-capacidades, los nive-
les del NAMO, NAMINO, NAME, detalles de la lámina de evaporación neta, los niveles
de desfogue, las curvas de elevación-volumen, las curvas elevación-área, principalmente.
Los registros históricos hasta la década de los años 90’s, sirvieron para agrupar
meses con caracterı́sticas estadı́sticas homogéneas, según la media y la desviación
estándar. Estos modelos probabilı́sticos sirvieron de pauta para las distribuciones de pro-
babilidad, ya que los ingresos son aleatorios. Además, se discretizó el problema a partir
de incrementos de volumen de ingreso, denotados por ∆V , de acuerdo al volumen útil
de la presa, creándose intervalos desde 1 hasta un máximo llamado N S [Arzate, 2000].
Se estudiaron los niveles del Grijalva en Tabasco, por efecto de descargas de la obra
de toma, los cuales no fueron afectados en forma significativa para gastos inferiores a
1000 m3 /s. Se realizó además un estudio del efecto de los niveles de Malpaso por des-
cargas de La Angostura. Esto llevó a la determinación de factores de ajuste en las des-
cargas de La Angostura, dependiendo si el volumen útil en Malpaso excedı́a a la media al
inicio de cada mes. Con estos nuevos cálculos, los derrames originales de 4474 millones
de m3 disminuyeron a 772 millones de m3 [Calva, 1993].
Para el año de 1998, ya se tenı́a una tabla de polı́ticas de operación para el siste-
ma trabajando en conjunto, que se desprendı́an de las originales de 1993. Sin embargo,
las avenidas extraordinarias de 1999, obligaron a replantear el problema, puesto que
los almacenamientos llegaron a 3050 millones de m3 y 356 millones de m3 para Malpa-
so y La Angostura respectivamente. Cabe mencionar, que hasta esta fecha, todos los
estudios y simulaciones eran teóricas, y aún no se ponı́an en marcha en las presas
[Fuentes, Sánchez, 1993].
Se observó que a partir de los valores propuestos de 10,000 para los costos de de-
rrames, éstos ya no disminuı́an significativamente, lo que representaba que la energı́a
no tenı́a un aumento considerable. Esto se debió a las autocorrelaciones no tomadas
en cuenta entre los meses y los ingresos. Por ello, se construyeron histogramas de pro-
babilidades desfasados, para tomar en cuenta correlaciones cruzadas en cada uno de
los meses, principalmente, cuando las lluvias aumentaban, esto es, a partir de agosto.
El problema de las autocorrelaciones también afectaba en los meses de transición entre
el estiaje y la época de grandes avenidas; puesto que se contaban con niveles bajos y
extracciones mı́nimas durante todo un mes, para luego pasar a extracciones máximas
con niveles aumentando. Esto condujo a reducir los intervalos de tiempo de un mes a tan
sólo 15 dı́as en las simulaciones [Domı́nguez y Mendoza, 2000].
Las corridas para este año se efectuaron por perı́odos de registros históricos de 1959
a 1999, y de 1977 a 1999, interpolando algunos valores desconocidos, y otros que no se
tenı́an los detalles. También se efectuaron corrimientos hacia la derecha de los histogra-
mas de probabilidades para los efectos de las autocorrelaciones. En total, se realizaron
13 corridas con los ajustes anteriores [Domı́nguez y Mendoza, 2000].
3 DATOS DEL LUGAR DE ESTUDIO 124
Figura 3.46: Corridas efectuadas por perı́odo [Domı́nguez, Mendoza y Arganis, 2001]
Se observa en la figura 46, que conforme se aumentan los costos por derrames, la
probabilidad por derrame disminuye. De manera concreta, con penalizaciones bajas de
10 unidades, los derrames son excesivos, y al aumentarse a 1000 unidades, estos se
disminuyen por la mitad. Cuando se tienen penalizaciones de 10 000 unidades, los de-
rrames disminuyen en forma considerable, y la energı́a eléctrica no se ve afectada en
forma significativa.
También se definió una polı́tica de operación inicial, de tal manera que no se presen-
taran derrames en el sistema, esta polı́tica servı́a para tener una base para los valores
máximos de energı́a generada. Para ello, se trazó una curva de volúmenes acumula-
dos en función del tiempo. Sobre esta curva se graficaron las capacidades posibles de
generación como funciones lineales. Las pendientes de estas rectas de capacidades
no debı́an ser excedidas, pues fı́sicamente representaba que los niveles sobrepasaban
al NAMO, y se presentaban derrames. A fin de evitar tales situaciones, se proyectaron
3 DATOS DEL LUGAR DE ESTUDIO 125
rectas tangentes sobre los ingresos quincenales a partir de las rectas de generaciones
máximas. Esto permitió que los niveles no excedieran al NAMO y en forma natural los
ingresos siempre quedaran por debajo de las extracciones, manteniendo el nivel de las
presas sin presentar derrames [Arganis, 2004].
Este proceso se hizo para el perı́odo 1990 al 2001. Con estas consideraciones, se
llegó a la conclusión de mantener extracciones iguales a los ingresos para las primeras
14 quincenas. A partir de la quincena 15, las extracciones serán mayores a los ingresos
durante 3 quincenas, para nuevamente tomar extracciones iguales a los ingresos. Esto
permitió que los niveles entre NAMO, NAMINO y NAME se mantuvieran en las elevacio-
nes de acuerdo con la CONAGUA, y cumpliendo con las demandas mı́nimas de energı́a.
Una vez tenidas las extracciones, se hicieron los cálculos de las simulaciones de las
energı́as obtenidas para cada una de las presas [Arganis, 2004].
Las simulaciones de estos perı́odos tuvieron las siguientes posibilidades: sin consi-
derar la energı́a de pico y de base, sin considerar restricciones en extracciones mı́nimas,
considerando la energı́a pico y de base, y considerando las restricciones en las extrac-
ciones mı́nimas. Como resultado, cuando no se consideran restricciones en extracciones
mı́nimas, no es necesario tomar en cuenta los picos, pero entonces sı́ existen derra-
mes. Al considerar extracciones mı́nimas, La Angostura permaneció prácticamente igual,
mientras que la energı́a de Malpaso disminuyó. Al retringir a Malpaso, ocurre lo contrario,
aumenta la energı́a de Malpaso, disminuyendo la de La Angostura, pero los derrames de
Malpaso aumentan. Si ambas presas no se restringı́an, y se consideraban los picos, se
mejoraba la situación de evitar derrames en ambas presas. Cada una de estas simulacio-
nes trajo en consecuencia, buscar el equilibrio entre mantener las demandas energéticas
de CFE cumplidas, y evitar las posibles inundaciones aguas abajo en Tabasco.
La curva guı́a también llamada curva ı́ndice es la curva de seguridad que sirve de
referencia para saber si una presa está en riesgo o no, y son establecidas por el Comité
Técnico de Operación de Obras Hidráulicas. La CONAGUA solicita no exceder los niveles
de almacenamiento que se establecen en esta curva. Los diferentes algoritmos de optimi-
zación empleados hasta la fecha toman en cuentan los efectos de la curva guı́a a través
de coeficientes de penalización cuando se rebase alguno de estos niveles [CNA, 1997].
En caso de que se rebase la curva guı́a en cada uno de los embalses, se agrega
al algoritmo de optimización la lectura de un coeficiente de penalización, para tomar
en cuenta a la curva guı́a expresada generalmente en forma discreta y por etapa. Si
el almacenamiento final de la presa es menor o igual a la curva guı́a, el beneficio por
generación se expresa de la siguiente manera:
donde se tienen los mismos parámetros que la ecuación anterior, además de que CGGpresa
es la penalización por exceder la curva guı́a [Domı́nguez, 1998].
Históricamente, la CFE ha ensayado distintas curvas guı́a conforme pasan los años
de registros hidrológicos, ası́ como las afectaciones de avenidas extraordinarias. En el
año 2008, CFE propuso la curva guı́a llamada POL CG2008 utilizando los registros de
1959-2007. Esta curva guı́a permite una generación de 513.06 GWh/quincena, y se ob-
tiene un derrame total de 19.4 millones de m3 para la presa Malpaso. Este derrame se
distribuyen como 7.22 millones de m3 en la primera quicena de octubre de 1970 y 12.18
millones de m3 en la primera quincena de 2003 [CONAGUA, 2009].
En 2010, se propuso la curva guı́a POL XIV. Con esta polı́tica se obtiene una gene-
ración de 486.82 GWh/quincena, y no existen ni derrames ni déficits en las presas. Cabe
3 DATOS DEL LUGAR DE ESTUDIO 127
mencionar que se consideraron en esta polı́tica las láminas de evaporación que se regis-
traron entre el 22 de octubre de 2009 y el 22 de octubre de 2010. Con esto, la energı́a
disminuyó 14.71 GWh/quincena en La Angostura y 11.53 GWh/quincena en Malpaso.
Los almacenamientos aumentaron a 1294.18 y 618.53 millones de m3 en La Angostu-
ra y Malpaso, respectivamente. Para esta polı́tica se tomaron en cuenta los registros de
1959-2008. En esta polı́tica hay mejorı́as respecto a la polı́tica POL CG 2008, sin embar-
go, los niveles quedan por arriba del NAME, representando riesgos importantes en caso
de avenidas extraordinarias [Alegrı́a, 2010].
La polı́tica POL XIV sirvió de base para modificar y mejorar los resultados. A partir del
análisis de los volúmenes históricos de ingreso por cuenca propia a las presas La An-
gostura y Malpaso, se observó que existe una relación entre el incremento de volumen
del almacenamiento en las presas durante la época de avenidas y la desviación estándar
de los volúmenes de ingresos. Por esto, en las polı́ticas posteriores se consideró es-
te parámetro estadı́stico, para poder definir las curvas guı́a mediante un procedimiento
sistemático. Para este fin se realizó un análisis de las elevaciones alcanzadas en forma
histórica, de las cuales se tenı́a registro y se buscó las tendencias del Sistema Hidro-
eléctrico Grijalva. Este registro se tomó desde 1998 hasta el 2008 [CONAGUA, 2009].
4. Metodologı́a
En este capı́tulo se describe la metodologı́a planteada para encontrar las polı́ticas
óptimas de operación de un sistema de embalses en serie como se observa en las figuras
4.47, 4.48, 4.49 y 4.50.
4.1.2. Histogramas
Posteriormente, se dividió el año en 6 etapas, tal como lo realiza la CONAGUA, de
acuerdo con los registros históricos. Para cada mes, se realizó un análisis de frecuen-
cias con distintas distribuciones continuas de probabilidad y se eligió la que menor error
hubiera presentado. Por cada etapa, se calculan los ingresos totales por año y se orde-
naron cronológicamente. Se obtuvo la frecuencia, el porcentaje acumulado, la frecuencia
relativa y se elaboran los histogramas, con la respectiva curva de ajuste de la distribución
obtenida.
4.2.2. Restricciones
Se añadieron las restricciones de los niveles de desfogue de cada presa para ser con-
siderado en el modelo de restricción estocástico. Además, para cada etapa y para cada
una de las presas, se calcularon los lı́mites máximos de las extracciones.
Por otro lado, se recopiló la información de las curvas guı́a dadas por la CONAGUA
presentada por quincenas en el año. La información se presentó tanto en unidades de
volumen, como en unidades de elevaciones para cada una de las presas. Se elabora-
ron los gráficos de las curvas guı́as en ambas unidades descritas con aterioridad. Cada
curva guı́a se modeló analı́ticamente, y se realizaron interpolaciones a fin de tener un
conjunto de ecuaciones polinomiales para cada una de las curvas presentadas.
con n, representando la etapa y An (x) la extracción que se realizará. Esto se hizo para
ambas presas.
pues, el kernel depende de la etapa en cuestión, a través de una función f (x1 , a), donde
x1 es el estado y a la extracción. La ley de transición está dada por la ecuación
fk (x1 , a) = x1 + m − a, (4.5)
Además, las mi son variables aleatorias distribuidas en forma normal, de acuerdo con
un análisis de diferentes distribuciones posibles, buscando la que minimizara el error.
Para cada una de ellas se obtuvo su respectiva desviación estándar, ası́ como su estan-
darización
X −µ
Z= , (4.6)
σ
donde N (µ, σ 2 ) representa la distribución de la variable aleatoria X con distribución ori-
ginal N (0, 1).
Las primeras condiciones representan los derrames de ambas presas, y las segundas
condiciones los déficits de ambas presas, también.
donde d2 es el derrame, x2 el estado final y el lı́mite máximo está dado por el espacio de
estados. Análogamente, para el caso de déficits, se calcularon mediante la expresión
d1 = |x2 |, (4.9)
si x2 < 0, con d1 el déficit y x2 el estado final. Además la extracción que pasa por las
turbinas se calculó como
ET = a − d1 , (4.10)
con a la extracción y d1 el déficit. Cuando existieron los derrames, estos formaron parte
de las extracciones de la presa Malpaso; además, el volumen turbinado se añadió a estos
derrames, y formó parte de las extracciones.
a1 = d2 + ET, (4.11)
El nivel de salida x2 se comparó con los niveles aceptados por las curvas guı́a de la
CONAGUA. En caso de no corresponder, se corrigieron agregando los niveles iniciales y
restando los niveles entre el NAMINO y el nivel medio de desfogue. Si existieron derra-
mes o déficits, se aplicaron las penalizaciones mencionadas anteriormente.
donde e.m. es el estado máximo y se calculó con las ecuaciones de estado iniciales.
Nótese que BEM es una función que depende del estado inicial x1 , de la extracción a2 y
del ingreso aleatorio m.
4.6. Beneficios
La función de beneficio esperado se derivó con respecto a la extracción m, mante-
niendo constante los demás valores, es decir si
Z 6567
BEM = [r(x1 , a2 ) − c(x1 , a2 )] ∗ P (m)dx1 , (4.17)
0
entonces
d
BEM (x1 , a2 , m) = 0, (4.18)
da2
luego Z 6567
d
[r(x1 , a2 ) − c(x1 , a2 )] ∗ P (m)dx1 = 0. (4.19)
da2 0
4 METODOLOGÍA 139
Para la presa 2, las ecuaciones son similares, con la condición de que el ingreso tuvo
los volúmenes de cuenca propia, mas el volumen de la presa 1, que se extrajo, es decir
m = m2 + a2 , (4.21)
Al variar las extracciones en cada una de las presas, a partir de las ecuaciones ante-
riores se calcularon los beneficios esperados para un horizonte de vida útil de 50 años.
Estos beneficios se presentaron en unidades de GWhr y las extracciones en unidades de
hm3 .
5. Resultados
5.1. Registros históricos de escurrimientos
Uno de los factores que distingue el análisis determinista de un análisis estocástico
es la incertidumbre de las mediciones. El segundo de ellos generalmente se modela a
través de distribuciones de probabilidad. Como se mencionó en los capı́tulos anterio-
res, las ecuaciones diferenciales ordinarias tienen limitantes cuando las variables que
intervienen son desconocidas y más aún, son de carácter aleatorio, dando lugar a las
ecuaciones diferenciales estocásticas.
Figura 5.64: Registro de ingresos quincenales de Malpaso de 1952 a 1971 (hm3 ), los
meses de enero a abril.
Figura 5.65: Registro de ingresos quincenales de Malpaso de 1952 a 1971 (hm3 ), los
meses de mayo a agosto.
5 RESULTADOS 148
Figura 5.66: Registro de ingresos quincenales de Malpaso de 1952 a 1971 (hm3 ), los
meses de septiembre a diciembre.
Figura 5.67: Registro de ingresos quincenales de Malpaso de 1972 a 1991 (hm3 ), los
meses de enero a abril.
5 RESULTADOS 149
Figura 5.68: Registro de ingresos quincenales de Malpaso de 1972 a 1991 (hm3 ), los
meses de mayo a agosto.
Figura 5.69: Registro de ingresos quincenales de Malpaso de 1972 a 1991 (hm3 ), los
meses de septiembre a diciembre.
5 RESULTADOS 150
Figura 5.70: Registro de ingresos quincenales de Malpaso de 1992 a 2011 (hm3 ), los
meses de enero a abril.
Figura 5.71: Registro quincenales de Malpaso de 1992 a 2011 (hm3 ), los meses de mayo
a agosto.
5 RESULTADOS 151
Figura 5.72: Registro de ingresos quincenales de Malpaso de 1992 a 2011 (hm3 ), los
meses de septiembre a diciembre.
Figura 5.73: Registro de ingresos quincenales de Malpaso de 2012 a 2017 (hm3 ), los
meses de enero a abril.
Figura 5.74: Registro de ingresos quincenales de Malpaso de 2012 a 2017 (hm3 ), los
meses de mayo a agosto.
5 RESULTADOS 152
Figura 5.75: Registro de ingresos quincenales de Malpaso de 2012 a 2017 (hm3 ), los
meses de septiembre a diciembre.
Los años de 1952 hasta 1958 se registraron en forma mensual y en forma quincenal
para cada una de las cuencas de cada presa. A partir del año 1959, los registros se
comenzaron a realizar en forma diaria, esto significa, que los datos que aparecen en las
tablas anteriores, fueron obtenidos realizando las sumas parciales por cada una de las
quincenas de cada mes. Los registros mensuales llegan hasta el año de 1999. Durante
la toma de registros, no se cuenta con los años 1975 y 1976. Estos datos fueron incor-
porados interpolando con los registros de datos que ya se contaban.
A partir del año 2000, los registros fueron hechos en forma diaria para cada una de
las presas, es decir, para cada una de las cuencas, tanto para la Presa La Angostura,
como para la Presa Malpaso, a partir de las cuencas propias y los escurrimientos de
Chicoasén, para el caso de Malpaso. Los datos del 2006 al 2017 fueron recopilados di-
rectamente con la Comisión Nacional del Agua, y por la Gerencia de Despacho Eléctrico
de la Comisión Federal de Electricidad, a través de registros diarios de escurrimientos,
evaporaciones e infiltraciones. Los datos que aquı́ se presentan ya son los escurrimientos
netos por cada una de las cuencas de manera individual, presentados en forma quince-
nal.
Para cada una de las etapas presentadas, se requiere el cálculo de las respectivas
distribuciones de probabilidades. Esto se hará mediante las respectivas frecuencias de
cada uno de los registros históricos por los meses que involucren a una etapa.
Figura 5.91: Ingresos totales para los meses de noviembre y diciembre. Años 1952-1990
5 RESULTADOS 159
Figura 5.92: Ingresos totales para los meses de noviembre y diciembre. Años 1991-2017
Se recuerda que las variables numéricas son discretas si su conjunto de valores po-
sibles es finito o se puede enumerar en una sucesión infinita; mientras que las variables
numéricas son continuas si sus valores posibles abarcan un intervalo completo sobre la
recta real.
A partir de estos datos se calcularon las frecuencias relativas para la etapa 1 como se
observa en la figura 5.93.
5 RESULTADOS 160
Las distribuciones de este tipo no son necesariamente simétricas. Además, los paráme-
tros se pueden establecer de la siguiente manera:
n
" n #1/2
X ln xi X (ln xi − α)2
α= , β= . (5.2)
i=1
n i=1
n
De acuerdo con los respectivos registros históricos se realizó una distribución normal,
ya que es la que mejor ajuste presentó, obteniendo la respectiva tabla de probabilidades
de la etapa que se muestra en la figura 5.96.
5 RESULTADOS 162
Se observa en la figura 5.97 un ingreso de 4734.20 hm3 en el registro del año 2005,
el cual fue ocasionado por el Huracán Stan a finales del mes de octubre. Este dato se
5 RESULTADOS 164
observa también en la figuras 5.98 con frecuencia relativa de 0.0154 y en la figura 5.99,
como un dato aislado sobre el histograma. Este hecho se ve reflejado en la figura 5.100,
donde existe cierta asimetrı́a hacia la derecha, en la tabla de probabilidades, la cual se
ve en la figura 100.
Figura 5.109: Ingresos totales para los meses de Junio y Julio. Años 1952-1990
Figura 5.110: Ingresos totales para los meses de Junio y Julio. Años 1991-2017
5 RESULTADOS 170
Figura 5.114: Ingresos totales para los meses de Enero a Mayo. Años 1952-1971
5 RESULTADOS 172
Figura 5.115: Ingresos totales para los meses de Enero a Mayo. Años 1972-1990
Figura 5.116: Ingresos totales para los meses de Enero a Mayo. Años 1991-2003
5 RESULTADOS 173
Figura 5.117: Ingresos totales para los meses de Enero a Mayo. Años 2004-2017
Figura 5.121: Ingresos totales para los meses de noviembre y diciembre. Años 1952-1990
Figura 5.122: Ingresos totales para los meses de noviembre y diciembre. Años 1991-2017
5 RESULTADOS 176
Figura 5.138: Ingresos totales para los meses de Junio y Julio. Años 1952-1990
5 RESULTADOS 184
Figura 5.139: Ingresos totales para los meses de Junio y Julio. Años 1991-2017
Figura 5.143: Ingresos totales para los meses de Enero a Mayo. Años 1952-1971
Figura 5.144: Ingresos totales para los meses de Enero a Mayo. Años 1972-1990
5 RESULTADOS 187
Figura 5.145: Ingresos totales para los meses de Enero a Mayo. Años 1991-2003
Figura 5.146: Ingresos totales para los meses de Enero a Mayo. Años 2004-2017
5 RESULTADOS 188
Este tipo de anomalı́as climatológicas está muy ligada a los fenómenos de cambio
climático que la Organización Meteorológica Mundial ha venido dando seguimiento. Los
cambios en la temperatura del océano provocan este tipo de manifestaciones hidrológi-
cas, tales como los fenómenos conocidos de El Niño y La Niña.
1999. Este año ocurrió el Frente Frı́o No. 4 y la entrada de humedad del Golfo y del Cari-
be, lo que originó lluvias extraordinarias en la región sureste los dı́as 1 y 2 de octubre, con
avenidas de 8102.86 m3 /s como gasto pico en el hidrograma de Peñitas. La operación de
5 RESULTADOS 192
las centrales revistieron singular importancia tanto en septiembre como en octubre, por
las intensas precipitaciones y escurrimientos que ocurrieron en la cuenca de influencia
de estos aprovechamientos. Lo mismo ocurrió en la planicie del estado de Tabasco. Los
almacenamientos fueron extraordinarios, como volúmenes no antes almacenados en los
embalses de estas centrales, al igual que la presencia de inundaciones en la zona baja
del rı́o Grijalva [La Jornada, 2012].
2007. El 21 de octubre del 2007 se tuvo la afectación del Frente Frı́o No. 4. Este se exten-
dió desde el oriente de los Estados Unidos hasta el sur del estado de Veracruz, con un
movimiento este-sureste. A partir del 26 de octubre se mantuvo estacionario y se exten-
dió al Atlántico norte, luego a las costas orientales de los Estados Unidos, pasando por el
norte del Mar Caribe, hasta llegar al noreste de la Penı́nsula de Yucatán. Mientras tanto,
al sureste del Mar Caribe, se desarrolló la Depresión Tropical No. 16, hasta convertirse
en tormenta tropical [OMM, 2018].
En los dı́as 28, 29 y 30, el frente frı́o estacionario interactuó con la tormenta tropical,
formando un flujo de humedad hacia el sureste del paı́s, en los estados de Tabasco y
Chiapas. La depresión tropical se intensificó y se convirtió en Tormenta Tropical Noel, del
Atlántico. Peñitas presentó un gasto máximo instantáneo de 5,005.7 m3 /s y un nivel de
91.32 msnm, con lo que las compuertas del vertedor descargaron un volumen total de
233.22 millones de m3 [Samuelson, 2006].
2010. En el año 2010 se presentó un fenómeno hidrometeorológico poco común, que los
especialistas definieron como monzón del sureste. Un monzón son vientos estacionales
basados en diferentes capacidades calorı́ficas que cambian de dirección, provenientes
de una baja de presión hacia una alta presión, debido al gradiente de temperatura. Esto
origina lluvias de hasta 3 meses de duración [USDA, 2013].
Por otra parte, el vaso de Chicoasén es más reducido que el vaso de La Angostura
y Malpaso. Las aportaciones de los escurrimientos de junio a septiembre fueron de tipo
húmedo. En los meses de julio y agosto se rebasaron los niveles máximos registrados,
llegando a los 1,239.65 millones de m3 ; nivel que no se excedı́a desde el año de 1973
[Samuelson, 2006].
El dı́a 20 de julio, se presentó un gasto pico del orden de 736 m3 /s, producto de las
precipitaciones generadas por el paso de la onda tropical No. 16. Para el mes de agosto,
a partir del dı́a 17 y hasta final de mes, se presentó un incremento en las aportaciones
medias diarias que en promedio fueron del orden de 588 m3 /s. Ası́ pues, el 19 de agosto,
se presentó un gasto máximo de 840.62 m3 /s [OCDE, 2013]. La figura 5.105 presenta
estos datos.
5 RESULTADOS 195
Además, se tienen registros de que México vivió una de las peores sequı́as en 70
años, en el año de 2012. Esto ocasionó pérdidas económicas en la producción agrı́cola
por encima de los 16,000 millones de pesos. Esta afectación se vivió en un 70 % del paı́s,
siendo los estados del norte los principales afectados, entre ellos, Sonora, Chihuahua,
Coahuila, Nuevo León, San Luis Potosı́, Tamaulipas, Durango, Zacatecas y Aguascalien-
tes [La Jornada, 2012].
El INEGI reportó las pérdidas en un 10 % del producto interno bruto nacional, inclu-
yendo 9,000 millones de pesos en cultivos de maı́z y 6,000 millones de pesos en cultivo
de frijol. Además, 48 millones de mexicanos padecieron las consecuencias de la sequı́a
en zonas áridas, semiáridas y subhúmedas secas, afectando también al sureste del paı́s
[La Crónica, 2012].
2013. En el año 2013, se presentó el huracán Manuel. Este fue el tercero de los ciclones
de la temporada del Pacı́fico que tocó tierra en México. El dı́a 15 de septiembre, cerca
de las 14.00 horas, la tormenta tropical Manuel, tocó la costa occidental del estado de
Colima, cerca de la ciudad de Manzanillo, con vientos máximos sostenidos de 100 km/h
y rachas de 130 km/h. A las 16.00 horas, ya se habı́a aproximado a 3 km al noroeste de
Lima y a 20 km al norte de Manzanillo. Una vez en tierra, perdió fuerza, y a las 22.00
horas, se habı́a convertido en depresión tropical, en el estado de Jalisco, en la población
de El Limón. En este lugar, los vientos máximos sostenidos fueron de 55km/h, con ra-
chas hasta de 75km/h y un desplazamiento hacia el noroeste de 13km/h [OMM, 2018].
La trayectoria de Manuel se observa en la figura 5.157.
Luego de degradarse a una baja de presión, Manuel salió al mar nuevamente, y ge-
neró fuerza, hasta volverse depresión tropical el 17 de septiembre por la tarde. Para el
dı́a 18 de septiembre, por la mañana, se convirtió en tormenta tropical, y por la tarde de
ese dı́a se intensificó a huracán, alcanzando vientos máximos sostenidos de 120km/h
con rachas de 150km/h. Para el dı́a 19 de septiembre, tocó la costa de Sinaloa, ocasio-
nando severos daños, y nuevamente inundaciones al sureste del paı́s, por las fuertes y
duraderas lluvias [OCDE, 2013].
Este mismo año se formó el huracán Ingrid. Fue el cuarto ciclón de la temporada e
impactó en el Océano Atlántico, al mismo tiempo que Manuel, cuando éste se encon-
5 RESULTADOS 197
Cabe resaltar que desde el año 1958, no se presentaban dos afectaciones climatológi-
cas en ambos océanos al mismo tiempo. Estas dos tormentas tropicales simultáneas tu-
vieron gran impacto en el sur del paı́s ocasionando grandes avenidas e inundaciones.
2014. En el año 2014 tuvo impacto el Huracán Odile. El ojo de este huracán tocó tierra el
dı́a 14 de septiembre a las 23:45 horas, a 10km del Cabo San Lucas, en Baja California
Sur. Alcanzó vientos máximos sostenidos de 205km/h, con rachas hasta de 250km/h.
Tuvo un desplazamiento hacia el noroeste a 28km/h. Este huracán llegó a estar en la
categorı́a 3. El dı́a 17 de septiembre, aproximadamente a las 11.30 horas, volvió a tocar
tierra en su trayectoria, bajando a tormenta tropical en la costa noroeste de Sonora a
75km al sureste de Puerto Peñasco con viento máximos sostenidos de 65km/h, alcanza-
do rachas hasta de 85km/h [WU, 2018]. La trayectoria de Odile se observa en la figura
5.159.
2015. En 2015 tocó tierra el Huracán Patricia a las 18:15 horas, el dı́a 23 de octubre, sien-
do de categorı́a 4. Tuvo vientos máximos sostenidos de 240 km/h y rachas de 295 km/h.
Se aproximó a las costas del paı́s en el estado de Tabasco, afectando directamente a los
municipios de Tenacatita, Cuestecomate y Navidad, y de manera indirecta a Cuixmala, el
Estrecho, La Manzanilla, Melaque, Huerta y Cihuatán. Este fenómeno ocasionó fuertes
lluvias en el centro y sur del paı́s, además de traer consigo bajas de presión y entradas
de humedad [WU, 2018]. La trayectoria de Patricia se observa en la figura 5.160.
5 RESULTADOS 198
Para llevar a cabo esto se consideraron los mismos porcentajes de la máxima extrac-
cin turbinable de 0.4 y 0.7 de acuerdo a trabajos previos del Instituto de Ingenierı́a de la
5 RESULTADOS 199
UNAM, y se variaron las pendientes de acuerdo con los rangos. Para cada pendiente, se
realizó una simulacin correspondiente, oscilando sus valores entre 1 y 2. A partir de las
expresiones
Con esta información, se calcularon los lı́mites máximos de extracciones por cada
una de las etapas y para cada una de las presas. Estos lı́mites determinaron las restric-
ciones consideradas en el modelo matemático de control estocástico con restricciones,
siendo el caso para restricciones constantes. Estos lı́mites se observan en la figura 5.163.
Para este trabajo se utilizarán los datos de la última curva guı́a dada por la CONA-
GUA, mismas que se muestra en la figura 5.164.
La figura 5.166 representa la curva guı́a, con los niveles máximos permitidos por cada
una de las quincenas, también para la Presa La Angostura.
Con los datos de la curva guı́a se realizaron interpolaciones por el método de La-
grange, a fin de aproximar las curvas guı́as a ecuaciones analı́ticas. Para ello se utilizó
el software libre R. A fin de tener una mejor aproximación a través de distintos ensayos,
se optó por seccionar la curva guı́a por etapas, para después trabajar con una función
general seccionada. Para las primeras 6 quincenas, que corresponden a los meses de
Enero, Febrero y Marzo, se presentan los códigos empleados en el cálculo; donde x es
el tiempo en quincenas y P (x) el volumen almacenado en hm3 .
5 RESULTADOS 203
line 1 require(PolynomF)
line 2 x = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24)
line 3 y = c(16157, 16157, 16157, 15943, 15350, 14775, 14200, 13625, 13050, 12475,
11900, 11900, 13992, 14153, 14297, 14482, 14789, 15077, 15407, 15826, 16157,
16157, 16157, 16157)
line 6 polyAjuste
line 8 curve(polyAjuste,add=T)
line 1 require(PolynomF)
line 2 x = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24)
line 3 y = c(16157, 16157, 16157, 15943, 15350, 14775, 14200, 13625, 13050, 12475,
11900, 11900, 13992, 14153, 14297, 14482, 14789, 15077, 15407, 15826, 16157,
16157, 16157, 16157)
line 6 polyAjuste
line 8 curve(polyAjuste,add=T)
line 2 x = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24)
line 3 y = c(9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9560,
9270, 8925, 8680, 8418, 8557, 8811, 9106, 9339, 9801, 9801, 9801, 9801)
line 6 polyAjuste
line 8 curve(polyAjuste,add=T)
5 RESULTADOS 207
line 1 require(PolynomF)
line 2 x = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24)
5 RESULTADOS 208
line 3 y = c(9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9560,
9270, 8925, 8680, 8418, 8557, 8811, 9106, 9339, 9801, 9801, 9801, 9801)
line 6 polyAjuste
line 8 curve(polyAjuste,add=T)
line 1 require(PolynomF)
line 2 x = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24)
line 3 y = c(9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9801, 9560,
9270, 8925, 8680, 8418, 8557, 8811, 9106, 9339, 9801, 9801, 9801, 9801)
line 6 polyAjuste
line 8 curve(polyAjuste,add=T)
Figura 5.174: Capacidad, NAME, NAMO y NAMINO de cada una de las presas
5 RESULTADOS 211
Dichos intervalos son espacios de Borel con su respectiva σ-álgebra de Borel, con lo
cual se cumple el primer requisito para formar parte del sexteto de los modelos de control
de Markov.
Con lo cual, los espacios de acciones para la Presa La Angostura son los intervalos
cerrados [0, 4769.28], [0, 2384.64], [0, 2384.64], [0, 2384.64], [0, 4769.28] y [0, 11923.20], para
las etapas 1, 2, 3, 4, 5 y 6, respectivamente. De igual manera, para la Presa Malpaso,
los espacios de acciones son los intervalos cerrados [0, 6220.80], [0, 3110.40], [0, 3110.40],
[0, 3110.40], [0, 6220.80] y [0, 15552.0], para las etapas 1, 2, 3, 4, 5 y 6, respectivamente;
todos en unidades de metros.
Por ser subconjuntos de los espacios de acciones anteriores, también son conjuntos
borelianos, mismos que trabajarán con la σ-álgebra de Borel.
Este kernel dependió de los ingresos que estuvieron relacionados con los registros
históricos dados por distribuciones de probabilidad, con lo cual, el kernel se estable-
ció en términos de la etapa en cuestión. Esta dependencia se indicó como una función
f (x1 , a), donde x1 es el estado y a la acción.
fk (x1 , a) = x1 + m − a. (5.23)
Para cada etapa y para cada presa, se plantearon un total de 12 funciones de transición,
a saber:
Y para la presa 2:
En las ecuaciones anteriores, se debe recordar que mi son variables aleatorias idénti-
camente distribuidas, más aún, todas presentan una distribución normal, donde para ca-
da etapa se tuvieron las siguientes medias, que se observan en la figura 5.179.
X − 1176.41 X − 1101.30
Z11 = , Z12 =
392.92 355.30
X − 1962.51 X − 1292.44
Z21 = , Z22 =
752.51 715.18
X − 2319.75 X − 1516.70
Z31 = , Z32 =
860.14 819.90
X − 1447.56 X − 810.05
Z41 = , Z42 =
738.05 432.26
X − 1955.92 X − 1290.15
Z51 = , Z52 =
762.23 644.13
X − 1107.45 X − 1332.31
Z61 = , Z62 =
252.59 444.10
Para Malpaso:
Para Malpaso:
En caso de existir un posible derrame, éste estuvo dado por d2 = x2 − 6567 para La
Angostura y x2 − 3400 para Malpaso. Y en caso de existir déficits, estuvieron dados por
d1 = |x2 | cuando x2 < 0. La extracción que pasa por las turbinas fue además
ET = a − d1 , (5.38)
con a la extracción y d1 el déficit. En caso de existir derrames, esto pasaron a ser parte
de las extracciones de la presa de Malpaso, y el volumen turbinado se añadió a estos
derrames, los cuales fueron parte de las extracciones.
El nivel x2 se comparó con los niveles de las curvas guı́a dadas por la CONAGUA, y
en caso de cumplir lo requisitos, se corrigieron los niveles del embalse agregando los ni-
veles iniciales de 500 y 144 m, respecitvamente para La Angostura y Malpaso, y quitando
los niveles entre el NAMINO y el nivel medio de desfogue, de 269.5 y 93.9m, respecti-
vamente para cada embalse. En caso de existir déficit o derrame, se penalizó con los
coeficientes dados por el Instituto de Ingenierı́a para cada presa.
El cálculo del beneficio esperado máximo se realizó integrando sobre todo el espacio
de estados [0, 6567] de la presa La Angostura los beneficios esperados para cada uno de
los estados, junto con su distribución de probabilidades
Z 6567
BEM = BEdx1 . (5.43)
0
5 RESULTADOS 219
Z 6567
BEM = [r(x1 , a2 ) − c(x1 , a2 )] ∗ P (m)dx1 , (5.44)
0
El beneficio esperado máximo dependió del estado actual del sistema, de la extracción
y del ingreso. En este caso, se pudo considerar como una función BEM (x1 , a2 , m), con
x1 , a2 y m como se han definido anteriormente.
De acuerdo con la CFE, los lı́mites superiores e inferiores que en cada etapa se deben
mantener se observan en la figura 5.182.
En este caso, x representa el tiempo medido en quincenas de acuerdo con las etapas
establecidas por el Instituto de Ingenierı́a. Ya que en principio, estas restricciones no
dependen del estado, sino del tiempo, las funciones se consideraron después del cálculo
respectivo de los beneficios, como una cota de comparación para cada estado que se
analizó.
5.5.7. Beneficios
Ya que se deseó conocer la extracción que presentara el mayor beneficio, la función
de beneficio esperado se derivó con respecto a la extracción m, manteniendo los demás
valores constantes, y se igualó a cero.
Si Z 6567
BEM = [r(x1 , a2 ) − c(x1 , a2 )] ∗ P (m)dx1 , (5.51)
0
donde
r(x1 , a2 ) = B ∗ 0.01, (5.52)
c(x1 , a2 ) = d1 ∗ coef 1 + d2 ∗ coef 2, (5.53)
y
−3 e1 + e2
B = 2.725 × 10 ∗ a2 ∗ , (5.54)
2
luego
d
BEM (x1 , a2 , m) = 0 (5.55)
da2
Z 6567
d
[r(x1 , a2 ) − c(x1 , a2 )] ∗ P (m)dx1 = 0 (5.56)
da2 0
En este caso, para el ingreso m ahora se consideró el volumen turbinado por la presa
2, es decir, se contaron ingresos por cuenca propia más los ingresos provenientes aguas
arriba del cauce. Con lo cual
m = m2 + a2 (5.59)
5 RESULTADOS 221
donde BET es el beneficio total, BEM1 y BEM2 son los beneficios esperados de la pre-
sa 1 y la presa 2, respectivamente.
A partir de los datos anteriores, las polı́ticas de operación representaron una función
de dos variables (los estados de cada presa) y como resultado nuevamente es una fun-
ción de dos variables (la extracción de cada presa). En este sentido
donde x1 , x2 pertenecen a los respectivos espacios de estados de cada una de las presas
(volúmenes útiles) y a1 , a2 pertenecen a los espacios de acciones admisibles. Matemáti-
camente, la gráfica de estas ecuaciones se representa en un espacio de 4 dimensiones,
las cuales se graficaron únicamente por medio de sus trazas para algunos de los estados
representativos.
Figura 5.219: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 600 hm3
Figura 5.220: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 1200 hm3
5 RESULTADOS 242
Figura 5.221: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 1800 hm3
Figura 5.222: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 2400 hm3
5 RESULTADOS 243
Figura 5.223: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 3000 hm3
Figura 5.224: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 3600 hm3
5 RESULTADOS 244
Figura 5.225: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 4200 hm3
Figura 5.226: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 4800 hm3
5 RESULTADOS 245
Figura 5.227: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 5400 hm3
Figura 5.228: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 6000 hm3
5 RESULTADOS 246
Figura 5.229: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 6600 hm3
Figura 5.230: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 7200 hm3
5 RESULTADOS 247
Figura 5.231: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 7800 hm3
Figura 5.232: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 8400 hm3
5 RESULTADOS 248
Figura 5.233: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 9000 hm3
Figura 5.234: Beneficio esperado máximo para la etapa 1, para una estado de Malpaso
de 9600 hm3
5 RESULTADOS 249
Figura 5.235: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 600, 1200 y 1800 hm3
Figura 5.236: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 2400 hm3
5 RESULTADOS 250
Figura 5.237: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 3000 hm3
Figura 5.238: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 3600 hm3
5 RESULTADOS 251
Figura 5.239: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 4200 hm3
Figura 5.240: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 4800 hm3
5 RESULTADOS 252
Figura 5.241: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 5400 hm3
Figura 5.242: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 6000 hm3
5 RESULTADOS 253
Figura 5.243: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 6600 hm3
Figura 5.244: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 7200 hm3
5 RESULTADOS 254
Figura 5.245: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 7800, 8400 y 9000 hm3
Figura 5.246: Beneficio esperado máximo para la etapa 2, para una estado de Malpaso
de 9600 hm3
5 RESULTADOS 255
Figura 5.247: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 600 hm3
Figura 5.248: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 1200 hm3
5 RESULTADOS 256
Figura 5.249: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 1800 hm3
Figura 5.250: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 2400 hm3
5 RESULTADOS 257
Figura 5.251: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 3000 hm3
Figura 5.252: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 3600 hm3
5 RESULTADOS 258
Figura 5.253: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 4200 hm3
Figura 5.254: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 4800 hm3
5 RESULTADOS 259
Figura 5.255: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 5400 hm3
Figura 5.256: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 6000 hm3
5 RESULTADOS 260
Figura 5.257: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 6600 hm3
Figura 5.258: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 7200 y 7800 hm3
5 RESULTADOS 261
Figura 5.259: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 8400 hm3
Figura 5.260: Beneficio esperado máximo para la etapa 3, para una estado de Malpaso
de 9000 y 9600 hm3
5 RESULTADOS 262
Figura 5.261: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 600 hm3
Figura 5.262: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 1200 hm3
5 RESULTADOS 263
Figura 5.263: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 2400 hm3
Figura 5.264: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 3600 hm3
5 RESULTADOS 264
Figura 5.265: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 4800 hm3
Figura 5.266: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 6000 hm3
5 RESULTADOS 265
Figura 5.267: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 7200 hm3
Figura 5.268: Beneficio esperado máximo para la etapa 4, para una estado de Malpaso
de 8400 hm3
5 RESULTADOS 266
Figura 5.269: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 1200 hm3
Figura 5.270: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 2400 hm3
5 RESULTADOS 267
Figura 5.271: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 3600 hm3
Figura 5.272: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 4800 hm3
5 RESULTADOS 268
Figura 5.273: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 6000 hm3
Figura 5.274: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 7200 hm3
5 RESULTADOS 269
Figura 5.275: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 8400 hm3
Figura 5.276: Beneficio esperado máximo para la etapa 5, para una estado de Malpaso
de 9600 hm3
5 RESULTADOS 270
Figura 5.277: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 600 hm3
Figura 5.278: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 2400 hm3
5 RESULTADOS 271
Figura 5.279: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 3600 hm3
Figura 5.280: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 4800 hm3
5 RESULTADOS 272
Figura 5.281: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 6000 hm3
Figura 5.282: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 7200 hm3
5 RESULTADOS 273
Figura 5.283: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 8400 hm3
Figura 5.284: Beneficio esperado máximo para la etapa 6, para una estado de Malpaso
de 9600 hm3
5 RESULTADOS 274
Para cada una de las etapas, los estados de Malpaso recorren desde los 600 hm3
en adelante. Para cada estado de Malpaso, se tienen dos grupos de funciones: para La
Angostura y para Malpaso, que representan las extracciones a realizar, donde la variable
x representa el estado en La Angostura, el cual va desde los 600 hm3 hasta el volumen
máximo de 13,200 hm3 .
Etapa 1.
Estado de Malpaso: 600 hm3 .
Para La Angostura.
600,
si 600 ≤ x ≤ 2400
x − 1800, si 2400 ≤ x ≤ 3600
1800, si 3600 ≤ x ≤ 6600
7
x − 13600, si 6600 ≤ x ≤ 7200
3
3600, si 7200 ≤ x ≤ 7800
Extracción = . (5.64)
34800 − 4x, si 7800 ≤ x ≤ 8400
x − 7200, si 8400 ≤ x ≤ 9000
1800, si 9000 ≤ x ≤ 10800
3x − 30600, si 10800 ≤ x ≤ 11400
3600, si 11400 ≤ x ≤ 13200
Para Malpaso.
600, si 600 ≤ x ≤ 2400
x − 1800, si 2400 ≤ x ≤ 3600
1800, si 3600 ≤ x ≤ 6000
3x − 16200, si 6000 ≤ x ≤ 6600
3600, si 6600 ≤ x ≤ 8400.
Extracción = (5.66)
37200 − 4x, si 8400 ≤ x ≤ 9000.
x − 7800, si 9000 ≤ x ≤ 9600.
1800, si 9600 ≤ x ≤ 10200.
3x − 28800, si 10200 ≤ x ≤ 10800.
3600, si 10800 ≤ x ≤ 13200
Para Malpaso.
600, si 600 ≤ x ≤ 2400
x − 1800, si 2400 ≤ x ≤ 3600
si 3600 ≤ x ≤ 6000
1800,
3x − 16200, si 6000 ≤ x ≤ 6600
Extracción = 3600, si 6600 ≤ x ≤ 9000 (5.70)
30600 − 3x, si 9000 ≤ x ≤ 9600
1800, si 9600 ≤ x ≤ 10200
3x − 28800, si 10200 ≤ x ≤ 10800
si 10800 ≤ x ≤ 13200
3600,
Para Malpaso.
n
Extracción = 1200, si 600 ≤ x ≤ 13200 (5.71)
Para Malpaso.
(
1200, si 600 ≤ x ≤ 12600
Extracción = (5.73)
x − 11400, si 12600 ≤ x ≤ 13200
Estado de Malpaso: 3600 hm3 .
Para La Angostura.
5 RESULTADOS 277
600, si 600 ≤ x ≤ 2400
x − 1800, si 2400 ≤ x ≤ 3600
si 3600 ≤ x ≤ 5400
1800,
3x − 14400, si 5400 ≤ x ≤ 6000
Extracción = 3600, si 6000 ≤ x ≤ 8400 (5.74)
28800 − 3x, si 8400 ≤ x ≤ 9000
1800, si 9000 ≤ x ≤ 9600
3x − 27000, si 9600 ≤ x ≤ 10200
si 10200 ≤ x ≤ 13200
3600,
Para Malpaso.
(
1200, si 600 ≤ x ≤ 12000
Extracción = (5.75)
x − 10800, si 1200 ≤ x ≤ 13200
Etapa 2.
Estado de Malpaso: 600 hm3 .
Para La Angostura.
600, si 600 ≤ x ≤ 3600
2x − 2200, si 3600 ≤ x ≤ 4200
Extracción = 1800, si 4200 ≤ x ≤ 11400 (5.76)
x − 9600, si 11400 ≤ x ≤ 12000
2400, si 12000 ≤ x ≤ 13200
Para Malpaso.
Para Malpaso.
(
1200, si 600 ≤ x ≤ 12000
Extracción = (5.83)
1/19x − 21600/19, si 12000 ≤ x ≤ 13200
Estado de Malpaso: 3600 hm3 .
Para La Angostura.
600, si 600 ≤ x ≤ 3600
2x − 6600, si 3600 ≤ x ≤ 4200
Extracción = 1800, si 4200 ≤ x ≤ 11400 (5.84)
x − 9600, si 11400 ≤ x ≤ 12000
2400, si 12000 ≤ x ≤ 13200
Para Malpaso.
(
1200, si 600 ≤ x ≤ 11400
Extracción = (5.85)
x − 10200, si 11400 ≤ x ≤ 13200
Estado de Malpaso: 4200 hm3 .
Para La Angostura.
5 RESULTADOS 279
600, si 600 ≤ x ≤ 3600
2x − 6600, si 3600 ≤ x ≤ 4200
1800, si 4200 ≤ x ≤ 6600
15000 − 2x, si 6600 ≤ x ≤ 7200
Extracción = (5.86)
2x − 13800, si 7200 ≤ x ≤ 7800
si 7800 ≤ x ≤ 11400
1800,
x − 9600, si 11400 ≤ x ≤ 12000
2400, si 12000 ≤ x ≤ 13200
Para Malpaso.
1200, si 600 ≤ x ≤ 10800
2x − 20400, si 10800 ≤ x ≤ 11400
Extracción = 11400, si 11400 ≤ x ≤ 12000 (5.87)
x − 9600, si 12000 ≤ x ≤ 12600
3000, si 12600 ≤ x ≤ 13200
Etapa 3.
Estado de Malpaso: 600 hm3 .
Para La Angostura.
600, si 600 ≤ x ≤ 2400
2x − 4200, si 2400 ≤ x ≤ 3000
Extracción = (5.88)
1800, si 3000 ≤ x ≤ 7800
si 7800 ≤ x ≤ 13200
2400,
Para Malpaso.
1200,
si 600 ≤ x ≤ 10200.
Extracción = x − 9000, si 10200 ≤ x ≤ 12000. (5.89)
3000, si 12000 ≤ x ≤ 13200.
Para Malpaso.
1200,
si 600 ≤ x ≤ 9600.
Extracción = x − 8400, si 9600 ≤ x ≤ 11400. (5.91)
3000, si 11400 ≤ x ≤ 13200.
Para Malpaso.
1200,
si 600 ≤ x ≤ 8400.
Extracción = x − 7200, si 8400 ≤ x ≤ 10200. (5.95)
3000, si 10200 ≤ x ≤ 13200.
600, si 600 ≤ x ≤ 6600.
2x − 12600, si 6600 ≤ x ≤ 7200.
Extracción = 1800, si 7200 ≤ x ≤ 7800. (5.96)
x − 6000, si 7800 ≤ x ≤ 8400.
2400, si 8400 ≤ x ≤ 13200.
Para Malpaso.
1200, si 600 ≤ x ≤ 7800.
x − 6600, si 7800 ≤ x ≤ 8400.
Extracción = 1800, si 8400 ≤ x ≤ 9000. (5.97)
2x − 16200, si 9000 ≤ x ≤ 9600.
3000, si 9600 ≤ x ≤ 13200.
Estado de Malpaso: 3600 hm3 .
Para La Angostura.
600, si 600 ≤ x ≤ 6600.
2x − 12600, si 6600 ≤ x ≤ 7200.
Extracción = 1800, si 7200 ≤ x ≤ 7800. (5.98)
x − 6000, si 7800 ≤ x ≤ 8400.
2400, si 8400 ≤ x ≤ 13200.
Para Malpaso.
1200,
si 600 ≤ x ≤ 7800.
Extracción = x − 6600, si 7800 ≤ x ≤ 8400. (5.99)
3000, si 8400 ≤ x ≤ 13200.
Etapa 4.
Estado de Malpaso: 600 hm3 .
Para La Angostura.
600,
si 600 ≤ x ≤ 4800
2x − 9000, si 4800 ≤ x ≤ 5400.
Extracción = (5.100)
1800, si 5400 ≤ x ≤ 6600
x − 4800, si 6600 ≤ x ≤ 13200
Para Malpaso.
5 RESULTADOS 282
1200, si 600 ≤ x ≤ 6600.
x − 5400, si 6600 ≤ x ≤ 7200.
1800, si 7200 ≤ x ≤ 9600.
Extracción = (5.101)
x − 7800, si 9600 ≤ x ≤ 10200.
si 10200 ≤ x ≤ 10800.
2400,
x − 8400, si 10800 ≤ x ≤ 13200.
600, si 600 ≤ x ≤ 6000.
x − 5400, si 6000 ≤ x ≤ 6600.
si 6600 ≤ x ≤ 8400.
1200,
x − 7200, si 8400 ≤ x ≤ 9000.
Extracción = 1800, si 9000 ≤ x ≤ 10200. (5.110)
12000 − x, si 10200 ≤ x ≤ 10800.
1200, si 10800 ≤ x ≤ 12000.
13200 − x, si 12000 ≤ x ≤ 12600.
si 12600 ≤ x ≤ 13200.
600,
Para Malpaso.
n
Extracción = 3000, si 600 ≤ x ≤ 13200. (5.111)
Etapa 5.
Estado de Malpaso: 1200 hm3 .
Para La Angostura.
1800 − x, si 600 ≤ x ≤ 1200.
x − 600, si 1200 ≤ x ≤ 1800.
si 1800 ≤ x ≤ 2400.
1200,
x − 1200, si 2400 ≤ x ≤ 3000.
Extracción = 1800, si 3000 ≤ x ≤ 4800. (5.112)
3x − 12600, si 4800 ≤ x ≤ 5400.
3600, si 5400 ≤ x ≤ 7200.
x − 3600, si 7200 ≤ x ≤ 7800.
si 7800 ≤ x ≤ 13200.
4200,
Para Malpaso.
1200,
si 600 ≤ x ≤ 3600.
x − 2400, si 3600 ≤ x ≤ 4800.
2400, si 4800 ≤ x ≤ 5400.
x − 3000, si 5400 ≤ x ≤ 6600.
2x − 9600, si 6600 ≤ x ≤ 7200.
Extracción = (5.113)
4800, si 7200 ≤ x ≤ 7800.
x − 3000, si 7800 ≤ x ≤ 8400.
5400, si 8400 ≤ x ≤ 9600.
x − 4200, si 9600 ≤ x ≤ 10200.
6000, si 10200 ≤ x ≤ 13200.
5 RESULTADOS 285
Para Malpaso.
5 RESULTADOS 286
1200, si 600 ≤ x ≤ 3600.
x − 2400, si 3600 ≤ x ≤ 4800.
si 4800 ≤ x ≤ 5400.
2400,
x − 3000, si 5400 ≤ x ≤ 7200.
Extracción = 4800, si 7200 ≤ x ≤ 7800. (5.117)
x − 3000, si 7800 ≤ x ≤ 8400.
5400, si 8400 ≤ x ≤ 9600.
x − 4200, si 9600 ≤ x ≤ 10200.
si 10200 ≤ x ≤ 13200.
6000,
Para Malpaso.
3000, si 600 ≤ x ≤ 2400.
5400 − x, si 2400 ≤ x ≤ 3000.
4x − 9600, si 3000 ≤ x ≤ 3600.
Extracción = (5.121)
4800, si 3600 ≤ x ≤ 4800.
si 4800 ≤ x ≤ 6000.
x,
6000, si 6000 ≤ x ≤ 13200.
Etapa 6.
Estado de Malpaso: 600 hm3 .
Para La Angostura.
5 RESULTADOS 288
600, si 600 ≤ x ≤ 3000.
2x − 3600, si 3000 ≤ x ≤ 3600.
x, si 3600 ≤ x ≤ 4200.
4200, si 4200 ≤ x ≤ 5400.
Extracción = (5.124)
x − 1200 si 5400 ≤ x ≤ 6000.
si 6000 ≤ x ≤ 7200.
4800
6x − 38400 si 7200 ≤ x ≤ 7800.
8400 si 7800 ≤ x ≤ 13200.
Para Malpaso.
3000 si 600 ≤ x ≤ 7200.
2x − 11400 si 7200 ≤ x ≤ 7800.
Extracción = 4200 si 7800 ≤ x ≤ 9000. (5.125)
x − 4800 si 9000 ≤ x ≤ 12600.
7800 si 7200 ≤ x ≤ 13200.
Estado de Malpaso: 2400 hm3 .
Para La Angostura.
2400 si 600 ≤ x ≤ 3600.
x − 600 si 3600 ≤ x ≤ 5400.
Extracción = 4800 si 5400 ≤ x ≤ 7200. (5.126)
6x − 38400 si 7200 ≤ x ≤ 7800.
8400 si 7800 ≤ x ≤ 13200.
Para Malpaso.
3000 si 600 ≤ x ≤ 6600.
x − 3600 si 6600 ≤ x ≤ 7200.
4x − 25300 si 7200 ≤ x ≤ 7800.
Extracción = (5.127)
6000 si 7800 ≤ x ≤ 9000.
x − 3000 si 9000 ≤ x ≤ 12600.
9600 si 12600 ≤ x ≤ 13200.
2400 si 600 ≤ x ≤ 3000.
x − 600 si 3000 ≤ x ≤ 5400.
Extracción = 4800 si 5400 ≤ x ≤ 7800. (5.128)
6x − 42000 si 7800 ≤ x ≤ 8400.
8400 si 8400 ≤ x ≤ 13200.
Para Malpaso.
3000 si 600 ≤ x ≤ 4800.
x − 1800 si 4800 ≤ x ≤ 6000.
4200 si 6000 ≤ x ≤ 6600.
Extracción = 4x − 24000 si 6600 ≤ x ≤ 8400. (5.129)
7200 si 8400 ≤ x ≤ 9600.
x − 2400
si 9600 ≤ x ≤ 10800.
5x − 45600 si 10800 ≤ x ≤ 13200.
Estado de Malpaso: 4800 hm3 .
Para La Angostura.
2400 si 600 ≤ x ≤ 3000.
x − 600 si 3000 ≤ x ≤ 5400.
Extracción = 4800 si 5400 ≤ x ≤ 7800. (5.130)
6x − 42000 si 7800 ≤ x ≤ 8400.
8400 si 8400 ≤ x ≤ 13200.
Para Malpaso.
x + 3000 si 600 ≤ x ≤ 1200.
6600 − 2x si 1200 ≤ x ≤ 1800.
3000 si 1800 ≤ x ≤ 3600.
x − 600 si 3600 ≤ x ≤ 5400.
5400 si 5400 ≤ x ≤ 6600.
Extracción = (5.131)
x − 1200 si 6600 ≤ x ≤ 7800.
3x − 16800 si 7800 ≤ x ≤ 8400.
8400 si 8400 ≤ x ≤ 9600.
5x − 39600 si 9600 ≤ x ≤ 10200.
11400 si 10200 ≤ x ≤ 13200.
2400 si 600 ≤ x ≤ 3000.
x − 600 si 3000 ≤ x ≤ 5400.
Extracción = 4800 si 5400 ≤ x ≤ 7200. (5.132)
6x − 38400 si 7200 ≤ x ≤ 7800.
8400 si 7200 ≤ x ≤ 13200.
Para Malpaso.
x + 2400 si 600 ≤ x ≤ 1200.
3600 si 1200 ≤ x ≤ 3000.
x + 600 si 3000 ≤ x ≤ 5400.
6000 si 5400 ≤ x ≤ 6000.
x si 6000 ≤ x ≤ 7200.
Extracción = (5.133)
14400 − x si 7200 ≤ x ≤ 7800.
8x − 55800 si 7800 ≤ x ≤ 8400.
11400 si 8400 ≤ x ≤ 11400.
x si 11400 ≤ x ≤ 12600.
12600 si 12600 ≤ x ≤ 13200.
Las ecuaciones anteriores conforman cada una de las polı́ticas de extracción depen-
diendo de la etapa y dependiendo del nivel de cada una de las presas. La manera de
interpretar estas ecuaciones se explica a continuación. Cada bloque de ecuaciones está
agrupado por etapas de la 1 a la 6. Por cada etapa, se tienen distintos estados de la
presa Malpaso que van desde los 600 hm3 hasta los 9600 hm3 . Para cada uno de estos
estados de Malpaso, se cuentan con dos ecuaciones seccionadas, una para La Angos-
tura y otra para Malpaso; estas ecuaciones son las polı́ticas óptimas, donde la variable
x representa el estado ahora de la Presa La Angostura, y dependiendo en dónde se en-
cuentre, se sustituirá en la ecuación adecuada. Ası́, por ejemplo, si se está en el mes de
noviembre, corresponden a la etapa 6. Si además la capacidad de Malpaso es de 600
hm3 y la capacidad de La Angostura es de 7500 hm3 , basta con sustituir x = 7500, en las
ecuaciones respectivas para La Angostura y Malpaso, y se obtienen las extracciones de
cada una de las presas en millones de m3 .
Para la Angostura:
6(7500) − 38400 = 6600. (5.134)
Para Malpaso:
2(7500) − 11400 = 3600. (5.135)
De esta manera, el operador puede calcular en forma sencilla el volumen de extrac-
ción para cada etapa del año y para cada estado de las presas, identificando la ecuación
y sustituyendo los respectivos estados.
5 RESULTADOS 291
Para cada año de registro se calculó a partir de las polı́ticas óptimas, los volúmenes
de extracción, y con ellos se calcularon con las curvas elevación-capacidad, las respec-
tivas alturas del embalse. Posteriormente se hizo un promedio de todos los años y esto
5 RESULTADOS 292
Figura 5.289: Energı́a simulada promedio de la presa La Angostura con respecto al tiem-
po
6. Conclusiones y recomendaciones
6.1. Conclusiones
El funcionamiento de una presa hidroeléctrica depende de su operación, que a su vez
está en función de los volúmenes a extraer con base en la etapa del año y el estado
de llenado en que se encuentre el embalse. Para ello se necesitan aplicar técnicas que
permitan alcanzar los objetivos propuestos dadas las restricciones planteadas.
Los modelos de control markovianos han demostrado alta precisión para ajustarse a
la realidad de diferentes sistemas dinámicos y fenómenos que ocurren, ya que precisa-
mente involucran los estados, las etapas, los objetivos y las restricciones de un problema.
La aplicación del correcto modelo y el planteamiento de las funciones adecuadas permi-
tirá que los resultados sean los adecuados.
El Sistema Hidroeléctrico Grijalva es uno de los más importantes del paı́s al contribuir
con la aportación de energı́a eléctrica en un porcentaje notable; su ubicación geográfica,
al igual que sus condiciones topográficas, ha permitido aprovechar al máximo el caudal
de rı́o Grijalva en cuanto a generación. Sin embargo, la misma topogafı́a conduce a la
implementación de un sistema en cascada para aprovechamiento eléctrico cuya com-
plejidad de su organización en serie, radica en el hecho de trabajar con escurrimientos
propios y escurrimientos provenientes del sistema aguas arriba, los cuales son las entra-
das del modelo de optimización. Esta complejidad ha podido ser disminuida simplificado
el número de embalses debido a la aportación y volumen de almacenamiento de manera
individual. Además, contar con un monitoreo constante de los niveles de cada embalse,
permite a los operarios conocer la información de llenado o descarga del sistema, para
generación de energı́a eléctrica, a fin de maximiar la producción y a la vez minimizar
posibles pérdidas por derrames o déficits. Este factor es determinante, sobre todo cono-
ciendo la etapa en que se encuentre el modelo (sea temporada de lluvias o de estiaje), y
tener las proyecciones futuras en cuanto a las descargas realizadas.
6 CONCLUSIONES Y RECOMENDACIONES 295
Los distintos estudios realizados han demostrado las graves repercusiones que ocu-
rren cuando se presenta una avenida máxima, y los daños ocasionados en poblados
aguas abajo de los cauces o en este caso, de los embalses. Los impactos van desde los
niveles económicos, los niveles sociales, hasta niveles de tipo sanitario. Es por ello que
se deben minimizar las posibilidades de desbordes en los cauces, manteniendo niveles
mı́nimos en los sistemas hidroeléctricos. Los derrames que se lleguen a producir a través
de las obras de excedencia serán responsabilidad de la operación del sistema. Es nece-
sario entonces, mantener un balance entre los niveles máximos energéticos y los niveles
mı́nimos de seguridad.
Más aún, el sistema debe satisfacer una energı́a mı́nima de base en todo el año, y
además garantizar que se cumplan las demandas por energı́as pico. En resumen, son
4 las condiciones que un sistema hidroeléctrico debe cumplir, y todas son consecuencia
de los niveles del embalse, que a su vez dependen directamente de las extracciones que
realice el operador; cumplir demandas máximas, demandas mı́nimas, la demanda base
y no poner en riesgo a la población.
Para cumplir lo anterior, se planteó el estudio del sistema mediante la modelación es-
tocástica. Las diferentes técnicas probabilı́sticas han evolucionado con el paso del tiem-
po, recreando condiciones cada vez más reales de los funcionamientos de un sistema.
La complejidad cada vez más es mayor, pero los resultados obtenidos motivan a estudiar
estos modelos con mayor detalle. La naturaleza aleatoria de los escurrimientos lleva a
pensar, que los modelos estocásticos pueden simular de una mejor manera el funciona-
miento de un embalse.
las variables que se estudian en un sistema hidroeléctrico; con él, se han modelado los
niveles de los embalses, como espacio de estados; las extracciones como espacio de
acciones; las extracciones permitidas como espacio de acciones admisibles; la ecuación
de continuidad como el kernel de transición; la energı́a generada como la función de
ganancia; los penalizaciones por derrames y déficits como la función de costo; y las deli-
mitaciones de la CFE y la CONAGUA como la función de restricción. Ya que los espacios
de estados y acciones, forman conjuntos compactos, existen para la función beneficio, un
máximo posible de alcanzar, sujeto a las restricciones impuestas. Las condiciones ma-
temáticas de las funciones anteriores ha garantizado que exista dentro de las polı́ticas
aleatorizadas un conjunto de polı́ticas deterministas, que puedan optimizar a la función
objetivo.
El presente trabajo ha logrado implementar las ténicas actuales del área de control
estocástico a la simulación de las leyes de transición que modelan un embalse, incor-
porando los conceptos de modelos markovianos, espacios borelianos, continuidad de
funciones, entre otros.
Además, el poder trabajar con funciones continuas sobre espacios de estados com-
pactos, garantizó la existencia de esta clase de polı́ticas en el conjunto dado. Los fun-
damentos matemáticos usados lograron que la función objetivo pudiera ser resuelta para
cada una de las extracciones, ası́ como las restricciones impuestas delimitaron la factibi-
lidad del problema.
Por otro lado, los escurrimientos históricos presentaron un comportamiento con dis-
tribución normal, entre otras posibles distribuciones de probabilidad. El contar con estas
distribuciones permitió que se cuenten con tablas de probabilidades para cualquier es-
currimiento, no únicamente para los casos discretos sino para todos los posibles escurri-
mientos que se presenten dentro de los intervalos de espacios admisibles.
Ası́ pues, los beneficios se calcularon como soluciones de integrales definidas con
respecto a la medida de probabilidad, de las distribuciones de los escurrimientos. Y los
óptimos fueron alcanzados por los procesos de diferenciación sobre esta clase de espa-
cios.
El haber trabajado con estados continuos permitió que las polı́ticas óptimas también
se expresaran en forma continua. Para las primeras etapas mateniendo los estados de
6 CONCLUSIONES Y RECOMENDACIONES 297
Malpaso fijos, se observó que las extracciones de Malpaso son mayormente constan-
tes, siendo las extracciones de La Angostura las que mayor variabilidad presentan; esto
cuando los niveles de Malpaso son inferiores a los 4,200 hm3 . A partir de este nivel,
las polı́ticas de Malpaso se asemejaron a las polı́ticas de La Angostura hasta niveles de
7,800 hm3 . A partir de estos volúmenes, las polı́ticas formaron curvas bastante simétri-
cas, es decir, salvo una constante, las polı́ticas de extracción son las mismas.
6.2. Recomendaciones
Pese a que el modelo propuesto en este trabajo realizó la caracterización del sistema
para estados continuos, se recomienda realizar simulaciones para estados discretos a fin
de comprender mejor la forma en la que opera el modelo.
Por otro lado, las curvas elevación-capacidad dependen de los estudios batimétricos
de las presas, los cuales no han sido actualizados desde su construcción. Se recomienda
un nuevo estudio a fin de generar nuevas ecuaciones.
Finalmente se recomienda una actualización de las curvas guı́a, sobre todo en la
época de lluvias debido a las afectaciones de cambio climático y los eventos particulares
como el Niño o la Niña, haciendo un estudio de las afectaciones a corto, mediano y largo
plazo.
discretizado por 6 etapas, sino serı́a para todos los dı́as del año.
Es posible también, dadas las condiciones climatológicas de los últimos 5 años, los
cambios climáticos y las recientes sequı́as en la zona, que ocasiona niveles del embalse
bajos, un nuevo estudio sobre las curvas guı́as para incorporarlas al modelo. Adicional-
mente en el funcionamiento del vaso es importante incluir la correlación existente entre
los volúmenes del ingreso de una quicena a otra, o de un mes a otro, mediante com-
paraciones con los volúmenes medios de ingreso, para incluir la posibilidad de hacer
adecuaciones a las polı́ticas. En este mismo sentido, se deben actualizar periódicamente
los datos históricos obteniendo nuevas funciones de distribución de probabilidades del in-
greso, tomando en cuenta el comportamiento estadı́stico ante nuevos eventos extremos.
Por otro lado, se pueden implementar diferentes métodos de optimización para rea-
lizar un contraste entre las variables discretas y las variables continuas; estos métodos
pueden ir desde anális teóricos y variantes de la Programación Dinámica Estocásticas y
el Control Estocástico, hasta los métodos robustos y computación evolutiva desarrollado
en los últimos años.
Finalmente las funciones de restricción, en este caso, los lı́mites máximos y mı́nimos
dados por la CFE, se compararon por fuera de la función objetivo, de hecho como cons-
tantes, más que como funciones. Es posible incorporarlas al modelo, mediante funciones
restrictivas y estudiar su influencia en el cálculo de las polı́ticas.
REFERENCIAS 299
Referencias
[Alegrı́a, 2010] A LEGR ÍA D., A., Polı́tica de operación óptima del sistema de presas del
rı́o Grijalva. Efectos de la curva guı́a., Tesis de Maestrı́a, México, UNAM, 2010.
[Altman, 1999] A LTMAN , E., Constrained Markov Decision Processes, Chapman
Hall/CRC, Boca Raton, FL., 1999.
[Arellano, 1997] A RELLANO M ONTERROSAS , J. L., Manejo de cuencas hidrológicas,
CNA, Simposio 4, IX Congreso Nacional de Irrigación, 1999.
[Arganis, 2004] A RGANIS J., M. L., Operación óptima de un sistema de presas en cas-
cada para generación hidroeléctrica tomando en cuenta condiciones reales de ope-
ración y el uso de muestras sintéticas para el pronóstico. Tesis Doctoral, UNAM,
2004.
[Arzate, 2000] A RZATE , R. S., Localizado en la cuenca del Grijalva el
mayor número de centrales hidroeléctricas, El Nacional en Internet.
www.unam.mx/servhem/nacional/1997/ene97//29ene97//29ec292.html.
[Ash, 2008] A SH , R OBERT B., Basic probability theory. Dover Publications, New York,
2008.
[Asmadi, Ahmed, 2014] A SMADI A HMAD, A HMED E L -S HAFI , S ITI FATIN, Reservoir Opti-
mization in Water Resources: A Review. Water Resources Managemente, Springer,
2014, Vol. 28, No. 11.
[Avilés, 1994] AVIL ÉS H ERRERA , R OBERTO, Optimización en lı́nea de Presas Hidro-
eléctricas. Tesis de Licenciatura. UNAM. 1994
[Bartle, 1964] B ARTLE , R. M., The elements of real analysis, Wiley International, Estados
Unidos, 1964.
[Basu, 2003] B ASU, A. K., Introduction to stochastic processes. Alpha Science Interna-
tional, Ltd, 2003.
[Bellman, 1957] B ELLMAN , R. E., Dynamic Programming. Princeton, NJ, 1957.
[Berezowsky, Domı́nguez, 1983] B EREZOWSKY V., M OIS ÉS , D OM ÍNGUEZ M., R AM ÓN ,
F UENTES . M., Ó SCAR. Manual de diseño de obras civiles. C.F.E. Hidrotecnia. A.2.8.
Planeación de sistemas de aprovechamiento hidroeléctrico. México, 1983.
[Bertsekas, 1995] B ERTSEKAS , D. P., Dynamic Programming and Optimal Control, Athe-
na Scientific, Belmont, Massachusetts, Prentice Hall, 1995.
[Borkar, 1990] B ORKAR , V. S., G HOSH , M. K., Controlled diffusions with constraints, J.
Math. Anal. Appl. 152 (1990), 88-108.
[Borkar, 1994] B ORKAR , V. S., Ergodic control of Markov chains with constraints- the
general case, SIAM, J. Control Optim. 32 (1994), 176-186.
REFERENCIAS 300
[Calva, 1993] C ALVA , S ALAZAR A., Bondad de una polı́tica de operación de una presa
de almacenamiento, VIII Congreso Nacional de Hidráulica. Asociación Mexicana de
Hidráulica, 1993.
[Campbell, 1992] C AMPBELL , J. C., R ADKE , J., G LESS , J.T., An Application of Linear
Programming and Geographic Infomation Systems: Cropland Allocation in Antigua.
Environment and Planning A. Abril, 1992.
[CFE, 1976] C OMISI ÓN F EDERAL DE E LECTRICIDAD, México construye: Proyectos Hidro-
eléctricos Chicoasén y La Angostura, México, 1976.
[Chen, 2007] C HEN , R. C., F EINBERG , E. A., Nonradomized policies for constrained
Markov decision processes, Maths. Methods Oper. Res. 66 (2007), No. 1, 165-179.
[CNA, 1997] CNA, SEMARNAP, Presas de México, 1982-1994, Volumen II, México,
1997.
[Cohn, 1980] C OHN , D ONALD, Measure Theory, Birkhauser Advanced, Second Edition,
1984. USA.
[Cohon, Marks, 1975] 17. C OHON , J. L., M ARKS , D. H., A review and Evaluation of Mul-
tiobjective Programming Tecniques. Water Resources Research 11(2), pp. 208-220,
1975.
[Comisión del Grijalva, 1964] C OMISI ÓN DEL G RIJALVA, Presa Netzahualcóyotl, Chiapas,
Secretarı́a de Recursos Hidráulicos, 1964.
REFERENCIAS 301
[CONAGUA, 1994] C OMISI ÓN N ACIONAL DEL AGUA, Presas de México, Vols. Secretarı́a
del Medio Ambiente y Recursos Naturales, 1994.
[Conant, 1948] C ONANT L, J. U., Use of storage water in hydroelectric system. Tesis
Doctoral. M.I.T., 1948.
[Costa, 2012] C OSTA , O. L. V., D UFOUR , F., Average control of Markov decision proces-
ses with Feller transition probabilities and general action spaces, J. Math. Anal. Appl.
396 (2012) 5869.
[Crichigno, Talavera, 2011] C RICHIGNO, J., TALAVERA , F., P RIETO, J., Enrutamiento
Multicast Utilizando Optimizacin Multiobjetivo. Ciencias y Tecnologa. Universidad
Catlica Nuestra Seora de Asuncin. Centro Nacional de Computacin. Universidad Na-
cional de Asuncin. 2012.
[Domı́nguez, 1998] D OM ÍNGUEZ , M. R., M ENDOZA R. R., Operación Integral del Siste-
ma Hidroeléctrico del Rı́o Grijalva, elaborado por la CFE por el Instituto de Ingenierı́a,
UNAM, México, 1998.
[Dufour, Stockbridge, 2010] D UFOUR , F., S TOCKBRIDGE , R. H., Existence of strict opti-
mal controls for discounted stochastic control problems, in Modern Trends in Con-
trolled Stochastic Processes: Theory and Applications, edite by A. B. Piunovskiy,
Luniver Press, Frome, U. ;., 2010, 12-22.
[Dutta, 1991] D UTTA , P. K., What do discounted optima converge to? A theory of discount
rate asymptotics in economic models, J. Econom. Theory 55 (1991) 6494.
[Feinberg, Shwartz, 1996] F EINBERG , E. A., S HWARTZ , A., Constrained discounted dy-
namic programming, Math. Oper. Res. 21 (1996), 922-945.
[Goicochea, Duckeistein, 1982] G OICOCHEA , A., D UCKEISTEIN , L., F OGEL , M. M., Mul-
tiobjective Decision Analysis with Engineering and Business Application. John Wiley
and Son, N. Y., USA, 1982.
[Gordienko, Hernandez, 1995] G ORDIENKO, E., H ERNANDEZ -L ERMA , O., Average cost
Markov control processes with weigthed norms: existence of canonical policies, Appl.
Math (Warsaw) 23 (1995), 199-218.
[Guo, Hernandez, 2009] G UO, X. P., H ERNANDEZ -L ERMA , O., Continuous-Time Markov
Decision Processes, Sprinver-Verlag, Berlin, Heidelberg, 2009.
[Guo, Huang, Song, 2012] G UO, X. P., H UANG , Y., S ONG , X., Linear programming and
constrained average optimality for general continuous-time Markov decision proces-
ses in history-depent policies, SIAM J. Control Optim. 50 (2012), 23-47.
[Guo, Quanxin, 2006] G UO X. P., Q UANXIN Z., Average optimality for Markov decision
processes in Borel spaces: a new condition and approach, J. Appl. Probab 43, 2006.
[Haviv, 1996] H AVIV M., On constrained Markov decision processes, Oper. Res. Lett. 10,
1996, pp. 25-28.
[Hernández, Laserre, 1991] H ERN ÁNDEZ -L ERMA , O, L ASERRE , J. B., 1991, Discrete-
Time Markov Control Processes. Basic Optimality Criteria, Springer, Scien-
te+Business Media New York, 1991.
REFERENCIAS 303
[Hernández, Laserre, 1991] H ERN ÁNDEZ -L ERMA , O, L ASERRE , J. B., 1991, Further To-
pics on Discrete-Time Markov Control Processes, Springer, Sciente+Business Media
New York, 1991.
[INE, 2005] I NSTITUTO N ACIONAL DE E COLOG ÍA, La cuenca de los rı́os Grijalva y Usu-
macinta, http://www.ine.gob.mx/uejaei/publicaciones/libros/402/cuencas.html, 2005.
[INEGI, 2007] I NSTITUTO N ACIONAL DE E STAD ÍSTICA Y G EOGRAF ÍA, Anuario Estadı́stico
de Chiapas, 2007, http://www.inegi.gob.mx.
[Karloff, 2009] K ARLOFF , H., Linear Programming, Modern Birkhausser Classic, Atlanta,
2009.
[Karterakis, 2007] K ARTERAKIS , S TEFANOS M., K ARATZAS G. P., N IKOLOS I. K., Appli-
cation of linear programming and differential evolutionary optimization methodologies
for the solution o coastal subsurface water management problems subject to envi-
ronmental criteria. Journal of Hydrology. Volumen 342. Issues 3-4, Septiembre 2007,
Pag. 270-280.
[Ko, Fontane, 1992] KO, S. K., F ONTANE , D. G., L ABADIE , J. W., Multiobjective Optimi-
zacion of Reservoir System Operation. Water Resources Bulletin. AWRA., Vol. 28,
No. 1, February, 1992.
[Labadie, 2000] L ABADIE , J. W., Reservoir System Optimization Models. Colorado State
University, Water Resources Update, University Council of Water Resources, 108.
Summer, 1997.
REFERENCIAS 304
[Loeve, 1977] P ROBABILITY T HEORY I, Loeve, M., Springer-Verlag, New York, 1977.
[Loganathan, 2002] L OGANATHAN , G.V., PARK , S., AND S HERALI , H.D., A Threshold
Break Rate for Pipeline Replacement in Water Distribution Systems, J. of Water Re-
sou. Plan. and Mgmt., ASCE, vol. 128, no. 4, July 2002, pp. 271279.
[Lund, Guzmán, 1999] L UND, J AY R., G UZM ÁN , J OEL, Derived operating rules for reser-
voirs in series or in parallel. Journal of waters resources planning and management.
Vol. 125, No. 3, Junio, 1999.
[Maidment, 1993] M AIDMENT, D. R., Handbook of Hydrology, Editor in Chief, USA, 2013.
[Marglin, 1967] M ARGLIN , S. A., Public Investment Criteria. The M.I.T. Press, Camgridge,
Massachusetts, USA. 1967.
[Marsden, 1974] M ARSDEN , J. E., Elementary Classical Analysis, Editorial W. II. Free-
man and Company, San Francisco, 1974.
[Mocholi, Sala, 1993] M OCHOLI , A., M., S ALA G., M., Programación Lineal. Metodologı́a
y problemas, Editorial Tebar Flores, Albacete, 1993.
[OCDE, 2013] OCDE, Estudio de la OCDE sobre el sistema nacional de Protección Civil
en México, México, 2013, p. 242.
[Peñuela, Granada, 2007] P EUELA M., C. A., G RANADA E., M., Optimizacin multiobjeti-
vo usando un algoritmo gentico y un operador elitista basado en un ordenamiento
no dominado (nsga-ii). Revista UTP. Universidad Tecnolgica de Pereira, Vol 1, Nmero
35, 2007.
[Rincón, 2012] R INC ÓN , L. Introducción a los procesos estocásticos. Facultad de Cien-
cias, UNAM. México, D.F., 2012.
[Rubio, Triana, 2006] RUBIO, G. H., T RIANA R. C., Gestión integrada de crecientes ca-
so de estudio México: Rı́o Grijalva, Programa Asociada de Gestión de Crecientes,
Organización Meteorológica Mundial (OMM), Global Water Partnership (GWP), sep-
tiembre, 2006, pp. 1-14.
[Salas, 2013] S ALAS M., J., Cadenas de Markov desde un punto de vista de Aplicacio-
nes, Tesis Licenciatura, Universidad Veracruzana, 2013.
[Samuelson, 2006] S AMUELSON , R. J., The Washington Post, Julio, 2006, USA.
[SJSU, 2001] SJSU, Pontryagins Maximum Principle. Notas San Jose University. www-
sjsu.edu/faculty/watkins7pontryag.html. Julio, 2001.
[SRH, 1971] S ECRETAR ÍA DE R ECURSOS H IDR ÁULICOS, Región Hidrológica No. 30
Grijalva-Usumacinta, México, 1971.
[Torres, Sandoval, 2015] TORRES , A., S ANDOVAL , A., Avances en hidrologı́a urbana,
Editorial Pontificia Universidad Javeriana, Colombia, 2015.