TEMÁTICA - revistapolitecnica.epn.edu.ec · Se autoriza la reproducción total o parcial de su...

Comisión de Apoyo Editorial

Ing. Javier Páez G.

Diseño de Portada

Sr. Diego Flores

TEMÁTICA Y ALCANCE La Revista Politécnica es una publicación periódica semestral, editada por la Escuela Politécnica Nacional del Ecuador, cuyo objetivo

es contribuir al conocimiento científico y tecnológico, mediante la publicación de estudios científicos relacionados a las áreas de

ciencias básicas (física, química y matemática) e ingenierías (agroindustria, ambiental, civil, eléctrica, electrónica, geología, mecánica,

petróleos, sistemas y química). La Revista Politécnica está dirigida a profesionales e investigadores que trabajan en estos campos del

conocimiento.

La Revista Politécnica está incluida en Latindex, catálogo y directorio: Sistema Regional de información en línea para Revistas

Científicas de América Latina, el Caribe, España y Portugal. Además la Revista Politécnica se encuentra en SciELO - Scientific

Electronic Library Online.

Se autoriza la reproducción total o parcial de su contenido siempre y cuando se cite la fuente. Los conceptos expresados son de

responsabilidad exclusiva de sus autores.

CONSEJO EDITORIAL

Coordinador Técnico Operativo

Ing. Javier Páez G.

[email protected]

Teléfono: (+593) 2976300 ext. 5220

EDITOR Oscar Eduardo Camacho, Ph.D.

Escuela Politécnica Nacional

[email protected]

Enio Da Silveira, Ph.D.

Universidad Católica de Río, Brasil.

Carlos Smith, Ph.D.

University of South Florida, Estados

Unidos

Gyimah-Brempong Kwabena, Ph.D.

University of South Florida, Estados

Unidos

Raymundo Forradelas, Ph.D.

Universidad Nacional del Cuyo,

Argentina

Ricardo Carelli, Ph.D.

Universidad Nacional de San Juan,

Argentina.

Vanderlei Bagnato, Ph.D.

Universidad de Sao Paulo, Brasil.

Rui Pedro Pinto de Carvalho, Ph.D.

University of Coimbra, Portugal

Vicenzo Vespri, Ph.D.

Università degli studi di Firenze, Italia

Oscar Ortiz, Ph.D.


Argentina

Gustavo Scaglia, Ph.D.


Argentina

Chen Ning, Ph.D.

Universidad de Mineralogía y

Tecnología de China, China.

Alex Ruiz Torres, Ph.D.

Universidad de Puerto Rico, Puerto

Rico.

CO-EDITORA Silvana Ivonne Hidalgo Trujillo, Ph.D.

Escuela Politécnica Nacional

[email protected]

Lizandro Solano, Ph.D.

Universidad de Cuenca, Ecuador

Romel Montufar, Ph.D.

Pontificia Universidad Católica,

Ecuador

Marcos Villacís, Ph.D.

Escuela Politécnica Nacional, Ecuador

Andrés Rosales, Ph.D.


Danilo Chávez, Ph.D.


Oscar Camacho, Ph.D.

Universidad de Los Andes, Venezuela

Carlos Ávila, Ph.D.


Rector

Jaime Calderón, MBA

Vicerrector de Investigación y

Proyección

Alberto Celi, Ph.D.

Vicerrector de

Docencia

Tarquino Sánchez, MBA

AUTORIDADES

ESCUELA POLITÉCNICA NACIONAL

mailto:[email protected]



PALABRAS DEL EDITOR

Esta edición de la Revista Politécnica muestra significativos trabajos realizados por importantes

investigadores sobre diversos temas. Estas diferentes contribuciones muestran un amplio alcance con un rico

contenido científico, algunos de ellos poseen una aplicación potencial en nuestro entorno.

El primer artículo es realizado por Mayorga y Rivera, quienes presentan una investigación asociada a los

Sistemas de Información Geográfica (GIS) como herramienta clave para la eficiente planeación y respuesta

ante la ocurrencia de derrames de petróleo en el Ecuador. Este trabajo pretende determinar las zonas proclives

a derrames de hidrocarburos, establecer las zonas vulnerables que podrían ser afectadas cuya conservación es

prioritaria y sugerir alternativas para reducir los tiempos de respuesta frente a emergencias.

En el segundo trabajo, Oscullo y Romero, muestran un análisis para el periodo 2011- 2015 mediante la

evolución de indicadores electro-energéticos del Sistema Nacional Interconectado (SNI) del Ecuador

considerando las diferentes fuentes disponibles de producción de energía eléctrica, así como también del

consumo. Seleccionando adecuadamente un conjunto de indicadores que son función de la expansión y

operación de los distintos recursos de generación, permiten, de una manera gráfica, determinar la robustez

para un sistema eléctrico. Esta disposición gráfica del valor en porcentaje de cada indicador en un eje, se lo

conoce como la herramienta denominada “Rosa de Robustez”; donde cada indicador a ser considerado se

obtiene en base a la información disponible en documentos oficiales de las instituciones del sector energético

En el tercer documento, Armijos-Abendaño y colaboradores, hacen un estudio de la curva de rotación de la

galaxia NGC 7331 situada a 14,7 Mpc de la Tierra. La curva de rotación se deriva usando observaciones

radioastronómicas del monóxido de carbono (CO. Se revela que el ensanchamiento de la línea de CO,

transición 2-1, está dominado por efectos turbulentos del gas de CO antes que por efectos térmicos. Asimismo,

se estudia el campo de velocidades de NGC 7331, lo que pone en evidencia la rotación de la galaxia en el

sentido de las agujas del reloj.

El cuarto trabajo es realizado por Quistial y colaboradores. Ellos presentan en su artículo modelos para el

cálculo de pérdidas de trayectoria en ambientes con línea de vista y sin línea de vista para microceldas en la

banda de 900MHz. Para obtener el modelo de pérdidas de trayectoria utilizan el método de aproximación por

mínimo cuadrado y lo determinan a partir de mediciones realizadas en la ciudad de Quito de potencia recibida

en diferentes puntos de una zona de cobertura. Las expresiones de los nuevos modelos obtenidos se basan en

los modelos existentes de espacio libre, Okumura, Okumura - Hata, COST - 231, Egli y Walfisch. Los métodos

de ajuste de los modelos existentes a las mediciones fueron el del error cuadrático medio y el método simple.

En el quinto artículo, Aguilar–Jaramillo y Aguinaga describen un modelo eléctrico equivalente de una celda

electrolítica simple, utilizada para la producción de hidrógeno a partir de vapor sobrecalentado de agua, con

la finalidad de brindar una base teórica que pueda ser aplicada al momento de elegir o diseñar equipos para

sistemas de electrólisis de este tipo. Para ello realizan el estudio de los componentes que conforman el

electrolizador, las leyes y ecuaciones que rigen el comportamiento del sistema, las variables que inciden en la

producción de hidrógeno, así como la espectroscopia de impedancia electroquímica, esta última útil para

interpretar de forma matemática los fenómenos que se presentan en el proceso de electrólisis. La

implementación y simulación del modelo análogo por circuito eléctrico equivalente (EEC) de la celda

electrolítica, utilizando las características de un electrolito de YSZ, un ánodo de LSM y un cátodo de Ni -

YSZ, permitió concluir que se pueden obtener resultados eficientes en la producción de hidrógeno a mayor

área transversal del electrolito y menor separación entre los electrodos, al trabajar con una señal de excitación

de pulsos a frecuencia de resonancia y a temperaturas de vapor sobrecalentado altas

Finalmente, Casa y colaboradores realizan un estudio que tiene como objetivo principal desarrollar la

modelación numérica del flujo rasante en una rápida escalonada aplicando el paquete comercial FLOW -3D.

En la actualidad el diseño de este tipo de estructuras se realiza con el uso de expresiones empíricas obtenidas

con base en la modelación física y estudios complementarios en la modelación numérica del flujo sobre la

rápida escalonada con apoyo de un código CFD. Con el modelo numérico obtenido se busca estimar la

velocidad del flujo en la región uniforme, y el coeficiente de fricción para cuatro caudales de operación de la

rápida escalonada.

Deseamos que el contenido de este volumen sea de interés para los lectores de la Revista Politécnica.

Oscar Eduardo Camacho Quintero, Ph.D.

EDITOR

CONTENIDO

Vol. 41, No. 2

MAYO – JULIO 2018

7

Mayorga, Henry Santiago; Rivera, José Luis

Environmental Risks for Oil Spills in Northeastern Ecuador using GIS: Response

Times and Vulnerable Areas

Riesgos Ambientales por Derrames de Petróleo en el Nororiente Ecuatoriano usando

GIS: Tiempos de Respuesta y Zonas Vulnerables

15

Oscullo, José; Romero, Luis

Determinación de la Rosa de Robustez para la Matriz Eléctrica del Ecuador

Determination of the Rose of Robustness for the Electrical Matrix of Ecuador

23

Armijos-Abendaño Jairo; López Ericson; Llerena Mario; Aldas Franklin

Cinemática y Masa dinámica de la Galaxia NGC 7331

Kinematics and Mass of the NGC 7331 Galaxy

29

Quistial, Alvin; Lupera Morillo, Pablo; Tipantuña, Christian; Carvajal, Jorge

Modelo Matemático Adaptado para el Cálculo de Pérdidas de Propagación en la

Banda de 900 MHz para Microceldas en la Ciudad de Quito

An 900 MHz Adapted Path Loss Model for Microcells in Quito

37

Aguilar-Jaramillo, Edwin; Aguinaga, Álvaro

Modelamiento y Simulación de la Producción de Hidrógeno en un Electrolizador a

Partir de Vapor Sobrecalentado de Agua

Modeling and Simulation of Hydrogen Production in an Electrolyzer from Superheated

Water Vapor

53

Casa E.; Hidalgo X.; Castro M.; Ortega P; Vera P

Modelación Numérica del Flujo Rasante en una Rápida Escalonada Aplicando la

Dinámica de Fluidos Computacional (CFD) Mediante el Uso de Flow-3D

Numerical Modeling of Flush Flow in a Rapid Step Applying Computational Fluid

Dynamics (CFD) Using Flow-3D

Environmental Risks for Oil Spills in Northeastern Ecuador using GIS: Response Times and Vulnerable Areas

Revista Politécnica - JULIO 2018, Vol. 41, No. 2

[email protected] Recibido: 17/08/2017

Aceptado: 20/07/2018

Publicado: 31/07/2018

11. INTRODUCTION

Ecuador is one of the most biodiverse countries in the world,

considering the number of species per unit area (Bravo, 2014).

Its territory has a set of ecological, geographical and climatic

conditions to denominate it as a megadiverse country (Finding

Species, 2010). This diversity has become a strategic axis for

Ecuador. Therefore, proposals are being implemented to

conserve their wealth and to exploit resources in a rational way

in order to ensure sustainable development (IGM, 2014).

Every industrial process implies the generation of a certain

impact on the environment that surrounds it (Escrig, 2008).

However, activities such as exploration, production and

transportation of oil are prone to cause major damage.

Throughout the history of Ecuador, there have been several

cases of environmental disasters related to oil industry,

especially oil spills.

Environmental Risks for Oil Spills in Northeastern Ecuador using

GIS: Response Times and Vulnerable Areas

Mayorga, Henry Santiago1; Rivera, José Luis1

1Escuela Politécnica Nacional, Facultad de Ingeniería en Geología y Petróleos, Quito, Ecuador

Abstract: In the present investigation Geographic Information Systems (GIS) are used as a key tool for an efficient

planning and response to the occurrence of oil spills in Ecuador; considered as a megadiverse country. It is intended

to identify areas prone to oil spills, to identify vulnerable areas which conservation is a priority, and to suggest

alternatives to reduce response times to these emergencies. The area of study corresponds to seven blocks of the

Northeastern Ecuador, an area of immense ecological value, in which oil activities are developed. GIS was used to

perform spatial analyzes of wells, platforms, production stations, oil pipelines and environmental factors such as

protected areas. In addition, a statistical analysis was carried out to determine if the proximity to roads, houses, towns,

rivers and flood zones are risk factors that could cause oil spills. An innovative proposal is presented to improve the

prior planning and remediation of oil spills by identifying patterns of occurrence and vulnerable zones with GIS; thus

minimizing environmental and social impact through rapid response. It also allows operating companies to update

their contingency plans to respond effectively to new emergencies and optimize their resources, providing a new

alternative to conserve the enormous biodiversity of Ecuador.

Keywords: GIS, Oil Spills, Ecuador.

Riesgos Ambientales por Derrames de Petróleo en el Nororiente

Ecuatoriano usando GIS: Tiempos de Respuesta y Zonas

Vulnerables

Resumen: En la presente investigación se utilizan Sistemas de Información Geográfica (GIS) como herramienta

clave para una eficiente planeación y respuesta ante la ocurrencia de derrames de petróleo en el Ecuador; considerado

como un país megadiverso. Se pretende determinar las zonas proclives a derrames de hidrocarburos, establecer las

zonas vulnerables que podrían ser afectadas cuya conservación es prioritaria y sugerir alternativas para reducir los

tiempos de respuesta frente a estas emergencias. La zona de estudio corresponde a siete bloques del Nororiente

ecuatoriano, área de inmenso valor ecológico a nivel mundial en la cual se desarrollan actividades petroleras. Se

utilizó GIS para realizar análisis espaciales de pozos, plataformas, estaciones de producción, oleoductos y factores

ambientales como áreas protegidas. Además, se realizó un análisis estadístico para determinar si la cercanía a vías,

casas, poblados, ríos y zonas inundables constituyen factores de riesgo que puedan ocasionar derrames. Se presenta

una propuesta innovadora para mejorar la planificación previa y remediación de derrames de petróleo al identificar

patrones de ocurrencia y zonas vulnerables mediante GIS; minimizando con ello el impacto ambiental y social

mediante una rápida respuesta. Además permite a las empresas operadoras actualizar sus planes de contingencias para

responder eficazmente ante nuevas emergencias y optimizar sus recursos, convirtiéndose en una nueva alternativa

para conservar la enorme biodiversidad del Ecuador.

Keywords: GIS, Derrames de Petróleo, Ecuador.

7



An oil spill represents the discharge (intentionally, by accident

or by incorrect work practices) of hydrocarbons into the

environment (Udoh, 2010). Oil spills are a major source of

human contamination, as they can spread rapidly depending on

the type of hydrocarbons and topographical and climatological

conditions in the area (Zambrano, 2015). Effects produced by

a spill can be persistent over time, putting lives of humans and

species at risk (ITOPF, 2011).

Insufficient investment in infrastructure and research, lack of

legislation and environmental education have been the main

factors for developing countries not having considered the

damaging effects generated by an oil spill (Lawler, Horst, &

Champan, 2011). Specifically, in Ecuador, there are several

historical causes that have caused oil spills such as natural

disasters, technical and human failures, lack of environmental

awareness of domestic and foreign companies, poor

maintenance of oil facilities, sabotage and vandalism

(Mohamadi, 2015).

Although we can not predict the occurrence of an oil spill, nor

the degree of involvement that could cause in the ecosystem,

we must be prepared by prevention and quick response tools

to reduce possible environmental risks (Udoh, 2010). The

updated contingency plans are the main procedures for dealing

with oil spills (Mohamadi, 2015). A contingency plan is a

program to achieve an immediate, organized and effective

response to an emergency within the different phases of the oil

industry; which details the actions, equipment and

responsibilities of personnel to control, remedy and minimize

the damage caused by an oil spill (Torres, 2011).

The areas where oil activities in Ecuador are developed,

particularly protected areas, represent zones of immense

ecological value worldwide. So, it is necessary to make

environmental sensitivity maps that could be of vital

importance for conservation of planet's biodiversity

(OILWATCH, 2004).

These maps can be made using Geographic Information

Systems (GIS), software used to manage geographic data that

allows the integration of different types of information and

communicate the results obtained visually (IPIECA, 2012).

The advantage of GIS is that it improves the pre-contingency

planning of an oil spill by establishing priority areas to be

conserved and facilitates communication through a rapid

information link (IPIECA, 2012). The use of GIS helps to

minimize the effects that a hydrocarbon spill may cause on the

environment and reduce the economic costs that should be

incurred by companies in their remediation.

Secretary of Hydrocarbons, entity of regulation of oil activities

in Ecuador, divided the country territorially into 83 areas

called “oil blocks” (Figure 1). This division was carried out in

accordance with Hydrocarbons Law, which establishes a space

of no more than 200 000 hectares for the awarding of oil

exploration and production contracts to operating companies.

These oil blocks are designated by a number, a name and its

operating company.

This research intends to identify areas prone to oil spills, to

determine vulnerable areas that could be affected and establish

optimal response times; analyzing patterns of occurrence of oil

spills based on historical data in seven oil blocks of the

Northeastern Ecuador. The results obtained in the present

project will be of great relevance for state and private

companies involved in the oil industry in Ecuador owing to the

lack of similar studies in the area of analysis.

2. METHODOLOGY

2.1 Study Area

The zone of analysis is located in Northeastern Ecuador

(Figure 1). The oil blocks studied were Block 12 (Eden-

Yuturi), Block 15 (Indillana), Block 58 (Cuyabeno-Tipishca)

and Block 59 (Vinita) operated by Petroamazonas EP, Block

16 (Iro) operated by Repsol, Block 53 (Singue) operated by

Gente Oil and Block 62 (Tarapoa) operated by Andes

Petroleum. These companies were in charge of the operation

of the blocks of analysis at the time of the study.

Figure 1. Oil Blocks in Ecuador and Study Area.

2.2 Data acquisition

The Environmental and Social Reparation Program of the

Ministry of Environment (PRAS-MAE) is the state agency that

keeps the official count of oil spills in Ecuador. This institution

makes a complete survey of information to determine its

causes in situ and propose remediation programs. It is a legal

obligation to operating companies to report oil spills and apply

contingency measures to reduce possible environmental and

social effects.

A total of 1168 oil spills have been registered from 1972 to

2015 in Ecuador by PRAS-MAE. Of this total, 61 oil spills are

located within the investigated blocks, as shown in Figure 2.

Wells, platforms and production stations were obtained from

PRAS-MAE database. Oil pipelines, roads, rivers, houses and

villages were downloaded from the website of Instituto

Geográfico Militar (IGM, 2015). The National Information

System provided the layers of mass movements and flood

zones (SNI, 2010). Finally, protected areas were drawn by

georeferencing the map of the National System of Protected

Areas of Ecuador. All these parameters will be used in later

analyzes.

8



Figure 2. Location of Oil Spills in the Study Area.

2.3 Data processing

Spatial analyzes were made using ArcGIS version 10.3. The

processes and tools used in each of the stages of study are

detailed below.

Causes of Oil Spills: Selections by attributes were used

to establish the causes of oil spills in the study area and

to determine what type of oil facility was associated

(well, platform or oil pipeline). A generalized linear

model for binary data was then performed using the

program SPSS version 23. A total of 74 possibilities (61

oil spills, 14 random control points and 1 case excluded

by the program) were taken, in which it is considered the

probability that an oil spill occurs as a dependent variable

and as categorical factors to its proximity to roads,

houses, towns, rivers and flood zones. This statistical

analysis allowed to determine if these environmental and

social variables represent risk factors that could directly

influence in the occurrence of oil spills.

Determination of Oil Spill Zone Prone (OSZP): The

maximum impact radius reported in the 61 oil spills is

approximately 100 meters, so that buffers of 300 meters

were created around the points of location of these

events. Then, areas with a pattern of occurrence prone to

oil spills were drawn.

Determination of Vulnerable Areas: Vulnerable areas

are the zones of primary protection when an oil spill

ocurrs. For their identification, environmental sensitivity

factors such as protected areas, flood zones, rivers,

houses and villages in each block were superimposed on

a map.

Response Times for Oil Spills: Information in the

contingency plans by block was analyzed. Then, stations

where it would be recommended to locate human and

material resources, as well as the type of transport

suitable to reduce response times against an oil spill were

proposed. For a vehicle traveling by land, routes were

drawn along the roads from the chosen station to the most

distant points in which a particular oil facility is located.

Also considering an average speed of 30 km/h (average

speed of heavy transport on the roads in the East of

Ecuador), maximum times of response were calculated.

Finally, to determine the time it would take to arrive by

helicopter at the most distant point of an oil facility in

each block of analysis, trajectories of rectilinear flight

were drawn and it was considered that in this type of

emergencies the average speed is 180 km/h (Exxon,

2008).

3. RESULTS

As can be seen in Figure 3, of the 61 oil spills occurred in the

study area, 50 of them are registered in Cuyabeno-Tipishca

Block, 5 in Tarapoa Block, 3 in Indillana Block, 2 in Iro Block

and 1 in Vinita Block. While in Eden-Yuturi and Singue Block

no oil spills have been recorded.

Figure 3. Number of Oil Spills per Block.

3.1 Causes of Oil Spills

Of the total oil spills registered, 19 were caused by corrosion

(31 %), 16 were caused by mechanical failures (26 %), 12

without determining their cause (20 %), 8 of them occurred by

attacks (13 %), 6 for human failures (10 %) and there was no

event that originated because of natural disasters. With respect

to the oil facility associated to each oil spill, it was verified that

20 occurred in production stations (33 %), 1 originated in flow

lines (2 %), 38 oil spills are associated to wells (62 %) and 2

of them occurred in platforms (3 %), as shown in Figure 4.

In the generalized linear model it was obtained that the

significance values for each categorical factor are greater than

0,05 (ρ> 0,05), so it is deduced that at least for the socio-

environmental factors considered there is no significant

influence that explains the reason why oil spills occur in the

study area. The ρ values are presented in Table 1.

9



Figure 4. Causes of Oil Spills and associated Petroleum Facility.

Table 1. ρ values per categorical factor studied.

Categorical Factor ρ

Houses 0,999

Towns 1,000

Roads 0,999

Rivers 0,764

Flood Zones 0,856

3.2 Determination of Oil Spill Zone Prone (OSZP)

An OSZP was identified in the study blocks, whose

approximate area is 52 square kilometers. This area is located

in the Cuyabeno-Tipishca Block as shown in Figure 5. In this

zone have occurred 43 of the 61 oil spills, that is, 70,5 % of

total oil spills are recorded in this zone.

Of these, 43 oil spills were determined that 13 were caused by

corrosion, 12 were caused by a mechanical failure, 7 spills

have not been established, 7 of them occurred by attacks and 4

by human failures. Considering the petroleum facility to which

they are related, 13 spills are associated with production

stations and 30 occurred in wells.

Figure 5. Zone prone to the occurrence of Oil Spills.

3.3 Determination of Vulnerable Areas

124 villages, around 6000 houses and 3 protected areas of

global ecological importance (Cuyabeno Fauna Production

Reserve, Yasuni National Park and Limoncocha Biological

Reserve) were identified as shown in Figure 6a. In addition,

flood zones and the rivers that cross the studied oil blocks

proposed in Figure 6b.

Figure 6. Vulnerable Areas.

3.4 Response Times for Oil Spills

It was determined that Sansahuari and Cuyabeno Station,

being within the OSZP and where 65 % of the events

associated with production stations occurred, represent the

best options for locating the equipment for oil spill response

throughout this area. The use of vehicles of heavy transport by

land constitutes the means of transport advisable for this zone.

The maximum time it would take to get from both stations to

the furthest point in this area is 13 minutes at an average speed

of 30 km/h (Figure 7).

In addition, two alternatives are proposed by air for an

effective response to all research blocks and areas without a

likely pattern of occurrence, using a helicopter containing spill

containment and remediation equipment (Figure 8a and 8b).

One option is to locate the helicopter in Tarapoa Base Camp,

which would achieve a specific point in an average time of

approximately 17 minutes. The other alternative is to locate it

in the Cuyabeno Station, taking on average 19 minutes to reach

a future spill. This is a reasonable time considering the most

distant point of the covered area and depending on the speed

of the helicopter, this time could be reduced.

Figure 7. Response Times for Oil Spills in the OSZP.

10



Figure 8. Response Times by Helicopter.

4. DISCUSSION

The area of study is a place of immense biological value

because it presents unique ecosystems in the world and a high

degree of endemism of species (Acción Ecológica, 2006), but

also represents a site of economic interest because of oil

reserves found. As a result, this research is of great importance

for responsible oil exploitation in Ecuador because areas with

the highest probability of occurrence of oil spills and

vulnerable zones that are of priority conservation have been

identified, thus improving response times.

The Cuyabeno - Tipishca Block registers the greatest number

of spills (82 % of the total). This block represents one of the

main mature fields of exploration and production of Ecuador,

with an approximate production of 25 000 BOPD. However, it

is verified that there is no relationship of dependence among

the production of a given block with respect to the number of

spills. An example is the Eden-Yuturi Block with a production

of around 50 000 BOPD but does not record the occurrence of

any emergencies caused by oil spills (ARCH, 2013).

4.1 Causes of Oil Spills

According to the BBC (2013), its study was to determine the

causes of oil spills in Ecuador, it is estimated that 539 spills

have been recorded from 2000 to 2010. The results showed

that 27 % of them occured through corrosion, 26 % for attacks,

18 % for mechanical failures, 15 % without data, 12 % for

human failures and 2 % for natural disasters. While in this

investigation, it was obtained that corrosion and mechanical

failures represent the causes with greater percentage.

Consequently, the two studies agree that corrosion is the main

cause of occurrence of oil spills in Ecuador.

In addition, the generalized linear model indicates that at least

for the categorical factors considered (houses, towns, roads,

rivers and flood zones) there is no direct risk relationship in

the occurrence of oil spills. Oil spills in the study area are

mostly caused by operational problems owing to the lack of

adequate and constant maintenance of the equipment and

installations by the operating companies.

Considering the oil facility in which oil spills have occurred,

the largest percentage occured in wells and production

stations. Contrary to what was presented in the analysis by

Achebe et al. (2012), which determines that oil pipelines

represent the facility with the highest number of associated oil

spills and despite the fact that main oil pipeline of Ecuador,

SOTE, has exceeded its optimal life (BBC, 2013), in the

investigated blocks there are only 2 % of them occurring in

flow lines.

The high percentage of oil spills for attacks is because of the

proximity of the study area with respect to the border with

Colombia, the possible presence of insurgent groups,

settlements of uncontacted indigenous communities, that the

operating companies dissociate themselves from any

responsibility and possible sabotage of people to obtain

compensation. In addition, most of the area presents lands of

low susceptibility to mass movements so there is a low

percentage of oil spills caused by natural disasters.

The number of oil spills related to human failures could be

significantly reduced through adequate staff training and

strong policies against alcohol and drug use during workdays

(DeCola & Sierra, 2006).

4.2 Determination of Oil Spill Zone Prone (OSZP)

The OSZP is located in the southern part of the Cuyabeno-

Tipishca Block, very close to the Singue Block. It represents

only 3 % of the total territory of Block 58. Although its area is

minimal compared to the total of the block, this zone has the

most oil facilities and human settlements. In this zone has

occurred a 70,5 % of the total events, which demonstrates an

area with a high probability of occurrence of future oil spills.

This identified area is the only one that presents a clear pattern

of occurrence to oil spills.

4.3 Determination of Vulnerable Areas

Ecuadorian Amazon is a very fragile habitat and events like oil

spills can cause ecological and economic consequences,

affecting fauna and flora, tourism and the communities that

inhabit it (BBC, 2013). This research found that the study area

is highly vulnerable resulting from the large number of houses

and villages near the oil facilities, the presence of 3 protected

areas of global ecological importance and rivers that could

cause emergencies of greater magnitude.

There are three fundamental problems that could lead to

catastrophic consequences in the study blocks. The first is the

fact that the oil spill zone prone is being limited to the

Cuyabeno Reserve and that oil spills have been recorded

within protected areas. The Yasuní National Park, with almost

one million hectares, is the most biodiverse place on the planet

and a symbol of preservation for its ecological and cultural

wealth (Larrea, 2011). Likewise, the Cuyabeno Fauna

Production Reserve and the Limoncocha Biological Reserve

are sensitive and priority conservation areas in case of the

occurrence of oil spills, as they are home to one of the largest

concentrations of wildlife, both in flora and fauna, many of

them in danger of extinction (Acción Ecológica, 2006).

The second problem is that there are 124 villages in the area

under investigation, including uncontacted indigenous

11



communities (such as Huaoranis, Tagaeris and Toromenanes),

who would be seriously affected by the occurrence of an oil

spill. In addition to causing diseases in people, depending on

the impact of the spill, communities that depend on fishing and

rivers as a source of drinking water would have great

implications for their survival (Chang et al., 2014).

Finally, it was observed that the Cuyabeno River crosses the

oil spill zone prone and flows directly into the Cuyabeno

Reserve. Also, most of the blocks contain permanently flooded

areas; this could generate a rapid dispersion of the spilled oil

and cause further environmental damage.

4.4 Response Times for Oil Spills

The great vulnerability of the study area altogether with

environmental conditions that favor a rapid spread of

hydrocarbons, it is necessary a prompt and effective response

to reduce the impacts that oil spills could generate. For this, all

information of the contingency plans, such as equipment and

personnel necessary to deal with these events, should be

carefully analyzed. But, despite being a requirement for the

award of the environmental license prior to any operation in

the oil industry in Ecuador (Constitution of Ecuador, 2008), it

is incomprehensible that contingency plans of all the analyzed

blocks are not available and they are not properly updated. For

an adequate ecological protection, there is a need for localized

GIS operations with local data to assist in immediate response.

This is the case of the Cuyabeno-Tipishca Block that is of

special interest since in its territory is the OSZP and whose

contingency plan establishes only the equipment required for

the containment of spills, but their location is not specified. So,

it has been recommended to locate them in strategic

production stations that allow obtaining appropriate response

times for this type of emergencies.

The best proposal to reduce the response time in the OSZP is

to place the equipment in the Sansahuari and Cuyabeno

Stations and to use vehicles of heavy transport by land. The

maximum estimated time is 13 minutes; this improves the

response time suggested in the contingency plan of the

Cuyabeno-Tipishca Block. While two alternatives were

proposed to contain oil spills by means of a helicopter, either

to locate it in the Cuyabeno Station or in the Tarapoa Base

Camp. In both cases response times are low (20 minutes on

average) and it could reach any point where an oil spill could

be generated. This approach was made in view of the fact that

there are not enough roads to reach all oil facilities by land and

it would take too long. The only limitation is that the blocks

are operated by different companies, which requires the

intervention of the national environmental authority to achieve

a mutual agreement among them. In addition, it is more

economically feasible to locate a helicopter in the whole area

than one per block.

Through these proposals, it is possible to deal with oil spills

effectively by reducing response times and involving all

companies in order to achieve sustainable oil exploitation in

Ecuador.

5. CONCLUSIONS

A total of 61 oil spills have been recorded in the study area.

The main causes are corrosion and mechanical failures; most

of the oil spills occurred in wells and production stations. It is

concluded that the oil spills in this area are caused by

operational problems because of the lack of adequate and

constant maintenance of the equipment and installations by the

operating companies. The best way to fight against the

environmental impact of oil spills is through prevention.

Therefore, government of Ecuador should invest more

resources on continuous monitoring and supervision on oil

infrastructure.

The OSZP is located in the southern part of the Cuyabeno-

Tipishca Block and represents an area with a high probability

of occurrence of oil spills. It occupies 3 % of the total territory

of Block 58 with an approximate area of 52 square kilometers.

There are very fragile habitats, highly vulnerable areas and

favorable conditions for a rapid spread of hydrocarbons, so an

immediate response is required to minimize possible

environmental risks that a spill may cause. In the OSZP, it is

proposed to locate the containment equipment in the

Sansahuari and Cuyabeno Stations and to use heavy transport

by land; this reduces to a maximum response time of 13

minutes. For the remaining study area, two alternatives are

proposed through the use of a helicopter. Placing it in the

Cuyabeno Station or the Tarapoa Base Camp, the response

time is 20 minutes on average for the entire investigated area.

This research should be extended to all oil blocks in Ecuador

because it allows operating companies to update their

contingency plans to respond effectively to the occurrence of

oil spills. For the correct use of GIS will be necessary that the

companies invest in a constant training of their personnel.

As a complement to this investigation, it is recommended to

make a numerical model for the prediction and dispersion of

oil spills in the study area.

REFERENCES

Acción Ecológica (2006). Derrame de petróleo en el Cuyabeno revela crítica

situación en que se encuentran las áreas protegidas. Retrieved from

http://www.accionecologica.org/.

Achebe, C., Nneke, U., & Anisiji, O. (2012). Analysis of Oil Pipeline Failure in the Oil and Gas Industries in the Niger Delta Area of Nigeria. IMECS,

2(1), 1-6.

Acosta, A., & Schuldt, J. (2014). Petróleo, Rentismo y Subdesarrollo: ¿Una Maldición sin Solución? Revista Latinoamericana de Comunicación

Chasqui, 94(1), 71-88.

ARCH (2013). Boletín Estadístico de Producción de petróleo por Bloiques del año 2013. Quito, Ecuador.

Astudillo, B. (2007). Actualización del Plan de Contigencias para las

Actividades Hidrocarburíferas en el Área Auca (Undergraduate Thesis). Escuela Politécnica Nacional.

BBC (2013). ¿Por qué hay un derrame petrolero por semana en Ecuador?.

Retrieved from

www.bbc.com/mundo/noticias/2013/06/130610_ciencia_ecuador_derra

me_rio_limpieza_ig.

12

http://www.accionecologica.org/



Biniwale, S. (2003). GIS, Effective Support for Oil Spill Contingency

Planning. Society of Petroleum Engineers. doi:10.2118/81206-MS

Bravo, E. (2014). La Biodiversidad en el Ecuador. Cuenca, Ecuador, Abya -

Yala. Cajamarca, O. (2014). Estudio Técnico - Económico del Cañoneo PERF STIM

en el Yacimiento "M2" del Campo Eden Yuturi (Undergraduate Thesis).

Universidad Tecnológica Equinoccial. Chang, S., Stone, K., Demes, K. & Piscitelli, M. (2014). Analysis of Oil

Pipeline Failure in the Oil and Gas Industries in the Niger Delta Area of

Nigeria. Ecology and Society, 19(2), 26-28. Constitución de la República del Ecuador (2008).

DeCola, E. & Sierra F. (2014). An Assessment of the Role of Human Factors

in Oil Spills from Vessels.. Alaska, Estados Unidos. Escrig, D. (2008). El Impacto Ambiental de las Actividades Industriales: El

Cambio Necesario. Andalucía, España: UBE Corporation Europe.

Exxon (2008). Oil Spill Response Manual. Texas, Estados Unidos, ExxonMobil.

Fernández , L., & Gutiérrez, M. (2013). Bienestar Social, Económico y

Ambiental para las Presentes y Futuras Generaciones. Información

Tecnológica, 24(2), 122-130.

Finding Species (2010). Armonía Ecuador. Quito, Ecuador, Ecuador Verde.

Guillaume, F. (2002). Sobre Bonanzas y Dependencia: Petróleo y Enfermedad Holandesa en el Ecuador. Iconos, 14(1), 1-9.

IGM (2014). El Medio Ambiente. Quito, Ecuador, Instituto Geográfico Militar.

IGM (2015). Geoportal. Retrieved from http://www.geoportaligm.gob.ec/portal/.

IPIECA (2012). Sensitivity Mapping for Oil Spill Response. Londres, Inglaterra, International Maritime Organization.

ITOPF (2011). Efectos de la Contaminación por Hidrocarburos en el Medio

Marino. Londres, Inglaterra, International Tanker Owners Pollution Federation Limited.

Jankilevich, S. (2003). Las Cumbres Mundiales sobre el Ambiente:

Estocolmo, Río y Johannesburgo. Buenos Aires, Argentina, Universidad de Belgrano.

Larrea, C. (2011). La Iniciativa Yasuní-ITT desde una perspectiva

Multicriterial. Quito, Ecuador, Universidad Simón Bolívar. Lawler, M., Horst, D., & Champan, L. (2011). Attacks on oil transport

pipelines in Nigeria: A quantitative exploration and possible explanation

of observed patterns. Applied Geography, 32(2), 636-651. Mohamadi, B. (2015). GIS Based Oil Spill Risk Assessment Model for the

Niger Delta´s Vegetation. Nature Environment and Pollution Technology,

14(3), 545-552. OILWATCH. (2004). Áreas Protegidas ¿Protegidas contra quién? Quito,

Ecuador, World Rainforest Movement.

SNI (2010). Información Geográfica. Retrieved from http://sni.gob.ec/coberturas.

Torres, C. (2011). Plan de Contingencias para Derrames de Hidrocarburos

en las Líneas de Flujo en el Campo Cuyabeno de Petroproducción (Undergraduate Thesis). Universidad de Guayaquil.

Udoh, J. (2010). A Gis Based Cost Distance Modeling for Oil Spill Hazard

Assessment in the Coastal Areas of South Eastern Nigeria. Research Journal of Applied Sciences, 5(4), 254-259.

Zambrano, J. (2015). Control de Derrames de Hidrocarburos. Quito, Ecuador,

Escuela Politécnica Nacional.

BIOGRAPHIES

Henry Santiago Mayorga Mayorga.

Nació en Quero-Ecuador el 5 de

Septiembre de 1993. Culminó sus

estudios secundarios en el Colegio

Bolívar de Ambato. Se graduó de

Ingeniero en Petróleos en la Escuela

Politécnica Nacional a inicios del 2017.

Actualmente labora como ingeniero de

reservorios en Grupo Synergy E&P. Está interesado en iniciar

sus estudios de posgrado.

José Luis Rivera Parra. Nació en Quito-

Ecuador el 22 de Julio de 1984. Se graduó

de Licenciado en Biología en el 2007 en

la Pontificia Universidad Católica del

Ecuador. Obtuvo su Maestría en Biología

en University of Missouri-St. Louis en el

2010. Y en el 2013 obtuvo su título de

Ph.D. en Ecología, Evolución y

Sistemática en University of Missouri-St.

Louis. Está interesado en modelamiento de riesgos

ambientales asociados a las actividades extractivas y cree que

se puede explotar los recursos naturales de una manera

responsable.

13






11. INTRODUCCIÓN

La característica principal de una matriz eléctrica es presentar

el nivel más alto posible de la oferta para asegurar el

suministro de energía eléctrica para el consumo de las

diferentes actividades de la sociedad. Esta característica del

sistema eléctrico depende de variables económicas y técnicas,

las mismas son derivadas de la planificación de la expansión

del sistema la cual está sujeta a decisiones del ámbito político,

económico y ambiental. Así, la seguridad de suministro de

energía eléctrica dada por un sistema eléctrico, es una de las

componentes, en búsqueda de la sostenibilidad de energía de

Determinación de la Rosa de Robustez para la Matriz Eléctrica del

Ecuador

Oscullo, José1; Romero, Luis1

1Escuela Politécnica Nacional, Facultad de Ingeniería Eléctrica y Electrónica, Quito, Ecuador

Resumen: El adecuado análisis en la expansión de la generación de un país, en un determinado periodo, permite

conocer el nivel de seguridad para el suministro de energía eléctrica; por medio de lo que se garantiza el

abastecimiento del consumo de la misma para las diferentes actividades requeridas por la sociedad. El desarrollo de

la matriz eléctrica se basa en la accesibilidad a las diferentes fuentes con fines de producción eléctrica, siempre y

cuando se cuente con un adecuado nivel de disponibilidad de la tecnología de generación considerada en la

planificación de la expansión del parque generador. Sin embargo, en cualquier país existen momentos, políticos,

económicos y ambientales que determinan de manera tácita el nivel de seguridad energético, el cual no es determinado

explícitamente en la planificación, y en la práctica es analizado considerando por lo general la dimensión técnica. Al

ser el abastecimiento un tema transversal es necesario tomar en cuenta las diferentes aristas de la expansión.

En el presente trabajo se muestra un análisis para el periodo 2011-2015 mediante la evolución de indicadores electro-

energéticos del Sistema Nacional Interconectado (SNI) del Ecuador considerando las diferentes fuentes disponibles

de producción de energía eléctrica y el consumo. Mediante un conjunto adecuadamente seleccionado de indicadores

que son función de la expansión y operación de los distintos recursos de generación; lo cuales, estructurados

adecuadamente y dispuestos gráficamente determinan la robustez para un sistema eléctrico, es decir, la adaptabilidad

del mismo ante variaciones del entorno. Esta disposición gráfica del valor en porcentaje de cada indicador en un eje,

se lo conoce como la herramienta denominada “Rosa de Robustez”; donde cada indicador a ser considerado se obtiene

en base a la información disponible en documentos oficiales de las instituciones del sector energético.

Palabras clave: Indicadores energéticos, gerenciamiento de la seguridad, economía de sistemas de potencia,

expansión de generación.

Determination of the Rose of Robustness for the Electrical Matrix

of Ecuador

Abstract: The adequate analysis in the expansion of generation for a country that allows knowing the level of security

for the supply of electrical energy, and that is guaranteed the supply of consumption for the different activities

required by society. The development of the electrical matrix is based on the accessibility to the different sources for

the purpose of electricity production that it is provided, there is an adequate level of availability of the generation

technology considered in the planning expansion of generator park. However, any country there are political,

economic and environmental moments, which tacitly determine level of energy security that is not explicitly

determined in planning, and in practice is analyzed considering technical dimension in general. As supply is a cross-

cutting issue, it is necessary to take into account the different edges of the expansion.

The present work an analysis is shown for 2011-2015 through evolution of electro-energy indicators for the National

Interconnected System (SNI) of Ecuador considering different available sources of electricity production and

consumption. Through a suitably selected set of indicators that are a function of expansion and operation different

generation resources; which, suitably structured and graphically arranged, determine the robustness of an electrical

system, its adaptability to environmental variations. This graphic layout of percentage value of each indicator on an

axis is known as the tool called " Rose of robustness "; where each indicator can be considered and must be obtained

in based on information available in official documents of institutions the energy sector.

Keywords: Power generation indicators, security management, power system economics, expansion of generation.

15



un país, la misma ha sido tratada de diversas formas, dadas las

coyunturas que cada ámbito puede ser analizado a través de

diferentes indicadores (DNETN, 2007; Blanco, 2015).

El término seguridad de suministro de energía eléctrica

representa la garantía de obtener energía eléctrica con calidad

a un nivel de precios que los consumidores puedan acceder.

Cada matriz eléctrica busca garantizar la oferta, la cual

depende de las centrales existentes, proyectos en construcción,

disponibilidad y diversidad de los combustibles de acuerdo a

las fuentes primarias que disponga o pueda acceder el país

(Retamales, 2005). En cada caso se considera una adecuada y

adaptada expansión del sistema de transmisión que permita

utilizar la energía de las distintas centrales; la cual sale del

alcance del trabajo propuesto.

Mientras que el término robustez del sistema eléctrico, en el

presente trabajo, se entiende como la adaptabilidad del mismo

ante variaciones del entorno; mediante el análisis de los

diferentes ámbitos de seguridad del suministro eléctrico.

(Molina,2005).

La matriz eléctrica de un país se planifica, aunque sea de forma

indirecta el nivel de seguridad del suministro de acuerdo a la

disponibilidad de recursos técnicos y económicos con los que

cuente el sector eléctrico, esta realidad determina que la

seguridad de suministro no trate de incrementar la

autosuficiencia del país por medio de solo un tipo de

tecnología para el suministro de energía y reducir al máximo

la dependencia energética de los recursos energéticos que no

cuente el país, sino que busque balancear la diversidad de

tecnologías de producción de energía eléctrica en la matriz

(Retamales, 2005).

Como se indicó anteriormente cada decisión ejecutada en la

matriz eléctrica puede ser evaluada mediante el

establecimiento de indicadores para determinar la seguridad en

el suministro, para lo cual es necesario contar con información

a detalle para determinar el indicador respectivo y que

adecuadamente engranados permita conocer el grado de

seguridad de suministro del sistema eléctrico de un país.

De acuerdo a la realidad del sector eléctrico de cada país, los

pocos estudios existentes respecto a este tema han planteado

un conjunto de indicadores para determinar el nivel de

seguridad de suministro eléctrico; pero solo se colocando

énfasis en uno de ellos (Retamales, 2005; Blanco, 2015).

Sin embargo, en un sistema hidrotérmico debido a que la

energía disponible depende de la aleatoriedad de los caudales

afluentes a las centrales hidroeléctricas esta situación

compromete la seguridad del abastecimiento de la demanda y

cuyo análisis debe ser ágil y rápido; estas características deben

se plasmadas en las herramientas que permitan ser observadas

de manera gráfica, esto se lo consigue a través de la ubicación

espacial de indicadores de las principales variables técnico-

económica de la expansión y operación de un sistema

eléctrico, siendo los principales: la diversidad de fuentes de

producción de energía, energía firme y el impacto económico

de contar con estas fuentes energéticas. Como segundo paso la

herramienta debe ser fácilmente replicable en cada escenario

analizado; mediante lo cual evaluar el nivel de seguridad en el

suministro eléctrico.

El conjunto de indicadores dispuestos gráficamente permite

determinar la robustez del sistema, por medio del valor de un

indicador ubicado en cada eje que en conjunto con otros

determinan un polinomio al cual se lo denomina “Rosa de

Robustez” (DNETN, 2007).

En este trabajo se presenta la aplicación de la herramienta Rosa

de Robustez para determinar el nivel de seguridad de

suministro de energía eléctrica del SNI del Ecuador en el

periodo 2011-2015 considerando los siguientes indicadores:

- La diversidad de fuentes de producción de energía,

- Energía firme, de las tecnologías de producción eléctrica

en base a fuentes autóctonas,

- Valor presente de los costos de operación e inversiones de

los proyectos en construcción, y;

- Generación del valor agregado del sector eléctrico en la

economía nacional;

Debido a que el acceso a la información requerida para

determinar los distintos indicadores debe basarse en análisis de

datos proveniente de documentos oficiales disponibles de las

instituciones del estado. Es así que en la Sección II se presenta

a la Rosa de Robustez como la herramienta que determina la

seguridad a través de un conjunto de indicadores. En la

Sección III, se introduce la metodología para la determinación

de la Rosa de Robustez. Dicha metodología es aplicada para el

periodo 2011-2015 y se obtiene la Rosa de Robustez del

sistema eléctrico ecuatoriano, comparándola con la Rosa de

Robustez del sistema eléctrico uruguayo y del sistema

eléctrico chileno para un horizonte que contenga el periodo de

estudio realizado. Finalmente, en la Sección IV, se presentan

las conclusiones y recomendaciones del trabajo, así como

futuras líneas de desarrollo.

2. INDICADOR DE ROBUSTEZ PARA UN SISTEMA

ELÉCTRICO MEDIANTE ROSA DE ROBUSTEZ

El indicador de robustez fue planteado por primera vez por la

Dirección Nacional de Energía y Tecnología Nuclear de

Uruguay para analizar diferentes estrategias de expansión de

la oferta de generación a través de la determinación de 5

indicadores para analizar las variables técnicas y económicas

de la misma (DNETN, 2007).

El indicador de cada variable analizada se coloca en un eje en

una escala de 0 a 1; cada nivel del indicador representa la

contribución a la robustez del sistema y la unión de cada uno

de los niveles alcanzados determina un polígono donde a

mayor área del mismo, es decir mayor perímetro determina la

mayor robustez, y si un indicador está ubicado internamente

afecta a la robustez del sistema; como se observa en la figura

1 (DNETN, 2007; Molina, 2005; Mosto, 2008).

16



Figura 1. Rosa de Robustez para un año t

A continuación, se presentan los indicadores de cada eje de la

rosa de robustez:

2.1 Diversidad de fuentes (IDF): Este indicador determina la

diversidad de las fuentes de producción de energía

eléctrica con que cuenta una matriz eléctrica, en la

ecuación 1 se muestra la determinación del indicador de

diversidad de fuentes.

t

t MaxIDF

1

))(

1(

(1)

En la ecuación 2 se determina la desviación estándar dado

el aporte de cada fuente.

i

i

NFPEFT

EFA 2)1

( (2)

La máxima desviación estándar se presenta en la ecuación

3.

22 )1

)(1()1

1()NFP

NFPNFP

Max(σ (3)

Siendo:

EFA: La energía firme anual de la generación de clase i

EFT: Energía firme total de la matriz eléctrica

NFP: Número de fuentes primarias existentes en la matriz

eléctrica

t : Periodo de análisis.

2.1.1 Energía Firme de una central se la entiende como la

mínima energía no interrumpible y garantizada de

acuerdo a la disponibilidad de la infraestructura y el

acceso al recurso para la producción de energía

eléctrica en todo momento.

La energía firme, presenta ciertas características

particulares según el tipo de generación. Para el caso de

las centrales hidroeléctricas se tiene principalmente dos

tipos de centrales, sin embalse de regulación también

conocidas como de “pasada” y con embalse de

regulación.

Para el caso de las centrales de pasada al no tener

embalse el cual sirva para el almacenamiento de agua,

situación que determina que se transforme la energía

potencial del agua en energía eléctrica inmediatamente;

lo cual determina que la energía firme corresponda a la

mínima producción del periodo de análisis.

Para el caso de las centrales con embalse es necesario

determinar el caudal mínimo histórico del periodo y

considerar el factor de planta de la misma a fin de

determinar la energía mínima del periodo en base a las

características físicas del embalse de la central.

Para las centrales térmicas el factor de planta del

periodo determina el consumo de combustible debido a

las condiciones operativas del sistema, las mismas que

son controlables y planificadas por el operador de la

central, es decir; la energía firme de este tipo de

tecnología depende de la potencia disponible en el

periodo en función de su factor de planta.

2.2 Fuentes autóctonas (FA): Este indicador determina la

proporción de la energía firme de fuentes locales respecto

del total de generación de la matriz eléctrica, como se

muestra en la ecuación 4.

i

i

EFT

EFAAFA (4)

Siendo:

EFAAi: La energía firme anual de la generación local de la

clase i

2.3 Energía firme del territorio nacional (EFTN): Este

indicador determina la energía firme, entendida como la

mínima energía posible de una central de generación para

la menor disponibilidad de la fuente del recurso de

producción energética, respecto al consumo de energía

eléctrica en un periodo de tiempo como se muestra en la

ecuación 5.

i

i

AnualEnergíaConsumo

FirmeEnergíaEFTN

__

_ (5)

2.4 Valor presente de los costos (IVPC): Este indicador

determina el costo total de la operación y del escenario de

inversión de los proyectos de generación realizada para

cada año considerando una tasa de descuento. Si el costo

es el menor posible permite a los consumidores acceder al

suministro de energía eléctrica, como se muestra en la

ecuación 6.

)( j

j

j

VPIVPCOMax

VPIVPCO

IVPC

(6)

Siendo:

VPCO: Valor presente de los costos de operación del periodo

de análisis.

VPCIj: Valor presente dela inversión del escenario j en el

periodo de análisis.

17



Para transferir una cantidad de dinero a valor presente es

necesario considerar el efecto de la inflación anual del periodo

en análisis.

2.5 Generación de valor agregado (IGVA): Este indicador

muestra el impacto del sector eléctrico respecto al valor

agregado nacional, como se presenta en la ecuación 7.

)(CVAGNSEMax

VPASEIGVA (7)

Siendo:

VPASE: Valor presente del agregado del sector eléctrico en

el periodo de análisis.

CVAGNSE: Valor agregado nacional codificado del sector

eléctrico

3. APLICACIÓN DE LA ROSA DE ROBUSTEZ

AL SISTEMA NACIONAL INTERCONECTADO

DEL ECUADOR

En el numeral anterior se presentó los indicadores que

conforman la rosa de robustez, mas la construcción de los

mismos requiere el procesamiento de información de acuerdo

a las características del sistema eléctrico, para el caso de un

sistema hidrotérmico. La información es obtenida de

documentos públicos de entidades del sector eléctrico y

económico del país.

Para el caso de la energía firme de la generación hidráulica se

determina la energía anual con probabilidad de excedencia del

90%, o el percentil del 10% de la fuente primaria de

producción de energía y a través de las funciones de

producción respectivas obtener la mínima energía que

garantiza anualmente suministrar la central.

Para el caso de la generación térmica se determina la potencia

efectiva considerando la tasa de salida forzada (FOR) que

disminuye la potencia nominal de la central o unidad, como se

indica en la ecuación 8, para las 8760 horas que posee un año.

8760*)1( FORPNEFGT (8)

Siendo:

EFGT : Energía firme de generación térmica.

PN : Potencia nominal de la central o unidad.

Para las centrales fotovoltaicas y eólicas disponibles en el

despacho centralizado del SNI al representar un porcentaje

reducido de la potencia instalada del sistema eléctrico, la

energía firme se la considera como si fuese una central

hidroeléctrica de pasada, es decir, como la menor energía

suministrada hacia el sistema eléctrico en el periodo de

análisis.

Dado el hecho de que se analiza la evolución de las inversiones

realizadas y los costos operativos de cada año es necesario

analizarlo al año más actual lo cual se utiliza la inflación para

determinar el valor presente de los montos de dinero.

Se calcula para cada año los diferentes indicadores con la

información de la disponibilidad de las unidades de generación

y energía producida, mediante lo cual se construye la rosa de

robustez para cada uno. Mientras que, para definir la rosa de

robustez del periodo total de análisis se considera el valor

promedio ponderado en base a la demanda del mismo.

La herramienta rosa de robustez se aplica al SNI del Ecuador,

considerando la operación de las centrales y unidades

disponibles obtenidas de los informes anuales del Operador

Nacional de Electricidad (CENACE) e inversiones y fechas de

los proyectos de generación realizadas en cada año en base a

información del operador del sistema eléctrico y de la agencia

de regulación y control eléctrico (ARCONEL).

A fin de contrastar la información económica se utilizó la

documentación pública del Ministerio de Electricidad y

Energía Renovable (MEER), Ministerio Coordinador de

Sectores Estratégicos (MICSE) y del Banco Central del

Ecuador (BCE).

Con la información de ARCONEL (n.d.) y CENACE (n.d.) en

la tabla 1, la energía firme anual de las centrales y unidades del

SNI donde se considera como energías renovables la central

eólica Villonaco de 15 MW y las centrales solares.

Tabla 1. Datos del SNI, Ecuador 2011-2015.

Año

DEMANDA

(GWh)

TÉRMICA

(GWh)

HIDRAÚLICA

(GWh)

ENERGÍA

RENOVABLE

(GWh)

2011 17747,80 3876,6 5156,88 0,12

2012 18605,91 4512,12 5145,48 22,8

2013 19458,95 3022,44 4335,36 209,88

2014 20882,55 2985,48 4658,52 406,08

2015 21934,39 3498,96 6301,44 594,96

Fuente: En base a datos de ARCONEL Y CENACE.

En la tabla 2, se presenta el consolidado de la inversión, costos

operativos y el valor agregado del suministro eléctrico para

cada uno de los años, así como el valor presente de los mismos

considerando la inflación de los años respectivos. Donde se

observa que en los años 2013 y 2014 debido a que los caudales

afluentes a las centrales hidroeléctricas disminuyeron los

costos operativos incrementaron debido a la operación de

centrales térmicas a fin de abastecer la demanda (MICSE, n.d.;

MEER, n.d.; BCE, n.d.; INEC, n.d.).

Tabla 2. Monto de la Inversión, Costos Operativos y valor Agregado del

Sector Eléctrico 2011-2015 Año

Inflación

(%)

Costo de Operación

e Inversión

(Miles de Dólares)

Valor Agregado de

Suministro Eléctrico

(Miles de Dólares)

2011 5,41 1 479 570,0 927 655,0

2012 4,16 1 747 370,0 1 046 322,0

2013 2,70 2 200 440,0 1 065 528,0

2014 3,67 2 258 970,0 1 301 923,0

2015 3,38 1 924 800,0 1 557 354,0

Fuente: En base a datos de BCE, INEC, MEER y MICSE.

18



En la tabla 3, se muestra los indicadores obtenidos mediante el

procesamiento de la información de los documentos oficiales,

para el caso de los valores monetarios los mismos se

encuentran trasladados al año 2015 por medio de la inflación.

Considerando que el SNI cuenta con siete fuentes primarias y

secundarias de energía; con la información de las tablas 1 y 2,

mediante la aplicación de las ecuaciones anteriores se

determina los indicadores de la rosa de robustez.

Tabla 3. Indicadores de la Rosa de Robustez del SNI, Ecuador 2011-2015

Año

IDF

(%)

EFTN

(%)

IVPC

(MM$)

FA

(%)

IGVA

(MM$)

2011 0,757 0,85 1267.42 0,442 805,13

2012 0,758 1,03 1582.44 0,387 947,58

2013 0,766 0,90 2048.04 0,392 991,76

2014 0,797 0,90 2182.62 0,389 1.257,92

2015 0,790 1,08 1924.8 0,450 1.557,35

Al normalizar los mismos respecto al máximo valor de cada

uno en el periodo de análisis; cuyos valores se muestran en la

tabla 4.

Tabla 4. Indicadores Normalizados Rosa de Robustez del SNI, Ecuador

2011-2015 Año

IDF

(%)

EFTN

(%)

IVPC

(%)

FA

(%)

IGVA

(%)

2011 0,95 0,79 0,59 0,98 0,52

2012 0,95 0,95 0,73 0,86 0,61

2013 0,96 0,83 0,94 0,87 0,64

2014 1,00 0,83 1,00 0,86 0,81

2015 0,99 1,00 0,88 1,00 1,00

La figura 2, presenta las rosas de robustez de las matrices

eléctricas de los sistemas ecuatoriano, uruguayo y chileno;

donde los indicadores IGVA y EFTN presenta un

comportamiento similar, mientras que los otros tres

indicadores son diferentes en las matrices eléctricas, las

siguientes figuras analizan de manera individual cada rosa de

robustez de cada sistema eléctrico.

Figura 2. Comparación de Rosa de Robustez.

En la figura 3, se presenta la rosa de robustez obtenida para el

SNI para el periodo 2011-2015, donde se observa que la misma

muestra una matriz eléctrica que cuenta con una alta diversidad

de fuentes ya que en el sistema existen dos tipos de centrales

hidráulicas (embalse y pasada); diferente tipo de tecnología de

centrales térmicas como son residuo, diésel, nafta; centrales en

base a fuentes de energía renovables e interconexiones

eléctricas. Debido a la alta participación de la generación

hidráulica del sistema los costos de energía son bajos, dado

que en el periodo ha existido una buena hidrología. Si bien en

el periodo la energía firme se ha incrementado; es necesario

tener presente que la misma se debe a la participación de

centrales hidráulicas de pasada, centrales eólicas y solares,

cuya energía depende de la disponibilidad estocástica del

recurso renovable (agua, viento y radiación solar

respectivamente). En el periodo el indicador presenta un alto

valor, es debido a la reducida demanda frente a la potencia

instalada existente, para el caso del indicador de fuentes

autóctonas este se ha incrementado debido al bajo despacho en

el periodo de la interconexión eléctrica internacional e

importación de diésel, lo cual es función de la hidrología que

depende del nivel de lluvias en los sitios donde están ubicados

las centrales hidroeléctricas.

A fin de mostrar la característica de análisis de esta

herramienta para determinar de manera gráfica la robustez de

una matriz eléctrica, en la figura 3 se muestra la rosa de

robustez para el sistema eléctrico del Uruguay establecida en

(DNETN, 2007) para el periodo 2007-2025 para la expansión

del sistema en base a centrales de gas natural, donde se ve que

una mejor robustez del sistema ecuatoriano debido a que la

rosa de robustez está en el perímetro externo, en especial los

indicadores de fuentes autóctonas y diversidad de fuentes del

territorio nacional, mientras que no existe diversidad de

fuentes. Los indicadores económicos del sistema eléctrico

uruguayo muestran un mejor desempeño lo cual es una medida

de que la expansión del sistema considero un mejor

aprovechamiento para el financiamiento de los proyectos.

Figura 3. Rosa de Robustez del SNI, Ecuador 2011-2015.

En la figura 4, se presenta la rosa de robustez para el sistema

eléctrico chileno establecida en Molina (2005) para 2007-2025

en base a la planificación de la expansión del sistema eléctrico

chileno, donde se presenta una rosa de robustez de

comportamiento muy semejante a la de la matriz eléctrica

ecuatoriana, en la cual dos de los cinco indicadores son

aproximadamente iguales, fuentes autóctonas y energía firme

19



del territorio nacional. Mientras que, para los indicadores

económicos presenta un mejor desempeño el sistema eléctrico

ecuatoriano lo que da una medida de que la expansión del

sistema considero un mejor aprovechamiento del

financiamiento de los proyectos.

Figura 4. Rosa de Robustez del Sistema Eléctrico Uruguayo 2007-2025.

El análisis de la rosa de robustez permite establecer

lineamientos para la expansión de la generación y para el caso

del sistema eléctrico ecuatoriano es posible determinar los

siguientes lineamientos:

- Las instituciones estatales deben invertir recursos a fin de

contar con un catálogo de proyectos de la cadena de

producción eléctrica con estudios actualizados y con una

base de datos de información pública sobre el potencial de

los recursos energéticos con fines de generación eléctrica

existentes en el país.

- Mantener y continuar en la senda de la diversificación de

fuentes de generación eléctrica, se recomienda incentivar

la expansión del parque generador mediante el uso energías

renovables de manera masiva por medio de generación

distribuida.

- Optimizar la planificación de expansión del sector eléctrico

considerando el uso eficiente de recursos energéticos y

económicos, a fin de incrementar el índice IDF.

- Revisar los esquemas de comercialización local e

internacional de energía eléctrica: analizar esquemas de

contrato a largo plazo, mediante lo cual puede mejorar el

valor agregado de cada kWh de energía eléctrica.

Figura 5. Rosa de Robustez del Sistema Eléctrico Chileno 2007-2025.

4. CONCLUSIONES

Se mostró que la herramienta Rosa de Robustez permite

determinar de manera gráfica la robustez del sistema eléctrico

respecto a variables técnicas y económicas; así como comparar

esta característica entre sistemas colocándole en una base que

no depende del nivel de consumo o de oferta de generación.

Mediante el conjunto de indicadores, estos permiten evaluar el

grado de robustez alcanzado en un periodo por un país y

mediante el cual es posible plantear lineamientos, políticas

para la planificación de la expansión del mismo, así como

determinar metas para el cumplimiento de las mismas.

La seguridad de suministro en el abastecimiento de la energía

eléctrica debe buscar un equilibrio entre las fuentes primarias

en base a recursos estocásticos y el nivel de energía térmica

dada por fuentes autóctonas que demandan la producción de

combustibles en refinerías locales. Dada la dinámica de la

demanda se requiere que el sistema lleve una continua

planificación, desarrollo y construcción de proyectos

eléctricos de bajo costo operativo y de un alto impacto en el

valor agregado del país.

Para trabajos futuros se recomienda colocar indicadores del

nivel de CO2 emitido por la operación de la matriz eléctrica;

así como incluir la elaboración de escenarios para analizar la

expansión en un tipo o mix de tecnologías de producción de

generación eléctrica.

REFERENCIAS

Agencia Nacional de Regulación y Control de Electricidad (n.d.).

Información estadística del sector eléctrico ecuatoriano, Obtenido de:

http://www.arconel.gob.ec (Mayo, 2017).

Banco Central del Ecuador (n.d.). Información de la Política Económica Obtenido de: http://www.bce.fin.ec (Mayo, 2017).

Blanco, N. (2015). Análisis de seguridad y productividad del suministro de

energía eléctrica en el sistema eléctrico de Nicaragua en el periodo comprendido desde el año 2010 hasta el 2018, Revista Iberoamericana

de Bioeconomía y Cambio Climático, Vol.1, Nº 2, pp. 20-53.

DNETN (2007). Robustez del Sistema Eléctrico Nacional: Aporte Metodológico y Ejercicio de Aplicación, Nota Técnica, Dirección

Nacional de Energía y Tecnología Nuclear, Ministerio de Industria,

Energía y Minería Uruguay.

Instituto Nacional de Estadísticas y Censos (n.d.). Información estadística del

país, Obtenido de: http://www.ecuadorencifras.gob.ec (Mayo, 2017).

Ministerio Coordinador de Sectores Estratégicos (n.d.). Información Energética de Sectores Estratégicos, Obtenido de:

http://www.sectoresestrategicos.gob.ec (Mayo, 2017).

Ministerio de Electricidad y Energía Renovables (n.d.). Información del sector eléctrico, Obtenido de:http://www.energia.gob.ec (Mayo, 2017).

Molina, J., Martínez, V., y Rudnick H. (2005). Indicadores de Seguridad

Energética: Aplicación al Sector Energético de Chile, PUC-Chile. Mosto, P. (2008). Robustez de Sistemas Eléctricos: Conclusiones del Aporte

Metodológico y ejercicio de aplicación para el Sistema Uruguayo,

Revista Energía - CEDDET, Nº 2, pp. 14-17. Operador Nacional de Electricidad (n.d.). Información sobre la operación del

SNI, Obtenido de: http://www.cenace.org.ec (Mayo, 2017).

Retamales, G. (2005). Indicadores de Seguridad de Suministro Eléctrico (SSE) en Chile, Tesis Universidad de Chile.

0.0

0.2

0.4

0.6

0.8

1.01

2

34

5

Uruguay

0.0

0.2

0.4

0.6

0.8

1.0IDF

FA

EFTNIVPC

IGVA

Chile

20

http://www.arconel.gob.ec/

http://www.bce.fin.ec/

http://www.ecuadorencifras.gob.ec/

http://www.sectoresestrategicos.gob.ec/

http://www.energia.gob.ec/

http://www.cenace.org.ec/



BIOGRAFÍAS

Oscullo Lala José. Nació en

Sangolquí, Ecuador, en 1971. Recibió

su título de ingeniero eléctrico en la

Escuela Politécnica Nacional en

1996, Master en ingeniería eléctrica

de la Universidad Estatal de

Campinas, Sao Paulo Brasil en 2002.

Actualmente se desempeña como

profesor titular del Departamento de Energía Eléctrica de la

Facultad de Ingeniería Eléctrica y Electrónica de la Escuela

Politécnica Nacional. Su campo de investigación se encuentra

relacionado a la aplicación de herramientas de inteligencia

computacional para la estabilidad de sistemas eléctricos de

potencia.

Romero Luis. Nació en Quito, Ecuador,

en 1990. Recibió el título de ingeniero

eléctrico en la Escuela Politécnica

Nacional en 2018. Actualmente se

encuentra ejerciendo el libre ejercicio

profesional en el campo de la eficiencia

energética. Sus áreas de interés son:

Economía de la Energía, Eficiencia

Energética y Optimización de la planificación de largo,

mediano y corto plazo de sistemas hidrotérmicos.

21


Revista Politécnica – JULIO 2018, Vol. 41, No. 2

*[email protected] Recibido: 01/12/2016 Aceptado: 20/07/2018 Publicado: 31/07/2018

11. INTRODUCCIÓN

NGC 7331 es una galaxia espiral cercana, con un tamaño aparente de 13.5 minutos de arco (Bosma, 1978) y situada a 14,7 Mpc de la Tierra (de Blok et al., 2008). La normal del plano de esta galaxia posee una inclinación (i) de 74.8 grados con respecto a la línea de mira del observador (de Blok et al., 2008). Gráficas de la velocidad de rotación del gas de hidrógeno atómico representadas en función del radio galáctico se han usado en un estudio del campo de velocidades de NGC 7331 (Bosma, 1981). Este tipo de gráficas se conocen como curvas de rotación. La masa de la materia visible (gas atómico y molecular, polvo cósmico, estrellas, etc.) observada en NGC 7331 y en otras galaxias deberían dar lugar a curvas de rotación con velocidades que disminuyen con el aumento del radio galáctico, contrario a lo que revelan las observaciones de NGC 7331 (Bosma, 1981) y de otras galaxias (Bosma, 1981; Sofue y Rubin, 2001). Un estudio llevado a cabo por de Blok et al., (2008) revela que las curvas de rotación de NGC 7331 obtenidas a partir del gas HI de las regiones que se acercan y se alejan del observador difieren en alrededor 25 km s-1. La curva de rotación en NGC 7331 se mantiene relativamente constante en regiones externas situadas a gran distancia de su centro (Rubin, 1965; Bosma, 1981). Este hecho es característico de las galaxias espirales

(Sofue et al., 1999). Se ha propuesto la presencia de un halo de materia oscura en las galaxia espirales, materia responsable de mantener constante la velocidad rotacional en regiones tan alejadas de los centros galácticos. Hasta el momento se desconoce la naturaleza de la materia oscura, pues ésta se puede estudiar únicamente por sus efectos gravitatorios. Sofue et al., (1999) estudiaron las curvas rotacionales de una muestra de galaxias espirales, hallando que las galaxias barradas muestran mayor dispersión de velocidad que aquellas de galaxias “normales”, ello posiblemente debido a movimientos no circulares del gas. Estos autores también descubrieron que estallidos de formación estelar nuclear y núcleos activos parecen no correlacionar con las propiedades de las curvas de rotación, sugiriendo así que estos procesos nucleares son causados por efectos locales y no por propiedades dinámicas globales. Teniendo en cuenta una aproximación esférica, la masa (M) de una galaxia se relaciona con la velocidad de rotación (Vrot) mediante la siguiente expresión (Sofue y Rubin, 2001):

𝑀 =R

(1)


Armijos-Abendaño Jairo1*; López Ericson1; Llerena Mario1; Aldas Franklin1

1Escuela Politécnica Nacional, Observatorio Astronómico, Quito, Ecuador

Resumen: Se llevó a cabo un estudio de la curva de rotación de la galaxia NGC 7331 situada a 14,7 Mpc de la Tierra. La curva de rotación se derivó usando observaciones radioastronómicas del monóxido de carbono (CO). La forma de la curva de rotación y las velocidades son muy similares a aquellas derivadas previamente, empleando datos de hidrógeno atómico, lo que sugiere la coexistencia de ambos elementos en las regiones estudiadas en NGC 7331. Se descubrió que el ensanchamiento de la línea de CO, transición 2-1, está dominado por efectos turbulentos del gas de CO antes que por efectos térmicos. Asimismo, se estudió el campo de velocidades de NGC 7331, lo que puso en evidencia la rotación de la galaxia en el sentido de las agujas del reloj. Finalmente, asumiendo una aproximación esférica para la forma de la galaxia se estimó una masa dinámica de 1,4E+10 masas solares para NGC 7331. Palabras clave: Radioastronomía, galaxia, curva de rotación, masa galáctica, gas molecular.

Kinematics and Mass of the NGC 7331 Galaxy

Abstract: The study of the rotation curve of the galaxy NGC 7331, located at 13,870 Mpc from the Earth, was carried out. The rotation curve is obtained using radio astronomy observations of the carbon monoxide (CO). The shape of the rotation curve and the velocities are very similar to those derived from atomic hydrogen data, suggesting the coexistence of both elements in the studied regions of NGC 7331. We found that the CO (2-1) line broadening is dominated by turbulent effects rather than by thermal effects. The velocity field of the galaxy is studied as well, showing that the galaxy is rotating clockwise. Finally, assuming a spherical approximation for the galaxy shape, a dynamical mass of 1,4E+10 solar masses is estimated for NGC 7331. Keywords: Radio astronomy, galaxy, rotation curve, galactic mass, molecular gas.

23



donde R es el radio de la galaxia y G es la constante gravitacional e igual a 6,67E-11 m3 kg-1 s-2. Uno de los objetivos del presente trabajo consiste en la estimación de la masa de la galaxia NGC 7331, ello a partir de su curva de rotación derivada a partir de datos de monóxido de carbono CO (2-1) (Leroy et al., 2009), donde (2-1) denota el hecho que se observó la transición rotacional 2-1 a 230,5 GHz. Estos datos se describen brevemente en la Sección 2. En la Sección 3 se muestra un mapa y espectros de CO (2-1) de NGC 7331. En la Sección 3.1 se lleva a cabo un estudio sobre el ensanchamiento de la línea de emisión de CO (2-1). La Sección 4 ofrece una explicación sobre la metodología empleada en la obtención de una curva de rotación, en base a los datos de monóxido de carbono. Un estudio sobre el campo de velocidades en la galaxia NGC 7331 se presenta en la Sección 5. El procedimiento usado en la derivación de la masa dinámica de la galaxia NGC 7331 se expone en la Sección 6. Finalmente, en la Sección 7 se presentan las conclusiones del presente trabajo.

2. DATOS

Para llevar a cabo el presente estudio se emplearon datos públicos de CO (2-1), los cuales se publicaron por primera vez en Leroy et al., (2009). Estos datos se obtuvieron con el radio telescopio IRAM de 30 metros de diámetro (http://www.iram-institute.org/EN/30-meter-telescope.php). A la frecuencia (230,5 GHz) de la transición de CO (2-1) el radiotelescopio proporcionó una resolución espacial de 13 segundos de arco. El autocorrelador WILMA se usó en las observaciones del radio telescopio IRAM, produciendo espectros con una resolución de 2,6 km s-1 y un ancho de banda de 1200 km s-1 a la frecuencia de la transición CO (2-1) (Leroy et al., 2009). En la observación de la galaxia NGC 7331 y de otras galaxias con el radio telescopio IRAM se usó el modo de observación On the Fly. La información de las observaciones se guardó en cubos de datos, donde las intensidades de los espectros de CO (2-1) están dadas en temperatura de antena y los espectros suavizados a una resolución de 10,4 km s-1. Los datos que se descargaron para nuestro trabajo ya están corregidos por la línea de base. Para el desarrollo del presente trabajo se descargó el cubo de datos de CO (2-1) correspondiente a la galaxia NGC 7331 que está disponible en la página web de IRAM (http://www.iram-institute.org/). 3. MAPA Y ESPECTROS DE LA GALAXIA NGC 7331

Se usó el software GILDAS (https://www.iram.fr/IRAMFR/GILDAS/) para la visualización del cubo de datos de NGC 7331. Seguidamente con este software se procedió a ejecutar una tarea que permitió obtener un mapa integrado de NGC 733 entre +550 y +1100 km s-1. Este mapa se muestra en la Figura 1. En esta figura los ejes X y Y corresponden a la ascensión recta (α) y la declinación (δ), respectivamente, en el sistema ecuatorial de coordenadas J2000. La cruz roja en la Figura 1 representa el

centro radio (fotométrico) de NGC 7331 con las siguientes coordenadas: 22h37m04,0s y 34º24’56.4’’ (Israel y Baas, 1999). En la Figura 1 también se ha dibujado una elipse de color azul con un eje mayor igual a 290 segundos de arco y el eje menor igual a 84 segundos de arco, y un ángulo de posición de 167 grados (Israel y Baas, 1999). El tamaño de la elipse ha sido elegido a ojo y aproxima bastante bien la forma de la galaxia. El centro de la elipse es el centro radio de NGC 7331. En esta figura se usan, asimismo, contornos que permiten resaltar la forma de la galaxia. Estos contornos inician a 3σ y continúan con pasos 10σ, con σ siendo igual a 0,33 Kelvin (K) km s-1. Por otro lado, en la Figura 2 se presenta el espectro de emisión de CO (2-1) que corresponde al centro radio de la galaxia NGC 7331, indicado con una cruz roja en la Figura 1. Como se puede apreciar en la Figura 2, la línea de emisión de CO (2-1) tiene una forma relativamente gaussiana y un pico de emisión a una velocidad de ~831 km s-1. En la Figura 3 se presenta el espectro que corresponde a una región alejada alrededor de 60 segundos de arco del centro de la galaxia NGC 7331. En este espectro no se detecta emisión de CO (2-1) a diferencia de aquella detectada en el centro de NGC 7331. Otro hecho importante que se infiere a partir de las Figuras 2 y 3 es que el ruido de los espectros es de aproximadamente 10 mK.

Figura 1. Mapa de intensidad integrada de CO (2-1) de la galaxia NGC 7331. En la figura se usó el sistema de coordenadas ecuatorial J2000. La cruz roja señala el centro radio de la galaxia, mientras que la cruz azul representa

la posición del espectro que se muestra en la Figura 3. Toda la región de emisión de la galaxia puede aproximarse por una elipse de color rojo. La

línea azul es el eje mayor de la elipse. La paleta de grises muestra la intensidad integrada del mapa dada en K km s-1. En esta figura se usan, además, contornos de intensidad que resaltan la forma de la galaxia. La

elipse rellena en rojo muestra el tamaño del haz del telescopio de 13 segundos de arco.

Se realizó, asimismo, un ajuste gaussiano a la línea de CO(2-1) mostrada en la Figura 2, con el objetivo de derivar los valores del pico de intensidad (T), la velocidad radial (V) y el ancho de la línea (ΔV). De esta manera se obtuvieron los siguientes valores: T=99,0 mK, V=831,2±2,3 km s-1 y ΔV=142,3±5,9 km s-1. El valor de ΔV permitirá llevar a cabo un breve estudio (Sección 3.1), sobre las contribuciones de los

24



efectos térmicos y turbulentos del gas molecular en el ensanchamiento de la línea de CO (2-1). La velocidad de 831 km s-1 se considera como la velocidad sistémica de NGC 7331 y ésta es similar a aquella dada en Israel y Baas (1999) y mayor que aquella de 818 km s-1 dada en de Blok et al., (2008).

Figura 2. Espectro de emisión de CO (2-1) de la galaxia NGC 7331. El

espectro corresponde al centro radio de NGC 7331, indicada en la Figura 1. La velocidad sistémica de NGC 7331 se considera igual a 831 km s-1, la que

se señala con una línea roja entrecortada. En el eje Y la intensidad de la radiación está dada en Kelvin. La línea de color azul muestra el ajuste

gaussiano a la línea de CO (2-1) (Véase la Sección 3).

Figura 3. Espectro observado en el mismo rango espectral, centrado a

230,538 GHz, al igual que el de la Figura 2. La posición, alejada alrededor de 60 segundos de arco del centro de NGC 7331, a la que corresponde este espectro, se indica con una cruz azul en la Figura 1. En este espectro no se

detecta emisión de CO (2-1) a diferencia de lo que sí se observa en la Figura 2.

3.1. Ensanchamiento De La Línea De CO (2-1)

Si se supone que la anchura de una línea molecular es consecuencia solamente de efectos térmicos y turbulentos, entonces el ancho de la línea ΔV de monóxido de carbono viene dado por la siguiente ecuación (Teague et al., 2016):

𝛥𝑉tur = ln− 𝛥𝑉 (2)

y

𝛥𝑉 =kTcin (3)

donde ΔVtur es el ensanchamiento de la línea molecular originada por turbulencia, ΔVter es el ensanchamiento térmico de la línea molecular, k=1,38E-23 J K-1 es la constante de Boltzmann, Tcin es la temperatura cinética, μ=28 es la masa molecular de CO, y mH es la masa molecular del hidrógeno atómico e igual a 1,67E-24 g. Si se reemplazan estos valores en las Ecuaciones (2) y (3) y considerando una Tcin=20 K (Israel y Baas, 1999), se halla una ΔVter=0,11 km s-1 y una ΔVtur=85,0 km s-1. Este resultado pone en evidencia que en NGC 7331 domina la turbulencia en el ensanchamiento de la línea de CO (2-1) ante los efectos térmicos. Sin embargo, el ensanchamiento de las líneas atómicas y moleculares no depende únicamente de los efectos térmicos y turbulentos, sino también que puede estar afectado, en cierto grado, por los movimientos sistémicos y relativos de unas partes de la región estudiada con respecto a otras. El tamaño del haz del telescopio IRAM de 13 segundos de arco a 230,5 GHz traza regiones de alrededor 900 pc de NGC 7331. Esto puede causar que movimientos relativos del gas afecte a la anchura de la línea. Sin embargo, incluso suponiendo que la anchura no afectada por movimientos relativos sea cuatro veces menor que lo que medimos, es decir un ancho de la línea de ~36 km s-1, se obtendría una ΔVtur=21.3 km s-1, que sigue siendo mucho mayor que aquella causada por los efectos térmicos. La resolución espectral de los datos usados en este trabajo de 10 km s-1 empezaría a afectar a la anchura de las líneas de CO si éstas tuvieran anchuras menores a 60 km s-1.

4. CURVA DE ROTACIÓN DE NGC 7331 En esta sección se procederá a derivar una curva de rotación de la galaxia NGC 7331, estudio previo a la determinación de la masa dinámica de NGC 7331. Dicha curva a su vez se obtiene a partir de una figura donde se representan las velocidades de la emisión de CO (2-1) correspondientes a diferentes posiciones. Para este estudio se han obtenido las velocidades de las posiciones situadas a lo largo del eje mayor de NGC 7331 mostrada en la Figura 1. Como se mencionó este eje tiene un ángulo de posición de 167 grados. Los resultados se muestran en la Figura 4, donde el centro de las posiciones corresponde al centro radio de NGC 7331. La figura resultante se conoce como diagrama posición-velocidad, donde se usan también contornos que ponen de manifiesto de una mejor manera la forma relativamente simétrica de las velocidades del gas con respecto a la velocidad sistémica de 831 km s-1. Dichos contornos inician a 3σ y continúan con pasos 4σ, con σ siendo igual a 10 mK. En la Figura 4 también se indica las curvas de rotación (en color rojo) que viene determinadas por la siguiente ecuación (Sofue y Rubin, 2001):

𝑉rot = 𝑉 − 𝑉sys /sin𝑖 − (𝜎 + 𝜎 ) / (4)

25



donde Vt es la velocidad terminal, Vsys es la velocidad sistémica de la galaxia con respecto al observador, sin i es un factor de corrección relacionado con la inclinación i de la galaxia (i=75,8 grados para NGC 7331), σISM es la dispersión de velocidades y σobs es la resolución en velocidades dada por las observaciones. Vt se define como la velocidad que corresponde a un 20% del pico de intensidad en una posición dada (Sofue y Rubin, 2001). En este trabajo se considera un valor de σISM=7 km s-1, que es un valor típico para el medio interestelar de nuestra Galaxia (Sofue y Rubin, 2001), entretanto, de las observaciones se conoce que σobs=10,4 km s-

1. De esta manera se obtienen las curvas de rotación indicadas con color rojo en la Figura 4. Promediando velocidades de rotación al rojo y al azul se obtienen los valores, a un determinado radio, que se listan en la Tabla 1.

Tabla 1. Velocidades de rotación derivadas para NGC 7331 Radio Vrot

(segundos de arco) (km s-1)

0 0 13,4 89,7 26,8 229,0 40,0 254,8 53,6 265,1 67,0 265,1 80,4 262,7 93,8 266,1

107,2 273,3 120,6 268,2 130,0 265,1 145,0 265,0

Figura 4. Diagrama posición-velocidad.

Las posiciones en el eje X inferior, dadas en segundos de arco, corresponden al semieje mayor indicado con la línea azul en la Figura 1. En la misma figura, el centro de posiciones, que corresponde al centro radio de NGC 7331, es indicado con la cruz roja. El eje X superior está dado en kpc. En la presente figura, la línea roja representa la curva de rotación que se obtiene siguiendo los pasos descritos en la Sección 4. En la figura se usan contornos que resaltan la estructura relativamente simétrica de la emisión de CO(2-1) con respecto a la velocidad sistémica de la galaxia de 831 km s-1, la misma que se indicada con un línea entrecortada.

Se aprecia en el diagrama posición-velocidad que la emisión de CO (2-1) llega hasta un radio de ~145 segundos de arco, valor que equivale a 10 kpc considerando una distancia para la fuente de 14,7 Mpc (de Blok et al., 2008). A esta distancia del centro de la galaxia NGC 7331, se infiere una velocidad de rotación de 265 km s-1. La forma de la curva de rotación y las velocidades derivadas en este trabajo son muy similares a aquellas derivadas en base a las observaciones de la línea de hidrógeno atómico (Bosma, 1981), sugiriendo la coexistencia del gas de monóxido de carbono e hidrógeno atómico dentro de un radio de 10 kpc en NGC 7331. Sin embargo, la emisión producida por el hidrógeno atómico es mucho más extensa que aquella de CO(2-1), debido a que en el trabajo llevado a cabo por Bosma (1981) se muestra la detección de la emisión del hidrogeno atómico en regiones alejadas 30 kpc del centro de NGC 7331.

5. CAMPO DE VELOCIDADES EN NGC 7331 En la Figura 5 se muestra un mapa de la velocidad promedio del gas de CO, ello con el objetivo de estudiar el campo de velocidades en la galaxia NGC 7331. Este mapa se ha obtenido usando una tarea del software GILDAS, donde se tiene que especificar el rango de velocidades donde se buscará la velocidad promedio en el cubo de datos descrito en la Sección 2. Como se observa en la Figura 5, existe una diferencia marcada entre las velocidades del gas de CO, para las regiones sur y norte de NGC 7331, cuya forma ha sido delineada con una elipse similar a aquella mostrada en la Figura 1. Las velocidades más altas, mayores a 831 km s-1, se encuentran en las regiones sur de la galaxia, mientras que las velocidades menores a 831 km s-1 se hallan en las regiones norte de NGC 7331. La Figura 5 proporciona información sobre la rotación de la galaxia; el gas de monóxido de carbono que posee velocidades mayores que la velocidad sistémica de la galaxia, indica que este gas se aleja del observador, entretanto que el gas molecular con velocidades menores a 831 km s-1, se acerca al observador, es decir, la galaxia NGC 7331 se encuentra rotando en el sentido de las agujas del reloj. Se puede usar el campo de velocidades presentado en la Figura 5 para obtener una curva de rotación. Usando el software Gildas se puede extraer información sobre las velocidades a lo largo del eje mayor mostrado en la Figura 5. Las velocidades extraídas de esta manera para el gas que se acerca al observador se muestran con una línea de color verde en la Figura 5. Como se puede apreciar en esta figura las velocidades de rotación obtenidas a partir del campo de velocidades son menores (en la mayoría de las posiciones) a aquellas estimadas a partir del método expuesto en la Sección 3. Este hecho es consecuencia de que el método usado en la Sección 3 considera la velocidad terminal y no la velocidad promedio, lo que no subestima la velocidad de rotación.

26



Figura 5. Campo de velocidades de la galaxia NGC 7331.

La cruz roja muestra el centro de la galaxia, mientras que la elipse es similar a aquella mostrada en la Figura 1. La línea azul indica el eje mayor de la elipse. La paleta de gris señala los valores de velocidades mostrados en esta figura. La elipse rellena en rojo muestra el tamaño del haz del telescopio de 13 segundos de arco.

6. MASA DINÁMICA DE NGC 7331

Asumiendo una aproximación esférica para la forma de la galaxia NGC7331, se puede estimar su masa dinámica mediante la ecuación (1). El valor del semieje mayor de 145 segundos de arco de la elipse, que aproxima bastante bien la forma de NGC 7331 (véase la Sección 3), se considerará como radio de la esfera. Este valor equivale a 10.3 kpc. En la Tabla 1 se halla que una velocidad de rotación de 265 km s-1 corresponde a un radio de 145 segundos de arco. Reemplazando valores en la ecuación (1), se obtiene una masa dinámica de 1,7E+11 masas solares para NGC 7331. Asumiendo una densidad constante de la materia, se halla que la masa calculada para una geometría esférica, es un factor 12 mayor que aquella de un elipsoide con dos de sus semiejes iguales a 42 segundos de arco y el tercer semieje igual a 145 segundos de arco. Como se dijo, una elipse con semiejes de 42 y 145 segundos se empleó para ajustar la forma de NGC 7331 en la Figura 1. Considerando ese factor 12, se encuentra una masa dinámica de 1,4E+10 masas solares para NGC 7331. Es importante recalcar que esta masa dinámica excluye la masa de las regiones externas galácticas que en esta galaxia llegan hasta radios de alrededor 30 kpc (Bosma, 1981). La línea de HI de 21 cm traza bien las regiones externas de galaxias (Sandstrom et al., 2013), donde se piensa que se halla la mayor cantidad materia oscura. También se halla una masa igual a 1.1E+11 masas solares para NGC 7331 a partir de la intensidad integrada de CO(2-1) y factores de conversión entre la masa de HII y las masas estelar y de HI dados en la Tabla 3 de Leroy et al., (2009). De esta manera la masa de 1.1E+11 masas solares, considera las masas

de HII, HI y de estrellas. Los valores dados en la Tabla 3 de Leroy et al., (2009) se calcularon para una elipse con un tamaño sólo un factor 1,2 mayor que la elipse usada en nuestros cálculos. La masa hallada a partir de la curva de rotación bajo la suposición de una densidad constante de la materia en NGC 7331 subestima la masa dinámica en un factor 7,9. La estimación de masas dinámicas a partir de curvas de rotación y una aproximación esférica puede ser más correcta para galaxias con formas proyectadas relativamente circulares. En estos casos las aproximaciones usadas en este trabajo pueden dar valores más cercanos a aquellos que se basan en la intensidad integrada de CO, HI y en la emisión del polvo cósmico.

7. CONCLUSIONES

El ensanchamiento de las líneas de CO (2-1), en la galaxia NGC 7331, está dominado por movimientos turbulentos del gas antes que por efectos térmicos. La forma de la curva de rotación y las velocidades de NGC 7331 derivadas a partir de observaciones de la transición 2-1 del monóxido de carbono son muy similares a aquellas derivadas usando el hidrógeno atómico (Bosma, 1981). Este resultado sugiere la coexistencia del gas de monóxido de carbono e hidrógeno atómico en las regiones estudiadas de NGC 7331. Se estudió el campo de velocidades del gas de monóxido de carbono en la galaxia NGC 7331, trabajo que muestra que esta galaxia rota en el sentido de las agujas del reloj. Usando una aproximación esférica para la forma de la galaxia, se deriva una masa dinámica de 1,4E+10 masas solares para la galaxia NGC 7331, subestimando en un factor 7,9 el valor obtenido empleando un método basado en la intensidad integrada de la línea de CO (2-1). El método usado en este trabajo para derivar la masa dinámica de una galaxia puede ser más correcta para galaxias con formas relativamente circulares.

REFERENCIAS

Bosma, A. (1978). The distribution and kinematics of neutral hydrogen in

spiral galaxies of various morphological types. Obtenido de: https://ned.ipac.caltech.edu/level5/March05/Bosma/Bosma_contents.html. (octubre, 2016).

Bosma, A. (1981). 21-cm line studies of spiral galaxies. I. Observations of the galaxies NGC 5033, 3198, 5055, 2841, and 7331. The Astronomical Journal, vol. 86, pp. 1791-1846.

de Blok, W. J. G., Walter, F., Brinks, E., Trachternach, C., Oh S. -H., Kennicutt Jr., R. C. (2008). High-resolution rotation curves and galaxy mass models from THINGS. The Astronomical Journal, vol. 136, pp. 2648.

Leroy, A. K., Walter, F., Bigiel, F., et al., (2009). Heracles: the HERA CO line extragalactic survey. The Astronomical Journal, vol. 137, pp. 4670-4696.

Israel, F. P., & Baas, F. (1999). Molecular gas in the bulge and ring of NGC 7331. Astronomy and Astrophysics, vol. 351, pp. 10-20.

Rubin, V. C., Burbidge, E. M., Burbidge, G. R., Crampin, D. J., & Predergast, K. H. (1965). The rotation and mass of NGC 7331. Astrophysical Journal, vol. 141, pp. 759.

Sandstrom, K. M., Leroy, A. K., Walter, F., et al. (2013). The CO-to-H2 conversion factor and dust-to-gas ratio on kiloparsec scales in nearby galaxies. The Astrophysical journal, vol. 777, pp. 5

27



Sofue, Y., Tutui, Y., Honma, M., Tomita, A., Takamiya, T., Koda, J., Takeda, Y. (1999). Central rotation curves of spiral galaxies. The Astrophysical Journal, vol. 523, pp. 136-146.

Sofue, Y., & Rubin, V. (2001). Rotation curves of spiral galaxies. Annual Review of Astronomy and Astrophysics, vol. 39, pp. 137-174.

Teague, R., Guilloteau, S., Semenov, D., et al., (2016). Measuring turbulence in TW Hydrae with ALMA: methods and limitations. Astronomy and Astrophysics, vol. 592, pp. A49.

BIOGRAFÍAS

Ericson López. Doctor en Astrofísica (PhD) por la Academia de Ciencias de Rusia y Físico Teórico por la Escuela Politécnica Nacional. Ha realizado investigaciones post doctorales en

Brasil y Estados Unidos. Es científico colaborador del Harvard-Smithsonian Center para la Astrofísica y profesor adjunto del Departamento de Astronomía de la Universidad de Sao Paulo. Ha realizado más de treinta publicaciones científicas y varias otras publicaciones relevantes. Es Director del Observatorio Astronómico de Quito desde 1997 y miembro de la Academia de Ciencias del Ecuador. Es profesor principal de la Facultad de Ciencias de la EPN por más de 25 años, en la que imparte cursos formales de Física y Astrofísica.

Jairo Armijos. Realizó sus estudios de pregrado en Astronomía en la Universidad Estatal de San Petersburgo, Rusia. Obtuvo su título de Maestría en Astrofísica y su título de Doctorado en Astrofísica en la Universidad Autónoma de Madrid, España. Actualmente es investigador del Observatorio Astronómico de Quito. Mario Llerena. Realizó sus estudios de pregrado en la Escuela Politécnica Nacional, obteniendo su título de Físico en el 2016. Es miembro del Observatorio Astronómico de Quito como parte de Unidades Científicas de investigación en Gravitación y Cosmología,

Radioastronomía y Astrofísica de Altas Energías. Franklin Aldás. Realizó sus estudios de pregrado en la Escuela Politécnica Nacional, obteniendo su título de Físico en el 2018. Es miembro del Observatorio Astronómico de Quito como parte de Unidades Científicas de investigación en Gravitación y Cosmología, Radioastronomía y Clima Espacial.

28

Modelo Matemático Adaptado para el Cálculo de Pérdidas de Propagación en la Banda de 900 MHz para Microceldas en la Ciudad de Quito





11. INTRODUCCIÓN

Para el dimensionamiento de un sistema de radiocomunicación

se pueden utilizar diferentes modelos de propagación, ya sean,

empíricos, teóricos o semi-empíricos.

Los modelos de propagación se utilizan para predecir las

pérdidas medias en la trayectoria en función de variables, tales

como altura de las antenas transmisora y receptora, frecuencia

y distancia. La predicción de las pérdidas de trayectoria es muy

importante cuando se trata de determinar la cobertura de redes

inalámbricas, ya que en este medio es donde se presenta mayor

cantidad de efectos negativos que provocan atenuación de la

señal transmitida (He R. 2010). El modelo de propagación log-

distancia es muy popular cuando se trata de modelación de las

pérdidas de trayectoria (Rappaport, 2002). En este modelo, el

exponente de pérdida de trayectoria depende del entorno de

propagación de las ondas de radio y se lo puede determinar en

base a los datos de mediciones realizadas, dicho exponente

indica que tan rápido aumentan las pérdidas en función de la

distancia (Rappaport, 2002), (Cox, 1984).

He (2011) realizó mediciones de banda estrecha a 930 MHz en

un sector de la línea férrea de alta velocidad de China. El

modelo obtenido se basa en el modelo de pérdidas de

trayectoria de shadowing log-normal considerando el efecto de

la distancia de referencia.

Aboul-Dahab (2010) tomó en cuenta mediciones de

propagación efectuadas en diferentes lugares del planeta en la

banda de frecuencia de 3.5 GHz considerando redes con

diferentes topologías de la tecnología WIMAX. En el proyecto

se utiliza Least Square Aproximation para el cálculo de los

límites superiores e inferiores de las pérdidas de trayectoria y

se realiza una comparación entre diferentes modelos.

Los estudios desarrollados se han direccionado en determinar

modelos de propagación a partir de mediciones. Debido a que

los modelos existentes se han obtenido en base a mediciones

realizadas en lugares específicos del planeta que no son lo

suficientemente precisos si se aplican en regiones distintas.

Para la determinación de las pérdidas de propagación en la

ciudad de Quito se requieren modelos que se ajusten de mejor

Modelo Matemático Adaptado para el Cálculo de Pérdidas de

Propagación en la Banda de 900 MHz para Microceldas en la

Ciudad de Quito

Quistial, Alvin1; Lupera Morillo, Pablo1; Tipantuña, Christian1; Carvajal, Jorge1

1Escuela Politécnica Nacional, Departamento de Electrónica, Telecomunicaciones y Redes de Información, Quito, Ecuador

Resumen: En este artículo se presentan modelos para el cálculo de pérdidas de trayectoria en ambientes con línea de

vista y sin línea de vista para microceldas en la banda de 900MHz. Para obtener el modelo de pérdidas de trayectoria

se utiliza el método de aproximación por mínimo cuadrado y se lo determinó a partir de mediciones realizadas en la

ciudad de Quito de potencia recibida en diferentes puntos de una zona de cobertura. Las expresiones de los nuevos

modelos obtenidos se basan en los modelos existentes de espacio libre, Okumura, Okumura-Hata, COST-231, Egli y

Walfisch. Los métodos de ajuste de los modelos existentes a las mediciones fueron el del error cuadrático medio y el

método simple.

Palabras clave: Modelos de pérdidas de trayectoria, mediciones de potencia recibida, banda de 900 MHz,

microceldas.

An 900 MHz Adapted Path Loss Model for Microcells in Quito

Abstract: This paper presents a path-loss model for predicting signal propagation in Line of Sight (LOS) and Non

Line of Sight (NLOS) environments for microcell from 928 MHz to 960 MHz. For to determine the path-loss model

is used the method of Least Square Approximation. The determination of model for predicting the path-loss was made

through received power measurements at different points in a coverage area. The new model is based on existing

propagation models such as Free space model, Okumura model, Okumura-Hata model, Cost-231 model, Egli model

and Walfisch model. The model was obtained by using a method of adjustment such as Mean Square Error (MSE)

and simple method.

Keywords: Path-loss models, received power measurements, 900MHz, microcells.

29



manera a las condiciones geográficas y a las características de

propagaciones intrínsecas y estocásticas de dicha región.

El presente artículo se encuentra estructurado de la siguiente

manera. En la sección II se incluye el fundamento matemático

utilizado para los posteriores análisis, ya que se presentan las

expresiones matemáticas de los modelos de propagación

utilizados en este proyecto. En la sección III se describe la

metodología de medición que se aplicó para obtener los

valores de pérdida de trayectoria. En la sección IV se propone

la obtención del modelo matemático de pérdidas de trayectoria

adaptado a las condiciones de medición y a las características

del entorno. En la sección V se muestra el análisis de los

modelos adaptados de pérdidas de trayectoria; la evaluación se

realiza con dos métodos estadísticos, que permiten determinar

si el modelo obtenido es adecuado. En la sección VI se

describen las conclusiones del trabajo desarrollado.

2. MODELOS DE PROPAGACIÓN

Para este estudio se toman como referencia algunos modelos

de propagación conocidos para su posterior adaptación a las

mediciones realizadas en la ciudad de Quito. Las expresiones

matemáticas para el cálculo en dB de las pérdidas de

trayectoria con los modelos de propagación, se incluyen en la

siguiente tabla y constituyen el fundamento matemático para

el análisis posterior.

Tabla 1. Modelos de propagación analizados

Modelos de

Propagación

Expresiones matemáticas de pérdidas

Espacio, IEEE

Vehicular

(1988)

32,44 + 20 𝑙𝑜𝑔 𝑑(𝐾𝑚) + 20 𝑙𝑜𝑔 𝑓(𝑀𝐻𝑧)

Egli, IEEE

Vehicular

(1988)

139,1 + 40 𝑙𝑜𝑔 𝑑(𝐾𝑚) − 20 𝑙𝑜𝑔 ℎ𝑡𝑒 (𝑚)

Okumura,

IEEE

Vehicular

(1988)

𝐿𝐹 + 𝐴𝑚𝑢(𝑓, 𝑑) − 𝐺(ℎ𝑡𝑒) − 𝐺(ℎ𝑟𝑒) − 𝐺𝐴𝑅𝐸𝐴

Okumura /

Hata,

European

Commission

(1999), UIT-R

(2009)

69.55 + 26.16 𝑙𝑜𝑔(𝑓(𝑀𝐻𝑧)) − 13.82 𝑙𝑜𝑔(ℎ𝑡𝑒(𝑚))

− 𝑎(ℎ𝑟𝑒)

+ [44.9 − 6.55 𝑙𝑜𝑔(ℎ𝑡𝑒(𝑚))]

𝑙𝑜𝑔(𝑑(𝐾𝑚))

Cost -231,

European

Commission

(1999)

46.3 + 33.9 𝑙𝑜𝑔(𝑓) − 13.82 𝑙𝑜𝑔(ℎ𝑡𝑒(𝑚)) − 𝑎(ℎ𝑟𝑒)

+ [44.9

− 6.55 𝑙𝑜𝑔(ℎ𝑡𝑒)] 𝑙𝑜𝑔(𝑑(𝐾𝑚))

+ 𝑐𝑚

Walfisch,

Giménez

(2011),

European

Commission

(1999)

42,6 + 26 𝑙𝑜𝑔 𝑑(𝐾𝑚) + 20 𝑙𝑜𝑔 𝑓(𝑀𝐻𝑧)

A continuación se describen algunos de los parámetros que se

incluyen en la tabla anterior:

ℎ𝑡𝑒 - la altura de la antena transmisora,

ℎ𝑟𝑒 - la altura de la antena receptora,

𝐿𝐹 - Pérdidas de propagación en el espacio libre,

𝐴𝑚𝑢 - Atenuación media relativa del espacio libre que se

determina mediante curvas.

𝐺(ℎ𝑡𝑒) - Factor de ganancia de la altura de la antena de

transmisión:

𝐺(ℎ𝑡𝑒) = 20 log (ℎ𝑡𝑒(𝑚)

200) 1000 𝑚 > ℎ𝑡𝑒 > 30 𝑚

𝐺(ℎ𝑟𝑒) - Factor de ganancia de la altura de la antena móvil:

𝐺(ℎ𝑟𝑒) = 10 log (ℎ𝑟𝑒(𝑚)

3) ℎ𝑟𝑒 ≤ 3 𝑚

𝐺(ℎ𝑟𝑒) = 20 log (ℎ𝑟𝑒(𝑚)

3) 10 𝑚 > ℎ𝑟𝑒 > 3 𝑚

𝐺𝐴𝑅𝐸𝐴 - Ganancia debido al tipo de entorno (Factor de

Corrección) que se determina mediante curvas.

𝑎(ℎ𝑟𝑒) - Factor de corrección por la altura de la antena móvil

o receptora, que depende del tamaño de la zona de cobertura.

Para una ciudad pequeña o de tamaño medio, el factor de

corrección de la antena móvil 𝑎(ℎ𝑟𝑒) está dada por:

𝑎(ℎ𝑟𝑒) = [1.1 log(𝑓(𝑀𝐻𝑧)) − 0.7]ℎ𝑟𝑒(𝑚)

− [1.56 log(𝑓(𝑀𝐻𝑧)) − 0.8]

Para una ciudad grande, el factor de corrección 𝑎(ℎ𝑟𝑒) está

dado por:

𝑎(ℎ𝑟𝑒) = 8.29[log(1.54hre(m))]2

− 1.1 𝑓𝑐 ≤ 300 𝑀𝐻𝑧

𝑎(ℎ𝑟𝑒) = 3.2[log(11.75hre(m))]2

− 4.97 𝑓𝑐 ≥ 300 𝑀𝐻𝑧

𝑐𝑚- Factor de corrección.

El factor 𝑐𝑚, se establece en Dalela (2012), tiene un valor de

0 dB para zonas suburbanas y zonas abiertas, y 3 dB para zonas

urbanas.

3. METODOLOGÍA DE MEDICIÓN

La metodología de medición se dividió en seis etapas que se

describen a continuación.

3.1 Definición de la zona de medición

Inicialmente se identificó la zona de cobertura donde se

realizarán las mediciones, ya que un factor muy importante

que influye en las pérdidas de trayectoria es la característica

geográfica del lugar. Esta zona se eligió para que permita tener

una referencia adecuada y sea posible generalizar los

resultados. La zona de cobertura seleccionada fué el campus

universitario de la Escuela Politécnica Nacional en la ciudad

de Quito.

Debido a que este proyecto se lo realiza para el cálculo de

pérdidas de trayectoria en microceldas, la distancia de

cobertura no sobrepasa los 350 metros.

3.2 Puntos de medición

Una vez identificada la zona de cobertura se seleccionó con

exactitud un número importante de puntos donde se realizaron

las mediciones de la potencia recibida. Y se identifican los

escenarios que se tienen entre el transmisor y receptor para

diferenciar la influencia de cada uno de ellos.

30



Debido a las características de la zona de cobertura, se

escogieron 39 puntos de medición, que se encuentran

distribuidos considerando ángulos de azimut de 15°, es decir

que desde el transmisor se tienen 10 rutas de medición. En

cada azimut las distancias de separación de la antena receptora

hacia la antena transmisora fué de 70, 140, 210, 280 y 350

metros. Para la ubicación precisa de los puntos de medición se

utilizó un GPS.

En la Figura 1 se observa la distribución de los puntos de

medición en la zona de cobertura planteada y la ubicación del

transmisor.

Figura 1. Ubicación de los puntos de medición

3.3 Escenarios de medición

En la zona seleccionada se identificaron tres escenarios: línea

de vista, sin línea de vista con obstrucción por edificios y sin

línea de vista con obstrucción por árboles. La agrupación de

las mediciones por el tipo de escenario permitió obtener

valores precisos de potencia de recepción, que se tomó en

cuenta para obtener modelos de pérdidas de trayectoria para

cada escenario.

3.4 Selección del rango de frecuencia de medición

El rango de frecuencias seleccionado se encuentra desde los

928 MHz hasta los 960 MHz. Utilizando un analizador de

espectros se verificó que dicho rango de frecuencias no se

encontraba ocupado durante la realización de las mediciones.

Las mediciones de los parámetros de antenas, potencia de

transmisión y recepción se ejecutaron cada 1 MHz.

3.5 Potencia real de transmisión y recepción

Con el uso de un medidor de potencia RF se estableció con

exactitud la potencia real transmitida a la salida del generador

de señales que se utilizó como transmisor.

Utilizando un analizador vectorial de redes se determinaron los

parámetros de VSWR de las antenas transmisora y receptora

en el rango de frecuencias de medición, con la finalidad de

establecer el nivel de potencia real transmitido aplicando el

parámetro de eficiencia de acoplamiento de las antenas.

En el lado del receptor se utilizó un analizador de espectro

configurado con las características de la antena receptora para

registrar la potencia real recibida.

Los equipos utilizados en este proyecto se enumeran en la

Tabla 2.

Tabla 2. Equipos utilizados en las mediciones

3.6 Especificación de las condiciones de medición

Las condiciones de medición se especifican en la Tabla 3 y por

tanto definen los rangos de aplicación de los modelos que se

obtuvieron en el presente trabajo.

Tabla 3. Condiciones de las mediciones

Condiciones de medición

Rango de Frecuencias 928 a 960 MHz

Distancia 70 a 350 metros (microceldas)

Altura de la antena

Transmisora 26 metros

Altura de la antena Receptora 2 metros

Escenarios

Línea de vista, sin línea de

vista (obstáculos: edificios y

árboles)

4. OBTENCIÓN DEL MODELO ADAPTADO DE

PÉRDIDAS DE TRAYECTORIA PARA LA CIUDAD

DE QUITO

El modelo matemático adaptado de pérdidas de trayectoria se

obtiene a partir de los modelos de propagación existentes, los

cuales se ajustan intencionalmente a las mediciones obtenidas.

En la figura 2 se observa la clasificación de los modelos de

propagación considerados (Saunders, 2007).

LISTA DE EQUIPOS UTILIZADOS

EQUIPO MARCA MODELO RANGO DE

FRECUENCIA

GENERADOR

DE SEÑALES Anritsu MG3691C 0,1 Hz a 10 GHZ

MEDIDOR DE

POTENCIA Tektronix PSM4120 10 MHz a 8 GHz

ANALIZADOR

VECTORIAL Keysight N9914A 30 kHz a 6,5 GHz

ANALIZADOR

DE ESPECTROS Anritsu MS2711E 9 kHz a 3 GHz

ANTENA LOG-

PERIÓDICA

(TRANSMISORA)

A. H.

Systems SAS-510-4 290 MHz a 4 GHz

ANTENA

MINI.MAG

(RECEPTORA)

Smarteq 1140.26SMA 824 a 960 MHz – 1710 a 2170 MHz

31



Figura 2. Clasificación de los modelos de propagación considerados

Las mediciones se realizaron en el campo, registrándose tres

valores de la potencia recibida en cada ubicación: potencia

mínima, potencia máxima y potencia promedio. Se observó

que la potencia recibida en los puntos de medición depende en

gran medida, además de la distancia, de las obstrucciones

presentes en el enlace entre el transmisor y el receptor.

Como ejemplo, en la Tabla 4 se presentan los valores medidos

de la potencia media recibida en los puntos de medición a la

frecuencia de 944 MHz, esta tabla se obtiene de una parte de

todas las mediciones realizadas.

En la figura 3 se gráfica la potencia media recibida en función

del logaritmo de la distancia.

Con los valores medidos de la potencia media recibida se

calcularon las pérdidas de trayectoria L, utilizando la siguiente

ecuación:

L = PTx + GTx + GRx − PRx (1)

Donde PTx y PRx representan la potencia real de transmisión y

la potencia real de recepción, GTx y GRx representan las

ganancias de transmisión y recepción respectivamente.

Con los resultados obtenidos, a continuación se procede a

determinar el modelo de pérdidas de trayectoria Logaritmo-

Distancia, que se expresa de la siguiente manera:

L̅(dB) = L̅(d0) + 10n log (d

d0) (2)

Donde d0 es la distancia de referencia de 1 metro, n es el

exponente de pérdida de trayectoria, L̅(d0) representa las

pérdidas inciales y d es la distancia de separación entre el

transmisor y el receptor en metros. El valor del exponente de

pérdida de trayectoria depende del entorno de propagación

específico, Aboul-Dahab (2010) indica que el exponente de

pérdida de trayectoria (n) es igual a 2 cuando la propagación

se da en el espacio libre.

Tabla 4. Potencia recibida a 944 MHz en [dBm]

Puntos Azimut

Distancia

70 m 140

m

210

m 280 m 350 m

1 0 -

42,52

-

47,69

-

57,87 -64,92

Fuera de

cobertura

2 15 -

44,11

-

55,37

-

56,99

Fuera de

cobertura

Fuera de

cobertura

3 30 -45,2 -

55,19

-

61,26

Fuera de

cobertura

Fuera de

cobertura

4 45 -

46,63

-

51,18

-

60,97

Fuera de

cobertura

Fuera de

cobertura

5 60 -

46,22

-

51,31

-

61,68

Fuera de

cobertura

Fuera de

cobertura

6 75 -

43,39

-

64,78

-

68,79 -68,5

Fuera de

cobertura

7 90 -

46,81

-

54,48

-

69,38 -63,72 -61,86

8 105 -

47,77

-

52,33

-

57,14 -57,7 -57,32

9 120 -

52,73

-

59,89

-

60,14 -62,03 -64,33

10 135 -

41,79

-

48,89

-

57,38 -56,8

Fuera de

cobertura

El modelo de pérdidas Logaritmo-Distancia permite establecer

una línea de tendencia de regresión lineal que relaciona las

pérdidas de trayectoria con la distancia. Los valores de los

coeficientes 10n y L̅(d0) se pueden determinar utilizando dos

métodos: el método de mínimos cuadrados o el método

matricial por eliminación de Gauss-Jordan (Alsayyari, 2015).

Obtenidos los valores de dichos coeficientes y reemplazando

en la ecuación (2) se puede obtener el modelo de pérdida

Logaritmo-Distancia para la zona de cobertura.

En la Figura 4 se aprecian las mediciones obtenidas en las

ubicaciones y la recta de regresión lineal obtenida.

Figura 3. Potencia media recibida a 944MHz en función del logaritmo de la

distancia.

32



Figura 4. Recta de regresión lineal de las pérdidas de trayectoria.

En este trabajo se propone adaptar los modelos de propagación

existentes a las mediciones mediante el método de ajuste del

Error Cuadrático Medio (MSE). Este método permite obtener

un valor numérico de ajuste para alcanzar el menor error

posible entre dos valores. La ecuación del Error Cuadrático

Medio (MSE) es la siguiente:

Donde Lm es el valor medio de las pérdidas de trayectoria

obtenido de las mediciones, Lp es el valor de pérdidas de

trayectoria calculados con los modelos de propagación

existentes y q es el número de puntos de medición.

El valor de ajuste se aumenta o disminuye a la ecuación de

pérdidas de trayectoria de los modelos de propagación

existentes para adaptarlos a las mediciones.

En la Figura 5 se observan las líneas de tendencia de las

mediciones para enlaces con línea de vista y de los modelos de

propagación de Okumura y Espacio Libre sin ajuste.

En la Figura 6 en cambio se aprecian los modelos de

propagación después del ajuste con el Método MSE.

Para enlaces sin línea de vista los modelos de propagación se

adaptaron a las mediciones con el método de ajuste simple. La

ecuación de ajuste simple es la siguiente:

𝐴𝑗𝑢𝑠𝑡𝑒 = Lm − Lp (4)

Este valor de ajuste, al igual que en el método MSE, se

aumenta o disminuye a la ecuación de pérdidas de trayectoria

de los modelos de propagación existentes.

Figura 5. Modelos de propagación antes del ajuste para enlaces con línea de

vista.

En la Figura 7 se visualizan las pérdidas de trayectoria

obtenidas con los modelos existentes ajustados para enlaces

sin línea de vista en escenarios con obstrucción por árboles.

Para enlaces sin línea de vista en un escenario con obstrucción

por edifícios, los modelos de pérdidas de trayectoria se ajustan

de acuerdo a lo que se muestra en la Figura 8.

Figura 6. Modelos de propagación ajustados con el método MSE para

enlaces con línea de vista.

MSE = √∑(Lm − Lp)

2

q − 1 (3)

33



Figura 7. Modelos ajustados de pérdidas de trayectoria con el método

simple en escenarios sin línea de vista por árboles.

Figura 8. Modelos ajustados de pérdidas de trayectoria con el método

simple en escenarios sin línea de vista por edificios.

En la Tabla 5 se incluyen las ecuaciones obtenidas de los

modelos adaptados de pérdidas de trayectoria obtenidas con el

método de ajuste MSE. Mientras que en la Tabla 6 se

encuentran los modelos adaptados con el método de ajuste

simple. En las expresiones adaptadas las distancias se indican

en metros, ya que la cobertura de las microceldas es menor a

los kilómetros.

Tabla 5. Ecuaciones de los modelos adaptados obtenidas con el método

MSE

Modelos de

Propagación

Modelo Adaptado para escenarios con línea de

vista

Espacio Libre −22,719 + 20 log d(m) + 20 log f(MHz)

Okumura −23,356 + 20 log d(m) + 20 log f(MHz)

Modelos de

Propagación

Modelo Adaptado para escenarios sin línea de vista

obstrucción por edificios

Espacio Libre −15,751 + 20 log d(m) + 20 log f(MHz)

Okumura −16,952 + 20 log d(m) + 20 log f(MHz)

Walfisch −28,787 + 26 log d(m) + 20 log f(MHz)

Modelos de Propagación

Modelo Adaptado para escenarios sin línea de vista obstrucción por árboles

Egli −2,561 + 40 log d(m)

Okumura / Hata −70,65 + 35,63 log d(m) + 24,42 log f(MHz)

Cost -231 −97,80 + 35,632 log d(m) + 33,9 log f(MHz)

Tabla 6. Ecuaciones de los modelos adaptados obtenidas con el método simple

Modelos de

Propagación Modelo Adaptado para escenarios con línea de vista

Espacio Libre 𝟑𝟏, 𝟔𝟓 + 𝟐𝟐, 𝟐𝟐 𝐥𝐨𝐠 𝐝(𝐦)

Egli 𝟓𝟗, 𝟗𝟓 + 𝟐𝟐, 𝟐𝟐 𝐥𝐨𝐠 𝐝(𝐦) − 𝟐𝟎 𝐥𝐨𝐠 𝐡𝐭(𝐦)

Okumura 𝟐𝟑, 𝟔𝟓 + 𝟐𝟐, 𝟐𝟐 𝐥𝐨𝐠 𝐝(𝐦) + 𝐀𝐦𝐮(𝐟,𝐝) − 𝐆𝐀𝐑𝐄𝐀

Okumura /

Hata 𝟑𝟏, 𝟔𝟓 + 𝟐𝟐, 𝟐𝟐 𝐥𝐨𝐠 𝐝(𝐦)

Cost -231 𝟐𝟖, 𝟔𝟓 + 𝟐𝟐, 𝟐𝟐 𝐥𝐨𝐠 𝐝(𝐦) + 𝐜𝐦

Walfisch 𝟑𝟏, 𝟔𝟓 + 𝟐𝟐, 𝟐𝟐 𝐥𝐨𝐠 𝐝(𝐦)

Modelo de

Propagación

Modelo Adaptado para escenarios sin línea de vista con

obstrucción por árboles

Espacio Libre 𝟑𝟖, 𝟑𝟓 + 𝟐𝟎, 𝟓 𝐥𝐨𝐠 𝒅(𝒎)

Egli 𝟔𝟔, 𝟔𝟓 + 𝟐𝟎, 𝟓 𝐥𝐨𝐠 𝒅(𝒎) − 𝟐𝟎 𝐥𝐨𝐠 𝒉𝒕(𝒎)

Okumura 𝟑𝟎, 𝟑𝟓 + 𝟐𝟎, 𝟓 𝐥𝐨𝐠 𝒅(𝒎) + 𝑨𝒎𝒖(𝒇,𝒅) − 𝑮𝑨𝑹𝑬𝑨

Okumura / Hata

𝟑𝟖, 𝟑𝟓 + 𝟐𝟎, 𝟓 𝐥𝐨𝐠 𝒅(𝒎)

Cost -231 𝟑𝟓, 𝟑𝟓 + 𝟐𝟎, 𝟓 𝐥𝐨𝐠 𝒅(𝒎) + 𝒄𝒎

Walfisch 𝟑𝟖, 𝟑𝟓 + 𝟐𝟎, 𝟓 𝐥𝐨𝐠 𝒅(𝒎)

Modelo de Propagación

Modelo Adaptado para escenarios sin línea de vista con obstrucción por edificios

Espacio Libre −𝟐𝟕, 𝟖𝟔 + 𝟓𝟎, 𝟒 𝐥𝐨𝐠 𝒅(𝒎)

Egli 𝟎, 𝟒𝟑𝟗 + 𝟓𝟎, 𝟒 𝐥𝐨𝐠 𝒅(𝒎) − 𝟐𝟎 𝐥𝐨𝐠 𝒉𝒕(𝒎)

Okumura −𝟑𝟓, 𝟖𝟔 + 𝟓𝟎, 𝟒 𝐥𝐨𝐠 𝒅(𝒎) + 𝑨𝒎𝒖(𝒇,𝒅) − 𝑮𝑨𝑹𝑬𝑨

Okumura / Hata

−𝟐𝟕, 𝟖𝟔 + 𝟓𝟎, 𝟒 𝐥𝐨𝐠 𝒅(𝒎)

Cost -231 −𝟑𝟎, 𝟖𝟔 + 𝟓𝟎, 𝟒 𝐥𝐨𝐠 𝒅(𝒎) + 𝒄𝒎

Walfisch −𝟐𝟕, 𝟖𝟔 + 𝟓𝟎, 𝟒 𝐥𝐨𝐠 𝒅(𝒎)

5. ANÁLISIS DE ERROR DE LOS MODELOS

ADAPTADOS

Para comprobar si las ecuaciones de los modelos adaptados de

pérdidas de trayectoria se ajustan a las mediciones obtenidas,

en la Tabla 7 se realiza un análisis de la exactitud del ajuste

mediante el contraste de regresión (𝐹1,𝑛−(𝑘+1)) y el coeficiente

de determinación (𝑅2).

De los resultados se observa que la utilización del método de

ajuste del Error Cuadrático Medio (MSE) permite obtener

modelos aceptables cuando se lo aplica para los enlaces de

transmisión con línea de vista, ya que en este escenario se

obtienen modelos adaptados de pérdidas de trayectoria con

77

82

87

92

97

2.13 2.18 2.23 2.28 2.33

Pér

did

a d

e tr

ayec

tori

a (d

B)

Distacia-log (d)Mediciones

Lineal (Mediciones)

Lineal (Modelo Ajustado Espacio Libre)

Lineal (Modelo Ajustado Egli)

Lineal (Modelo Ajustado Okumura)

Lineal (Modelo Ajustado Okumura/Hata)

Lineal (Modelo Ajustado Cost-231)

Lineal (Modelo Ajustado Walfisch)

34



valores promedio de contraste de regresión de 2609,851 y

coeficiente de determinación de 0,85.

El método de ajuste simple permite obtener modelos adaptados

aceptables para todos los escenarios de aplicación.

Tabla 7. Contraste de regresión y coeficiente de determinación de los

modelos adaptados

6. CONCLUSIONES

Las mediciones efectuadas en la zona de cobertura mostraron

que se requiere un ajuste a los modelos de propagación

conocidos.

En este proyecto el ajuste de los modelos de propagación se

realizó mediante el método del Error Cuadrático Medio y el

ajuste simple.

Del presente proyecto se han obtenido modelos de

propagación adaptados que pueden ser utilizados como

referencia para el diseño de microceldas en la banda de 900

MHz en entornos urbanos de la ciudad de Quito considerando

condiciones similares a las establecidas en las mediciones que

se efectuaron.

REFERENCIAS

Aboul-Dahab M. A., Kamel H. M. (2010), Methodology for Calculating Path

Loss Upper and Lower Bounds for WiMAX; IEEE.

Alsayyari A., Kostanic I., Otero C. E., (2015), An Empirical Path Loss Model for Wireless Sensor Network Deployment in a Concrete Surface

Enviroment, IEEE.

Cox D. C., Murray R. R., and Norris A. W., (1984), 800 MHz attenuation measured in and around suburban houses, AT&T Bel1 Lab. Tech.J., vol.

63, pp. 922-951.

Dalela Ch., Prasad M., Dalela P.; (2012), Tuning of Cost-231 Hata Model for

Radio Wave Propagation Predictions; CS & IT, pp. 255-267.

European Commission; (1999), COST Action 231, Digital mobile radio

towards future generation systems, Final report, pp. 115-197. Giménez J, López J., Gómez-Barquero D., Cardona D., (2011), Modelos de

propagación radio para redes de TDT móvil en la banda UHF, Revista

S&T Universidad ICESI, Colombia, pp. 9-27. He R., Zhong Z., and Ai B., (2010), Path loss measurements and analysis for

high-speed railway viaduct scene, in Proc. 4th Int. Wirel. Commun.

Mob. Comput. Conf., France, pp. 265-270. He R., Zhong Z., Ai B., Xiong L., and Ding J., (2011), The effect of Reference

Distance on Path Loss Prediction Based on the Measurements in High-

Speed Railway Viaduct Scenarios, IEEE, Chinacom. IEEE Vehicular Technology Society Committee on Radio Propagation;

(1988), Coverage Prediction for Mobile Radio Systems Operating in the

800/900 MHz Frequency Range, IEEE Transactions on Vehicular Technology, vol. 37, No.1, pp. 20-25.

Rappaport T. S., (2002), Wireless Communications: Principle and Practice,

2nd Edition ed., Upper Saddle River, NJ, USA: Prentice Hall.

Saunders S. R., Aragón-Zavala A.; (2007), Antennas and Propagation for

Wireless Communication Sysrtems, 2nd Edition, John Wiley & Sons ltd,

The Atrium, Southern Gate, Chichester, West Sussex, England, p.p. 89-91, 97-99, 163-165, 167-169, 178-180, 182.

UIT-R; (2009), P.1546-4 Métodos de Predicción de Punto a Zona para

Servicios Terrenales en la Gama de Frecuencias de 30 a 3000 MHz.

BIOGRAFÍAS

Alvin Patricio Quistial, realizó sus estudios

secundarios en el Instituto Nacional Mejía,

obteniendo el título de bachiller Físico-

Matemático. Posteriormente obtiene el título

de Ingeniero en Electrónica y

Telecomunicaciones en la Escuela Politécnica

Nacional, en 2016.

Pablo Lupera Morillo, obtuvo el título de

ingeniero en Electrónica y

Telecomunicaciones en la Escuela Politécnica

Nacional en el año 2002 y el título de Ph.D. en

ciencias técnicas en la Universidad Estatal de

Telecomunicaciones de San Petersburgo en

Rusia en el año 2009. Sus áreas de

investigación son el comportamiento del canal inalámbrico,

técnicas de transmisión aplicadas a la capa física y el diseño

de redes de comunicación móvil.

Christian José Tipantuña Tenelema,

graduado en Ingeniería en Electrónica y

Telecomunicaciones por la Escuela

Politécnica Nacional, en 2011. Sus estudios

de maestría los realizó en el Politécnico di

Torino, Turin-Italia, en 2013. Actualmente es

profesor asistente en la Escuela Politécnica

Nacional.

Jorge Eduardo Carvajal, realizó sus estudios

en la Escuela Politécnica Nacional,

obteniendo el título de ingeniero en

Electrónica y telecomunicaciones, sus

estudios de cuarto nivel los realizó en la

universidad de ciencias aplicadas Mannheim

obteniendo el título de máster en Tecnologías

de la Información. Actualmente trabaja en la

Escuela Politécnica Nacional como profesor auxiliar.

Modelos de

Propagación para

escenarios con

línea de vista

Método MSE Método simple

𝑭𝟏,𝒏−(𝒌+𝟏) 𝑹𝟐 𝑭𝟏,𝒏−(𝒌+𝟏) 𝑹𝟐

Espacio Libre 2554,53 0,847 1,09x10^11 0,999

Egli No se aplica 6,88x10^10 0,999 Okumura 2665,17 0,853 1,23x10^11 0,999

Okumura / Hata No se aplica 1,79x10^9 0,999

Cost-231 No se aplica 1,22x10^11 0,999

Walfisch No se aplica 1187164,46 0,999 Modelos de

Propagación para

escenarios sin línea de vista obstrucción

por edificios


𝐹1,𝑛−(𝑘+1) 𝑅2 𝐹1,𝑛−(𝑘+1) 𝑅2

Espacio Libre 393,696 0,413 9,31x10^10 0,999

Egli No se aplica 1,03x10^11 0,999

Okumura 391,095 0,412 1,02x10^11 0,999 Okumura/Hata No se aplica 1,69x10^9 0,999

Cost-231 No se aplica 1,16x10^11 0,999

Walfisch 634,827 0,532 55202,158 0,999 Modelos de

Propagación para

escenarios sin línea de vista obstrucción

por árboles


𝐹1,𝑛−(𝑘+1) 𝑅2 𝐹1,𝑛−(𝑘+1) 𝑅2

Espacio Libre No se aplica 3,58x10^10 0,999 Egli 393,696 0,413 3,83x10^10 0,999

Okumura No se aplica 3,82x10^10 0,999

Okumura/Hata 391,095 0,412 6,62 x10^8 0,999 Cost-231 634,827 0,532 4,63x10^10 0,999

Walfisch No se aplica 137606,8 0,999

35

Modelamiento y Simulación de la Producción de Hidrógeno en un Electrolizador a Partir de Vapor Sobrecalentado de Agua





11. INTRODUCCIÓN

Desde hace varias décadas, el hombre ha planteado el uso

combinado de la tecnología y de los combustibles fósiles para

dar paso a otras tecnologías que facilitan y permiten acceder

en la actualidad a diferentes servicios, pero lamentablemente

tanto el petróleo como el gas natural, entre las fuentes de

energía más representativas, son recursos finitos que han

llegado a provocar altos niveles de contaminación que afectan

al ecosistema y por ende a la humanidad. Ante tales

circunstancias, resulta necesario cambiar la base energética

actual por nuevas formas de energía que sean accesibles,

abundantes, de bajo costo, tratables y que reduzcan la

contaminación ambiental. Una de las nuevas fuentes

innovadoras de energía es el hidrógeno, caracterizado por ser

eficiente, inagotable, seguro, almacenable, transportable y

económico de producir. Como elemento químico en

condiciones normales de presión y temperatura es un gas

diatómico (H_2) liviano, incoloro, inodoro, insípido e

inflamable. Su átomo se caracteriza por tener un protón y un

electrón, además que en la naturaleza se lo encuentra como

parte de una gran cantidad de compuestos como el agua,

hidrocarburos, proteínas o ácidos. Su molécula es pequeña,

ligera, simple y relativamente estable, posee el más alto

contenido de energía por unidad de peso que cualquier otro

combustible y al combinarse con el oxígeno produce

electricidad en procesos electroquímicos mientras que al

quemarse con oxígeno produce vapor de agua lo que resulta en

una reducción de la contaminación. Finalmente, presenta un

Modelamiento y Simulación de la Producción de Hidrógeno en un

Electrolizador a Partir de Vapor Sobrecalentado de Agua

Aguilar-Jaramillo, Edwin1; Aguinaga, Álvaro2

1 Universidad de la Fuerzas Armadas ESPE, Departamento de Eléctrica y Electrónica, Quito, Ecuador

2 Escuela Politécnica Nacional EPN, Departamento de Ingeniería Mecánica, Quito, Ecuador

Resumen: La investigación presenta un modelo eléctrico equivalente de una celda electrolítica simple, utilizada para

la producción de hidrógeno a partir de vapor sobrecalentado de agua, con la finalidad de brindar una base teórica que

pueda ser aplicada al momento de elegir o diseñar equipos para sistemas de electrólisis de este tipo. Para ello se

realiza el estudio de los componentes que conforman el electrolizador, las leyes y ecuaciones que rigen el

comportamiento del sistema, las variables que inciden en la producción de hidrógeno, así como la espectroscopia de

impedancia electroquímica, esta última útil para interpretar de forma matemática los fenómenos que se presentan en

el proceso de electrólisis. La implementación y simulación del modelo análogo por circuito eléctrico equivalente

(EEC) de la celda electrolítica, utilizando las características de un electrolito de YSZ, un ánodo de LSM y un cátodo

de Ni-YSZ, permitió concluir que se pueden obtener resultados eficientes en la producción de hidrógeno a mayor

área transversal del electrolito y menor separación entre los electrodos, al trabajar con una señal de excitación de

pulsos a frecuencia de resonancia y a temperaturas de vapor sobrecalentado altas.

Palabras clave: electrólisis, producción de hidrógeno, espectroscopia de impedancia electroquímica, circuito

eléctrico equivalente.

Modeling and Simulation of Hydrogen Production in an

Electrolyzer from Superheated Water Vapor

Abstract: The research presents an electrical model equivalent of a simple electrolytic cell, used for the production

of hydrogen from overheated water steam, in order to provide a theoretical basis that can be applied when choosing

or designing equipment for electrolysis systems of this type. To do that, will be carried out the component study of

electrolyzer, laws and equations governing the behavior of the system; the variables that affect the production of

hydrogen as well as the electrochemical impedance spectroscopy, this last one useful for mathematically interpret the

phenomena that occur in the process of electrolysis. The implementation and simulation of the analog model, of the

equivalent electric circuit (EEC) of the electrolytic cell, using the characteristics of a YSZ electrolyte, LSM anode

and a Ni-YSZ cathode, allowed to conclude that optimal results in the hydrogen production can be obtained to greater

area of the electrolyte and smaller spacing between electrodes, working with an excitation signal of pulses at

resonance frequency and a high overheated water steam temperatures.

Keywords: electrolysis, hydrogen production, solid oxide electrolyzer cell, electrochemical impedance spectroscopy,

equivalent electrical circuit.

37



margen de explosión bajo, respecto a otros combustibles

debido a su densidad.

Uno de los métodos para la obtención de este gas ha sido por

mucho tiempo la electrólisis convencional cuyo inconveniente

ha llegado a ser su alto consumo de electricidad. En vista de

esto, producto de varios estudios, se llega a la conclusión de

que la descomposición del agua a altas temperaturas resulta ser

más ventajosa, pues con esto parte de la energía total, usada

para producir la electrólisis, ahora será suministrada o añadida

en forma de calor, para así obtener como resultado una notable

mejora en el índice de eficiencia total de energía.

Considerando lo anteriormente mencionado, el presente

trabajo busca modelar matemáticamente la producción de

hidrógeno en un electrolizador mediante el uso de vapor

sobrecalentado de agua a fin de brindar una base teórica que

pueda ser aplicada al momento de elegir o diseñar equipos para

sistemas de electrólisis que permitan producir hidrógeno de

forma óptima con un consumo de energía eficiente; a la vez

que se pretende incentivar todo proceso de investigación

encaminada al aprovechamiento del hidrógeno como vector

energético en el país.

2. MARCO TEÓRICO/METODOLOGÍA

2.1. Características del electrolizador

El electrolizador es el dispositivo físico en el cuál se produce

la descomposición del vapor de agua, en este caso, en

hidrógeno y oxígeno al hacer circular una corriente eléctrica

de forma conveniente. La celda simple, configuración del

electrolizador a utilizarse, está compuesta por dos electrodos,

ánodo y cátodo con un electrólito en medio operado en un

ambiente dual. Dicha configuración se presenta a continuación

en la Figura 1.

Figura 1. Celda simple de electrolito a base de óxido sólido. Fuente:

Annabelle Brisse, J. S. (2008).

En base a las características físicas y químicas con las que

deben cumplir los componentes del electrolizador para el

trabajo a altas temperaturas, de acuerdo con Alvarado (2013)

y Cachadiña, et al. (2001), los materiales del electrolito, cátodo

y ánodo seleccionados son YSZ (circonia estabilizada con

itrya), Ni-YSZ (níquel circonio estabilizado con escandio.) y

LSM (manganito de lantano dopado con estroncio),

respectivamente.

2.2. Reacciones químicas

Siendo el agua una de las sustancias más estables de la

naturaleza, su disociación, para la producción de hidrógeno y

oxígeno moleculares, no es espontánea y posee un alto grado

de consumo de energía (Pierre y Capitaine, 2006).

La reacción neta de electrólisis del agua se muestra en la

Ecuación (1).

H2O+ calor + electricidad → H2+ 1

2O2 ( 1)

Que indica que, por cada mol de agua ingresada, se obtiene

una mol de hidrógeno y medio mol de oxígeno.

En una celda electrolítica la energía eléctrica se suministra

aplicando una diferencia de potencial entre dos electrodos

acoplados con un electrolito. La conversión de la energía

eléctrica en energía química tiene lugar en la interfase

electrodo-electrolito a través de reacciones de transferencia de

carga.

Las reacciones que se producen en cada parte del

electrolizador de acuerdo a Brisse et al. (2008) y Nechache et

al. (2014) son:

En el cátodo se produce la reacción como se muestra en la

Ecuación (2).

2H2O + 4e− → 2H2 + 2O2− ( 2)

En el ánodo la reacción producida se muestra en la Ecuación

(3).

2O2− → O2 + 4e− ( 3)

En el electrolito la reacción es: 2O2−; circulando de forma

constante. Todas estas reacciones se resumen en la Figura 2.

Figura 2. Celda de electrólisis por vapor sobrecalentado. Fuente: Annabelle

Brisse, J. S. (2008).

38



Cabe destacar, que todo el proceso comprende un estado de

oxidación mientras esté en funcionamiento, es decir que el

compuesto inicial se disgrega por la energía aplicada, en

compuestos que tienden a combinarse debido a su

inestabilidad (como el oxígeno que circula por el electrolito,

hasta combinarse en 2 moléculas de sí mismo para formar el

gas de oxígeno).

2.3. Determinación de parámetros y variables que intervienen

en el proceso de generación de hidrógeno por electrólisis a

partir de vapor sobrecalentado.

Tomando en cuenta lo descrito por Mazloomi et al. (2012), se

consideran a continuación las variables que pueden influir

directamente en la obtención de hidrógeno a partir de vapor

sobrecalentado en una celda electrolítica simple. Dichas

variables son:

2.3.1. Temperatura

La temperatura (T) es la variable determinante en la

efectividad del proceso pues las características

termodinámicas del agua afirman un aumento de la energía

cinética de la sustancia que ayuda a reducir drásticamente el

uso de la energía eléctrica según Zoulias et al. (2010).

Adicionalmente, se conoce también que el calor puede reducir

el potencial reversible del agua (conocido también como el

equilibrio de voltaje). Un rango de temperatura con el cual se

han obtenido buenos resultados de electrólisis al utilizar YSZ

como electrolito es según el análisis de Brisse et al. (2008),

Cachadiña et al. (2001) y Nechache et al. (2014):

T = [773; 1 873] (K)

Lo cual indica que el proceso depende principalmente de altas

temperaturas, pues es en este rango en donde ocurren las

transferencias de electrones óptimos en los electrodos y

electrolito.

El valor de la temperatura con la que se trabajará inicialmente

para posteriormente proceder a variarla será de 1 123,15 (K).

2.3.2. Área transversal, alineación y espacio entre los

electrodos.

Se debe resaltar que el área transversal del electrolito es

equivalente al área transversal de los electrodos A (m2) por la

estructura de la celda, por lo que en este trabajo el área

calculado (A) en las fórmulas del electrolito corresponderá

recíprocamente al área transversal de los electrodos, parámetro

a considerar y variar dentro del proceso de la producción de

hidrógeno.

Se debe tomar en cuenta el área, ya que la impedancia de un

elemento depende en gran medida de ella, obteniendo una

menor resistencia y un flujo más rápido de los electrones.

Se han observado efectos del vacío entre interfaces, debido a

la formación de burbujas de gas (Mazloomi et al., 2012), pero

que no alteran de forma significativa el proceso de electrólisis

por la naturaleza del vapor. La resistencia a la corriente se

reduce si el área sobre la cual actúa es mayor en comparación

con su longitud, por lo que la altura de un electrodo (parte

superior del electrodo) es un factor poco determinante en el

proceso. Además, el poner el electrodo en posición vertical

aumenta la eficiencia que puede ser obtenida.

Esto último es causado por la reducida resistencia óhmica que

existe en esta zona (Mazloomi et al., 2012).

De acuerdo con Brisse et al. (2008) y Nechache et al. (2014)

los valores iniciales de espesor del electrolito (que están

relacionados con la separación de los electrodos d) y área

transversal de los electrodos para la interfaz cátodo-electrolito

ánodo electrólito “A” (m2) correspondientes para una

temperatura de 1 123,15 (K) serán:

d = 5 · 10−3 (m)

A = 1,6 · 10−3 (m2)

2.3.3. Voltaje suministrado

Los estudios para una electrólisis a partir de vapor

sobrecalentado de agua estable descritos por Brisse et al.

(2008) y Nechache et al. (2014) han demostrado una tendencia

a los valores pequeños en cuanto se refiere al voltaje

suministrado; es así que se tienen los siguientes valores para

un rango de 1 073 (K) a 1 272 (K):

V = [1,64 1,23 1,5] (V)

A un valor de flujo de carga o densidad de corriente;

j = [16 000 3 000 10 000] (A/m2).

De estos valores, el voltaje de entrada de alimentación de la

celda VS (V) correspondiente para las condiciones iniciales de

la celda electrolítica, que ha permitido lograr una mejor

eficiencia del proceso es:

VS = 1,33 (V)

2.3.4. Ciclo de trabajo

Al trabajar con una fuente tipo pulsos a frecuencia de

resonancia, se observa que mejora la eficiencia energética del

sistema, en lugar de trabajar solo a un valor de fuente DC

constante.

Esta frecuencia de resonancia, para este caso, será calculada a

partir del modelo del circuito eléctrico equivalente EEC

definido para celda electrolítica simple. Adjunto con ello se

analiza la variación del ciclo de trabajo dentro del siguiente

rango:

ciclo de trabajo = [0 → 100] (%)

2.4. Modelo de una celda electrolítica de óxido sólido simple

mediante un ECC

La Espectroscopia de Impedancia Electroquímica, (EIS), es un

método ampliamente utilizado para entender los fenómenos

que se presentan no sólo en la electrólisis por vapor

39



sobrecalentado, sino, en muchos sistemas en los que no se

puede determinar su impedancia por métodos convencionales.

Este método es utilizado en las celdas de electrólisis como una

herramienta muy útil, tanto al momento de caracterizar los

materiales de los electrodos y del electrolito como al analizar

la influencia de parámetros como la temperatura, voltaje, el

área transversal y la distancia entre electrodos en el

comportamiento de la celda.

De acuerdo con el análisis de Escobedo y Zamora (2006),

Nechache et al. (2014), Orlik. et al. (2012), se pudo encontrar

la topología más adecuada para el sistema, es decir, aquella

que logra representar cada uno de los fenómenos que en él se

presentan. En la ¡Error! No se encuentra el origen de la

referencia. se muestra la configuración resultante para el

análisis de la producción de hidrógeno a partir de vapor

sobrecalentado. Como modelo de circuito eléctrico

equivalente EEC, en la Figura 3 se observa lo siguiente:

Figura 3. Modelo de circuito eléctrico equivalente EEC para la generación

de hidrógeno por electrólisis a partir de vapor sobrecalentado.

La definición de cada elemento se presenta a continuación:

LH inductancia asociada al hardware del sistema (H); RH

resistencia asociada al hardware del sistema (Ω); Rce

resistencia de la interfaz cátodo-electrolito (Ω); Cce

capacitancia de doble capa de la interfaz cátodo-electrolito

(F); Raeresistencia de la interfaz ánodo-electrolito (Ω); Cae

capacitancia de doble capa de la interfaz ánodo-

electrolito (F); Re resistencia del electrolito (Ω) y Ce

capacitancia del electrolito (F).

Se conoce que en una celda de combustible de óxido sólido

SOFC, cuyo comportamiento y modelo EEC es semejante al

de la celda de electrólisis de óxido sólido SOEC, conforme

aumenta la temperatura la impedancia generada por el ánodo,

es menor a la generada por el cátodo el cual tiene mayor

influencia sobre las pérdidas en la celda.

Además, de acuerdo al análisis en el dominio de la frecuencia

hecho en Escobedo y Zamora (2006); se dice que para el

presente trabajo la impedancia del cátodo tiene una influencia

importante sobre las pérdidas de la celda mientras que la

impedancia del ánodo no, razón por la que la impedancia de

este último se considerará despreciable conforme aumenta la

temperatura, lo que quiere decir que Rae~0 y Cae~0.

Al ser Ra y Ca insignificantes con respecto a Rc y Cc y más aún

a mayor temperatura, es evidente entonces que a 1 123,15 (K)

se puede despreciar la impedancia de la interfaz ánodo-

electrolito. Con esto finalmente, en la Figura 4, se muestra el

EEC que se utilizará para la investigación.

Figura 4. Modelo de circuito eléctrico equivalente EEC propuesto para la generación de hidrógeno por electrólisis a partir de vapor sobrecalentado.

2.5. Análisis de electrolito

De 473 (K) a 1 873 (K) el Y2O3 estabilizado con ZrO2 sigue

una ley de Arrhenius no lineal. Para estudiar el

comportamiento de la temperatura en este material es

necesario encajar el rango completo de temperaturas

contenidas en una función teórica que permita entender los

parámetros físicos implicados tales como: entropía, entalpía,

energías contenidas y liberadas (Energía libre de Gibbs), entre

otros.

Conociendo la fórmula empírica de la entalpía ∆H (J/mol)

dependiente de valores locales de temperatura (absoluta)

Cachadiña et al. (1998), se puede afirmar como se muestra en

la Figura 4 que:

∆H(T) = Ef + Eb tanh [β (1

T−

1

Tf

)] ( 4)

Donde, Ef, Eb, β y Tf son parámetros ajustables,

independientes de la temperatura. A la Ecuación 4 se la puede

conocer también como función energía de activación. En

Cachadiña et al. (2001), se ha logrado simplificar una ecuación

que comprende el estudio de la conductividad iónica y

parámetros energéticos junto con la entalpía, usando una

exponencial para aproximar la curva de energías respecto la

temperatura como se muestra en la Figura 5, modelo que ha

sido aceptado para diferentes modelos de materiales

cerámicos.

Figura 5. Datos de temperatura vs entalpía obtenidas de las experiencias

con el 𝒀𝟐𝑶𝟑estabilizado con 𝒁𝒓𝑶𝟐. Fuente: Cachadiña, I., & Solier, J. (1993).

Con lo que, la ley de Arrhenius para la conductividad iónica se

puede escribir como la Ecuación 5.

40



σT

= A(T) exp {−(Ef + Eb tanh [β (

1

T−

1

Tf)])

kT}

( 5)

Que a su vez corresponde a la conductividad total para el, en

donde A (T) es un factor pre-exponencial dependiente de la

temperatura y k es la constante de Boltzmann igual a 1,381 ·10−23(J/K). Si la admitancia compleja total de la celda de

acuerdo con Cachadiña et al. (1995) está representada como en

la Ecuación 6.

Y∗(ω) = G + iωC∗(ω) + Ye∗ ( 6)

Donde Ye∗ es la admitancia en la interfaz electrodo-electrolito,

C∗(ω) representa la respuesta de la masa dieléctrica incluida

la capacitancia geométrica debido a la constante dieléctrica de

alta frecuencia ϵo, i como número imaginario √−1 y G =1

Re

es la conductancia másica (es decir, de la masa activa que

genera la capacitancia).

El modelo de circuito permite separar la conductividad, el

dieléctrico y las contribuciones de la interfase a la admitancia

total. Esta se usa para conocer la capacidad de conducción que

tiene el electrolito. Ahora es necesario conocer los valores de

Re y Ce del modelo descrito para lo cual se tiene que la

respuesta dieléctrica másica, permitividad compleja de

Havriliak-Nagami ϵ∗en (F/m), es como se muestra en la

Ecuación 7.

ϵ∗ =C∗(ω)

Lr ϵo

( 7)

Con Lr como factor que describe la relación entre el área del

electrodo y su espesor. Esta ecuación puede ser representada

por la ecuación de Havriliak-Negami como la Ecuación 8.

ϵ∗ =ϵs − ϵα

[1 + (jω

ωp)1−αe

]

βe+ ϵα

( 8)

Siendo: 0 ≤ αe ≤ 1, 0 ≤ βe ≤ 1, ϵs la permisividad estática

dieléctrica y ωp la frecuencia de pico de pérdida (que depende

de la frecuencia de resonancia) y, como indica su nombre,

representa la frecuencia con la que varían los picos de acuerdo

a su ancho de pulso; esto es así porque el YSZ tiene influencia

directa sobre la frecuencia de resonancia, y que a su vez

depende de la temperatura y otros parámetros como se puede

apreciar, (Cachadiña, Solier, Fatuarte, Sánchez, y Domínguez,

2001), en la Ecuación 9.

ωp = ωpo exp (−Ed

kT) ( 9)

Donde Ed es la energía de activación para la relajación

dieléctrica en (eV) y ωpo un factor pre-exponencial obtenido

de forma gráfica en (Hz).

Conforme a Escobedo y Zamora (2006); Mazloomi et al.

(2012), Macdonald (1992), Cachadiña, et al. (2001) y la

topología de celda planteada anteriormente se pueden

representar a los elementos internos del proceso electrolítico

como una resistencia Re y una capacitancia compleja Ce en

paralelo, cuyos valores expresados en términos de las variables

que inciden en la producción de hidrógeno, de acuerdo con

Cachadiña et al. (2001), como las Ecuaciones 10 y 11 son:

Re =d

σTA ( 10)

Ce = ϵ∗ε0

A

d ( 11)

2.6. Análisis de la interface electrodo-electrolito

2.6.1. Sobre-potencial de activación

Como se dijo en un principio, al existir dentro de la celda

electrolítica dos tipos de electrodos, ánodo y cátodo, se llegó a

tener dos interfaces electrodo-electrolito respectivamente.

Conocer detalladamente los procesos que se llevan a cabo

dentro de la celda es esencial para el estudio de la producción

del gas hidrógeno en un electrolizador a partir de vapor

sobrecalentado de agua. Dado que los electrodos son hechos

de un material metálico (conductor de electrones) y que el

electrolito es un conductor de iones, se produce una separación

y acumulación de cargas en ambos lados de cada interfaz, lo

cual se conoce como el fenómeno de capacitancia de doble

capa; así como también se producen al mismo tiempo, en

dichas interfaces, reacciones electroquímicas (Escobedo y

Zamora, 2006).

Los electrodos poseen una capa catalítica de tal forma que

cuando un electrón cruza esta capa, ya sea para pasar del

electrodo a un ión o viceversa, se dice que se produce la

transferencia de carga, proceso que genera una corriente en el

sistema.

Debe recalcarse, de acuerdo con lo anterior, que tanto si un

electrón pasa del electrodo al electrolito en una interfaz, como

si lo hace del electrolito al electrodo en la otra, la corriente se

produce en el mismo sentido. Este fenómeno se llega a ilustrar

en la ¡Error! No se encuentra el origen de la referencia. en

la que por ejemplo, si los iones al lado del electrolito aceptan

o ceden electrones al electrodo dentro de la reacción, este

intercambio modifica la valencia y por tanto se produce un

flujo de corriente.

41



Figura 6. Transferencia de carga en una interfase. Fuente: Escobedo, E., &

Zamora, L. (Agosto de 2006).

Se puede distinguir dos tipos en el proceso de transferencia de

carga, uno de electronización en el cual un ión acepta un

electrón del electrodo, por lo tanto, su estado de valencia

disminuye en uno o se reduce (reducción), y un proceso de

deselectronización, caso en el que un ión cede un electrón del

electrodo, por lo tanto, su estado de valencia aumenta en uno

o se aumenta (oxidación).

Cuando un ión se está desplazando a través de la doble capa

hasta el electrodo entonces en algún punto de este recorrido se

produce la transferencia de carga; en la ¡Error! No se

encuentra el origen de la referencia. se muestra la variación

de la energía del ión con respecto al recorrido.

Figura 7. Energía Potencial – Distancia en un ión que atraviesa la doble

capa. Fuente: Escobedo, E., & Zamora, L. (Agosto de 2006).

Se puede apreciar entonces que se necesita cierta cantidad de

energía para que se produzca la transferencia de carga, a esta

energía se le conoce como energía de activación ηact, la cual

está representada, de acuerdo con el análisis descrito en

Escobedo y Zamora et al. (2006), por la Ecuación 12.

ηact =2,3RT

αzFlog (

j

jo) ( 12)

Siendo z el número de iones, α el coeficiente de transferencia

de carga, jo densidad de corriente de intercambio, (A/m2), R

la constante universal de los gases 8,3143 (J · K−1 · mol−1) y

F la constante de Faraday 96 485 (C/mol). Dicha energía de

activación, conforme con la ecuación anteriormente descrita,

también puede ser manifestada como un voltaje requerido para

que se produzca el proceso de transferencia de carga o

activación el cual le produce pérdidas al potencial de Nernst o

potencial reversible (Atkins, 2008).

2.6.2. Sobre-potencial de concentración

A más de las pérdidas debido a la energía necesaria que se

tiene que dar a un ión para que ocurra la transferencia de carga

o activación (sobrepotencial de activación), existen también

fenómenos de transferencia de masa que producen pérdidas

adicionales dentro del proceso de electrólisis.

La velocidad de las reacciones químicas que se producen en el

ánodo y en el cátodo, es afectada por la concentración de

suministro presente en los mismos; debido a que no se puede

consumir inmediatamente todo el volumen de suministro,

entonces es necesaria una cantidad de energía para mantener

esa concentración de suministro.

Para estudiar los efectos de transporte de masa en la presente

investigación, se utilizó la aproximación de película de Nernst

o Nernstiana, que afirma que la superficie del electrodo está

cubierta por una película inactiva uniforme, esto produce una

convección suficientemente adecuada para que la

concentración se mantenga uniforme en el exterior, por otra

parte, se asume que la difusión se da en dirección normal a la

superficie de concentración, por lo tanto se puede expresar ∇ci,

gradiente de concentración (mol), como una función lineal,

esto se ilustra en la ¡Error! No se encuentra el origen de la

referencia..

Figura 8. Gradiente de concentración tomado como una función lineal.

Fuente: Escobedo, E., & Zamora, L. (Agosto de 2006).

Análogo al sobrepotencial de activación, la energía que

permite mantener la concentración de suministro se manifiesta

en un voltaje que el sistema toma de sí mismo (y por lo tanto

le produce pérdidas), conocido como sobrepotencial de

concentración ηcon Escobedo y Zamora (2006) como se

muestra en la Ecuación 13.

ηcon = −RT

zFln (1 −

j

jmax

) ( 13)

Con jmax igual a la máxima densidad de corriente (A/m2). De

todo lo hasta aquí expuesto se define finalmente una

resistencia de la interfaz electrodo-electrolito Rce que resulta

42



de la combinación de los potenciales de activación y

concentración, producidos en la interfase, divididos para la

densidad de corriente como se muestra en la Ecuación 14.

Rce =ηact + ηcon

j ( 14)

El fenómeno de capacitancia de doble capa producido en la

interfaz electrodo-electrolito por razón de separación y

acumulación de cargas, explicado anteriormente, puede ser

expresado de acuerdo con Escobedo y Zamora (2006),

Macdonald (1984), Macdonald (1992) y Fernández et al.

(1997) como en la Ecuación 16.

Rce =A [

2,3RT

αzFlog10 (

I

A jo) −

RT

zFln (1 −

I

A jmax)]

I

( 15)

Siendo Yo una constante que depende de la temperatura

medida en [FsnCPE −1]; y nCPE una constante que depende de

la geometría de la superficie del electrodo y de la temperatura,

que conforme el análisis realizado en Jorcin et al. (2005),

asume un valor de 1 por lo que Cce pasa a ser un elemento

puramente capacitivo.

2.6.3. Sobre-potencial óhmico

Los electrodos utilizados están hechos de un material

conductor cubierto con un electrocatalizador (material que

incrementa la velocidad de la reacción que participa en una

reacción electroquímica); los cuales presentan una pequeña

oposición al flujo de electrones que, como en cualquier

conductor, dependerá de la longitud, área y resistividad a su

vez dependiente de la temperatura.

Se definirá a esta resistencia inherente, a cualquier conductor,

como resistencia electrónica Re− Para la determinación de una

resistencia óhmica total se consideran otros elementos de la

celda tales como: placas colectoras, capa de difusión y cables

de conexión, que en conjunto dan como resultado una

resistencia de celda denominada RCELDA (Escobedo y Zamora,

2006), (Mazloomi, et al., 2012).

Así, se define una resistencia total del hardware RH como la

de la Ecuación 17.

Cce = Yo(ω)nCPE ( 16)

Por lo que el sobrepotencial óhmico (ηohm) estará determinado

por la ecuación 18.

RH = Re− + RCELDA ( 17)

Por lo que el sobrepotencial óhmico (ηohm) estará determinado

por la Ecuación 18.

ηohm = jRH ( 18)

Es necesario concluir en este punto, que en la investigación se

presentaron tres clases de pérdidas, a saber, una pérdida debida

a la energía necesaria para transferir carga, una pérdida debida

a la energía necesaria para transferir masa y por último una

pérdida óhmica debida a la resistencia propia del hardware del

sistema.

En este punto se resume a continuación, en la Tabla 1, los

valores iniciales óptimos de las variables descritas

previamente y que inciden directamente en la producción de

hidrógeno. Se asume dichos valores conforme a Brisse et al.

(2008), Cachadiña et al. (2001), Mazloomi et al. (2012) y

Nechache et al. (2014).

Tabla 1. Rangos y valores iniciales óptimos asumidos para las variables que

inciden en la producción de hidrógeno.

Variable Descripción Unidad Rango Valor

d Espesor

(separación

electrodos).

(m) [0 ; 0,01] 5 · 10−3

A Área

transversal

electrodos.

(m2) [0 ; 4,5 ·10−3]

1,6· 10−3

Vs Voltaje de

entrada.

(V) [0 ; 2] 1,33

D Ciclo de trabajo (Duty Cycle).

(%) [0 ; 100] 50

T Temperatura (K) [773,15 ; 1

873,15]

1 123,15

2.7. Implementación del modelo de ECC.

La implementación del modelo de EEC de la celda electrolítica

se lo hizo utilizando la herramienta PowerSystemsToolbox de

Matlab.

Las celdas para la electrólisis se ilustran esquemáticamente en

la ¡Error! No se encuentra el origen de la referencia. para

cuando se utiliza una fuente DC y para cuando se utiliza una

fuente de voltaje aplicada en forma de pulsos a frecuencia de

resonancia.

Figura 9. Implementación del modelo de EEC utilizando una fuente DC.

En este modelo, la inductancia y resistencia de cableado y

conexiones así como el equivalente eléctrico del cátodo y del

electrolito de la celda asumiendo que la impedancia del ánodo

es despreciable.

.

Para esta última implementación, esquematizada en la ¡Error!

No se encuentra el origen de la referencia., se obtiene una

fuente de voltaje en forma de pulsos utilizando para ello un

MOSFET de ultra baja resistencia, con respecto a los

componentes del sistema, que es colocado entre el cátodo y

tierra (Mazloomi, et al. 2012). Este dispositivo funciona como

un interruptor que conduce y corta el paso de la corriente a

través de la celda en un rango determinado de frecuencias

1(Hz) a 1(MHz).

43



Figura 10. Implementación del modelo de EEC aplicando una fuente de

voltaje en forma de pulsos a frecuencia de resonancia.

Las características ideales seleccionadas para este elemento se

destacan a continuación: resistencia FET: 1 · 10−4(Ω);

resistencia interna del diodo: 1 · 10−4(Ω); resistencia

Snubber: 1 · 105(Ω); capacitancia Snubber: ∞; inductancia

interna del diodo: 0 (H). Se utiliza finalmente, en esta

configuración, un generador de onda para la conducción del

MOSFET de potencia. Su incorporación permite generar una

forma de onda cuadrada con un ciclo de trabajo, que para este

caso convenientemente corresponde al 50% como se muestra

en la ¡Error! No se encuentra el origen de la referencia.. El

rango de frecuencia de salida de este dispositivo podría ser

sintonizado entre 1(Hz) y 1 (Mhz). La amplitud de la señal

generada se fijó en -0,5 (V) para el apagado y 18 (V) para el

encendido, valores usados con el fin de garantizar las mejores

condiciones de funcionamiento para el interruptor

semiconductor (Mazloomi et al., 2012).

Figura 11. Corriente de celda vs [%] de ciclo de trabajo.

2.8. Caudal de hidrógeno

Según las leyes de Faraday para electrólisis (López et al.

2010), se define la función para conocer la producción de

hidrógeno sobre la cual intervienen todas las funciones del

sistema, con masa molar, M = 4,032 (g/mol) y número de

valencia (electrones por ion) n = 4 como se muestra en la

Ecuación 19.

∅H2 =M

nF∫ I(t)dt

t

0

= (1,04466 · 10−5) · ∫ I(t) dtt

0

( 19)

En donde ∅H2(g/s) es el caudal de hidrógeno (flujo másico)

que depende de la ley de Faraday con la corriente I (A) en el

tiempo de operación t (s).

Finalmente, el desarrollo de las expresiones matemáticas, a

partir del modelo por EEC, para la obtención de la resistencia

de la interfaz cátodo-electrolito, resistencia del electrolito,

capacitancia del electrolito, capacitancia de la interfaz cátodo

electrolito, inductancia asociada al hardware del sistema,

frecuencia de resonancia y caudal de hidrógeno en función de

las variables temperatura, espesor y área de los electrodos,

fuente de voltaje y ciclo de trabajo se presenta a continuación:

Partiendo de la ley de Arrhenius para la conductividad iónica

(Cachadiña et al. 2001) presentada en la Ecuación 5, la

ecuación de Havriliak-Negami, expuesta por la Ecuación 8, así

como también la frecuencia de pico de pérdida mostrada en la

Ecuación 9, se puede obtener la resistencia de electrolito y la

capacitancia de electrolito de la siguiente forma; combinando

las Ecuaciones 5 y 10 se obtiene la resistencia del electrolito

como se indica en la Ecuación 20.

Re =d

A(T) exp {−(Ef+Eb tanh[β(

1

T−

1

Tf)])

kT} · A

( 20)

La combinación de las Ecuaciones 8 y 11 permiten encontrar

la capacitancia del electrolito (Ce), mostrada en la Ecuación

21.

Ce =

[

ϵs − ϵα

[1 + (iω

ωpo exp(−EdkT

))

1−αe

]

βe+ ϵα

]

· ε0

·A

d

( 21)

Donde: ε0 = 8,854 · 10−12 (F/m) siendo ésta la

permitividad del vacío.

Conforme Mazloomi et al., 2012) y Savova et al. (1987) se

asume una inductancia de cables y conexiones Lcc mostrada

en la Ecuación 22.

Lcc = LH = 1 (μH) ( 22)

Las ecuaciones de la resistencia de la interfaz cátodo-

electrolito (Rce) y de la capacitancia de la interfaz cátodo-

electrolito (Cce) ya se han tratado previamente en las

Ecuaciones 15 y 16 respectivamente.

De los valores tomados dentro de rangos de operación

considerables (Yong y Moon, 2010), (Cachadiña et al.,2001),

(Nechache et al., 2014), se seleccionan los presentados en las

Tablas 2, 3, 4, 5, 6 y 7 para las constantes asociadas:

44



Tabla 2. Rangos y valores asumidos para los parámetros constantes de la Ecuación 8.

𝝐∗ Descripción Unidad Rango Valor

αe Constantes de

simetría.

[0 ; 1] 0,5

βe

𝜖𝑠 Permitividad estática

dieléctrica. (F/m) 3,529478

· 108

3,529478 · 108

𝜖∞ Constante dieléctrica de alta frecuencia.

(F/m) 176,471 176,471

𝜔𝑝 Frecuencia del pico

perdido.

(Hz) Depende de

otras constantes.

Tabla 3. Rangos y valores asumidos para los parámetros constantes de la

Ecuación 9.

𝝎𝒑 Descripción Unidad Rango Valor

𝐸𝑑 Energía de activación para

la relajación dieléctrica (J) 1,0748

· 10−20

1,0748· 10−20

K Constante de Botlzmann (J/K) 1,3806· 10−23

1,3806· 10−23

T Temperatura variable del

proceso (℃) variable

𝜔𝑝𝑜 Factor pre-exponencial (Hz) 10 · 106


Ecuación 5.

𝝈𝑻 Descripción Unidad Rango Valor

𝐸𝑓

Parámetros variables

independientes de la

temperatura.

(eV) 1,01

𝐸𝑏 (eV) 0,29

𝛽 (K) 3 278

𝑇𝑓 (K) 1 328

T Temperatura del proceso. (K)

Depende de

otras constantes.


Ecuación 12.

𝛈𝐚𝐜𝐭 Descripción Unidad Rango Valor

R Constante universal de

los gases. (J · K−1

· mol) 8,3144 8,3144

T Temperatura de

operación absoluta. (K) 1 073,15

𝛼 Coeficiente de

transferencia de carga. [0 ; 1] 0,5

𝑧

Número de electrones implicados en la

reacción.

4

F Constante de Faraday. (C/mol) 96 485 96 485

j Densidad de corriente del electrodo.

(A/m2) Variable

jo Densidad de corriente

de intercambio (A/m2) 0,35


Ecuación 13.

𝛈𝐜𝐨𝐧 Descripción Unidad Rango Valor

R Constante universal de los gases.

(J · K−1

· mol) 8,3144 8,3144

T Temperatura de operación absoluta.

(K) 1 073,15

𝑧

Número de electrones

implicados en la reacción.

4

F Constante de Faraday. (C/mol) 96 485 96 485

jmax Densidad de corriente

máxima (A/m2) 2


Ecuación 16. 𝐂𝐜𝐞 Descripción Unidad Rango Valor

Yo Constante dependiente de

la temperatura. (mF) 1

nCPE Constante geométrica de

superficie del electrodo. 1

En cuanto a la resistencia de hardware, RH, ha resultado poco

práctico suponer que esta sea mayor que 1 ohmio de acuerdo a

Nechache et al. (2014) y Mazloomi et al. (2012), razón por la

que se optó por darle el valor presentado en la Ecuación 23.

RH = 0,1 (Ω) ( 23)

Para encontrar la frecuencia de resonancia en función de los

parámetros conocidos es importante considerar inicialmente la

Ecuación 24.

ZTOTAL = ZH + Zce + Ze ( 24)

Tras obtener las impedancias individuales del circuito EEC

mostrado en la Figura 4, al determinar la impedancia total se

obtiene la expresión dada en la Ecuación 25.

ZTOTAL = (RH +Rce

1 + ω2Cce2Rce

2

+Re

1 + ω2Ce2Re

2)

+ i (ωLH

−ωCceRce

2

1 + ω2Cce2Rce

2

−ωCeRe

2

1 + ω2Ce2Re

2)

( 25)

Sea la parte imaginaria de la impedancia total de la forma

expresada por la Ecuación 26.

𝐼𝑚{𝑍𝑇𝑂𝑇𝐴𝐿} = ωLH −ωCceRce

2

1 + ω2Cce2Rce

2

−ωCeRe

2

1 + ω2Ce2Re

2

( 26)

Para obtener la frecuencia de resonancia para esta

configuración de dos etapas se procede como se muestra en la

Ecuación 27.

ωLH −ωCceRce

2

1 + ω2Cce2Rce

2 −ωCeRe

2

1 + ω2Ce2Re

2

= 0

( 27)

Siendo entonces la Ecuación 28 la que nos permite obtener la

frecuencia de resonancia del sistema.

45



ωLH(1 + ω2Cce2Rce

2)(1 + ω2Ce2Re

2)

− ωCceRce2(1

+ ω2Ce2Re

2)

− ωCeRe2(1

+ ω2Cce2Rce

2) = 0

( 28)

2.8.1. Caudal de hidrógeno respecto de la temperatura

Dado el valor de Re mostrado en la Ecuación 29:

Re = (1,5625 · 10−2) exp{𝑋1} ( 29)

Dónde: 𝑋1 = −(1,01−0,29 tanh[3 278(

1

T−

1

1 328)])

4,4·10−4·T

Del modelo EEC se tendría la Ecuación 30:

I(s)

=Vs

0,1 + s(1 · 10−6) +Rce

1+sCceRce+

Re

1+sCeRe

( 30)

De la Tabla 1 se tiene que Vs = 1,33 (V).

Al reemplazar los valores de Rce, Cce, Re se obtiene la

Ecuación 31:

I(s) =1,33

X + Y ( 31)

Empleando las Ecuaciones 31, 32, 33, 34 y 35 se obtiene el

valor del caudal de hidrógeno en función de la temperatura,

presentado en la Ecuación 36.

X = 0,1 + s · (1 · 10−6)

+X2

1 + s · 1 · 10−6 · X2

( 32)

X2

=105 856 · 10−3 · Tlog10 (

I

5,6) − 3,446919 · 10−4 · ln (1 −

I

32)

I

(

33

)

Y =(1,5625 · 10−2) exp{X3}

1 + s · 1 · 10−6 · (1,5625 · 10−2) exp{X3}

(

34

)

𝑋3

= −(1,01 − 0,29 tanh [3 278 (

1

T−

1

1 328)])

4,4 · 10−4 · T

( 35)

∅H2(s) = (1,04466 · 10−5) ·I(s)

s ( 36)

2.8.2. Caudal de hidrógeno respecto del área transversal de

los electrodos

Considerando las Ecuaciones 30, 37, 38, 39, 40, 41 y 42 se

puede obtener el valor del caudal de hidrógeno en función del

área expresado por la Ecuación 43. Considerando un área

variante:

Rce

=A [0,113 log10 (

I

0,35A) − 2,41963 · 10−2 ln (1 −

I

2A)]

I

(

37

)

Re =0,5

0,3319 ∗ A=

1,5065

A ( 38)

I(s) =1,33

0,1 + s · (1 · 10−6) + G + H ( 39)

𝑋4

=A [0,113 log10 (

I

0,35A) − 2,41963 ∗ 10−2 ln (1 −

I

2A)]

I

(

40

)

G =X4

1 + s · 1 · 10−6 · X4

(

41

)

H =

1,5065

A

1 + s · 1 · 10−6 ·1,5065

A

( 42)

∅H2(s) = (1,04466 · 10−5)I(s)

s ( 43)

2.8.3. Caudal de hidrógeno respecto del voltaje de entrada

Considerando la Ecuación 44, Ley de Ohm:

𝑉𝑠 = 𝑍 · 𝐼 ( 44)

Se tiene el valor de caudal de hidrógeno en función del voltaje,

como se indica en la Ecuación 45.

∅H2(s)

VS(s)=

1,04466 · 10−5

Z(s)∗

1

s ( 45)

Donde la impedancia puede expresarse mediante las

Ecuaciones 46, 47 y 48 como:

Z(s) = 0,1 + s · 1 · 10−6 + X5 + X6 ( 46)

𝑋5 =6,81488 · 10−2

1 + s · 6,81488 · 10−5 ( 47)

𝑋6 =0,15

1 + s · 1,5 · 10−4 ( 48)

Entonces el caudal de hidrógeno en función del voltaje de

entrada se puede expresar mediante la Ecuación 49 como:

∅H2(s) =1,04466 ∗ 10−5 · Vs

0,1 + s · 1 · 10−6 + X5 + X6·1

s ( 49)

2.8.4. Caudal de hidrógeno respecto del ancho de pulso

Considerando las Ecuaciones 50, 51 y 52 para posteriormente

reemplazarlas en la Ecuación 53, con B representada mediante

46



la Ecuación 54, se puede entonces obtener el caudal de

hidrógeno como se muestra en la Ecuación 57 partiendo de la

Ecuación 55 y empleando la Ecuación 56.

𝑋7 =1,764710 · 108

[1 + (2πf

5·106)0,8

]0,5 + 1,764710 · 108

( 50)

Ce = [X7] · 2,833 · 10−12 ( 51)

X8 =6,81488 · 10−2

1 + iω · (6,81488 · 10−5)

( 52)

Z(iω) = 0,1 + iω · (1 · 10−6) + X8 + B ( 53)

B =0,15

1 + iω · 4,2499 · 10−13 · [X7] ( 54)

∅H2(iω) = (1,04466 · 10−5) ·Vs(iω)

Z(iω)·

1

iω ( 55)

𝑋9 =1,33

0,1 + iω(1 · 10−6) + X8 + B ( 56)

∅H2(iω) = (1,04466 ∗ 10−5) · X9 ·1

iω ( 57)

2.8.5. Caudal de hidrógeno respecto del espesor de electrolito

o separación de los electrodos.

Cuando el espesor es variante aplicamos las Ecuaciones 58 y

59.

Re =d

5,3104 ( 58)

Ce = Ce =5 · 10−4

d ( 59)

La impedancia se puede expresar mediante la Ecuación 61,

donde X5 y X10 están representados por las Ecuaciones 47 y

60 respectivamente.

𝑋10 =

d

5,3104

1 + s · 9,415 · 10−5 ( 60)

Z(s) = 0,1 + s · (1 · 10−6) + X5 + X10 ( 61)

Partiendo de la Ecuación 55, tomando en cuenta que s = iω,

se puede calcular el caudal de hidrógeno respecto del espesor,

Ecuación 63, empleando la Ecuación 62.

𝑋11 =1,33

0,1 + s · (1 · 10−6) + X5 + X10 ( 62)

∅H2(s) = (1,04466 · 10−5) · X11 ·1

s ( 63)

Finalmente, los valores de resistencias y capacitancias

obtenidos se muestran en la Tabla 8.

Tabla 8. Valores obtenidos de resistencias y capacitancias del EEC. 𝐑𝐞 𝟎. 𝟏𝟓 (𝛀)

Rce 6,81488 ∗ 10−2(Ω) Ce 1 (mF) Cce 1 (mF)

3. RESULTADOS Y DISCUSIÓN

Se presenta a continuación de forma gráfica la influencia de

algunas de las variables que influyen en la producción de

hidrógeno para luego realizar su respectivo análisis.

Para esto, se han seleccionado los valores comunes de éstas

variables y otros parámetros en base a las referencias

bibliográficas usadas a lo largo de este documento, de tal

forma que al analizar la influencia de una de ellas en la

producción de hidrógeno el resto permanezca sin modificarse.

3.1. Aplicación de una señal de excitación tipo DC

En la ¡Error! No se encuentra el origen de la referencia., la

respuesta del sistema de electrólisis por vapor sobrecalentado

para una señal de excitación tipo DC de valor 1,33 (V).

Figura 12. Respuesta del sistema para una señal de excitación tipo DC.

La respuesta del sistema en el tiempo, es constante, alcanzando

un valor de corriente de aproximadamente 5,12 (A), como se

observa en la Figura 12.

3.2. Aplicación de una señal de excitación en forma de pulsos

a frecuencia de resonancia

Obtenida la frecuencia de resonancia: fR = 6,82131 (kHz), la

respuesta del sistema ante una señal de pulsos con frecuencia

fRse observa en la Figura 13.

47



Figura 13. Respuesta del sistema para una señal en forma de pulsos a

frecuencia de resonancia.

En la Figura 13 se observa que la oscilación de la señal de

respuesta del sistema es mínima y por lo tanto se la puede

considerar como una señal DC con una amplitud de corriente

de 5,5449 (A).

Esto permite demostrar que para una señal de excitación en

forma de pulsos a frecuencia de resonancia aumenta la

eficiencia del sistema en un 8.2988% como se muestra en la

Ecuación 64.

Aumento de corriente (%)

=5,5449 − 5,12

5,12· 100 %

= 8,2988

( 64)

Debido al aumento de eficiencia de la celda, se procede a

trabajar con una señal en forma de pulsos a frecuencia de

resonancia, lo cual incide de forma positiva en la producción

del gas hidrógeno a partir de vapor sobrecalentado de agua.

3.3. Influencia de la temperatura en la producción de

hidrógeno

Figura 14. Influencia de la temperatura en la corriente del sistema.

Figura 15. Caudal de gas H2 vs Temperatura.

Análisis de resultados de las Figuras 14 y 15:

El comportamiento de la celda electrolítica simple

con respecto a la temperatura, tal y como se puede

apreciar en ambas figuras, es del tipo no lineal.

La eficiencia máxima para este análisis se logra

cuando la temperatura es aproximadamente

1123,15 (K).

3.4. Influencia del área transversal de los electrodos en la

producción de hidrógeno.

Figura 16. Influencia del área transversal en la corriente del sistema.

Figura 17. Caudal de gas H2 vs área transversal de los electrodos [m2].


El comportamiento de la celda con respecto al área transversal

se puede aproximar a una recta con pendiente positiva, esto

48



implica que, a mayor área, mayor corriente. Sin embargo,

aumentar el área transversal de la celda implica un aumento de

costos. Ante aquello, se presenta una alternativa que consiste

en aumentar el área transversal de la celda dependiendo la

aplicación en la que vaya a ser utilizada. Con esto, se calcula

el área que se puede ocupar maximizando el rendimiento y

evitando las pérdidas tanto monetarias como de recursos.

3.5. Influencia del espesor del electrolito (equivalente a la

separación de los electrodos) en la producción de hidrógeno.

Figura 18. Influencia de la separación de los electrodos en la corriente de

la celda.

Figura 19. Caudal de gas H2 vs voltaje de alimentación.


La caída de la corriente en una celda aumenta a

medida que crece la separación de los electrodos, lo

que sugiere que para que la electrólisis maximice su

producción, las placas deben encontrarse muy cerca

una de otra. Esto debido al hecho de que para la

corriente eléctrica es más fácil y ventajoso atravesar

distancias cortas.

Se recomienda finalmente, tener cuidado en cuanto a

la separación electrodo-electrolito se refiere, pues a

pesar de que estos deben encontrarse cerca, no

pueden estar en contacto, entonces, es necesario usar

instrumentos de precisión para acercar dichos

materiales.

3.6. Influencia del voltaje de entrada en la producción de

hidrógeno

Figura 20. Influencia del voltaje de entrada en la corriente de la celda.

Figura 21. Caudal de gas H2 vs separación de los electrodos.


A medida que aumenta el voltaje de alimentación de

la celda, aumenta la corriente y, por lo tanto, aumenta

el caudal de hidrógeno.

Para aumentar la eficiencia del sistema, se debe

aumentar el voltaje de alimentación del mismo, sin

embargo, la pérdida de energía que se produce debido

a la ley de Joule, es un limitante, que produce pérdida

de eficiencia. Por encima de cierto valor del voltaje

de alimentación la pérdida de eficiencia aumenta

progresivamente. Por lo tanto, el voltaje de una celda

debe oscilar entre 0,8 y 2,6 (V).

4. CONCLUSIONES

Se ha logrado cumplir con el objetivo de la presente

investigación, esto es, modelar la producción de hidrógeno en

una celda electrolítica simple usando como entrada vapor

sobrecalentado de agua con un consumo de energía eficiente,

alcanzada por medio del uso de una señal de excitación en

forma de pulsos a frecuencia de resonancia; con el fin de

permitir el diseño y dimensionamiento de equipos de

electrólisis.

El proceso de producción de hidrógeno mediante electrólisis

con vapor sobrecalentado puede ser representado mediante un

sencillo circuito eléctrico, en el cual, los valores de los

49



elementos dependen de los fenómenos electroquímicos que se

presentan en el sistema.

La representación del proceso mediante un sistema eléctrico

facilita significativamente su análisis, puesto que para su

comprensión sólo es necesaria la teoría de circuitos.

Algunos de los fenómenos electroquímicos que se presentan

en el proceso pueden ser representados de forma aproximada

mediante elementos eléctricos no ideales como el capacitor de

doble capa.

La eficiencia del sistema aumenta aplicando una señal de

excitación de pulsos a frecuencia de resonancia, esto se debe a

que a esta frecuencia la parte reactiva de la impedancia se

anula y se vuelve puramente resistiva.

La variación de la eficiencia del sistema con respecto a la

temperatura es no lineal y presenta su máximo valor en el

rango de 773 a 1 873 (K).

La eficiencia del proceso varía de forma aproximadamente

lineal respecto el área transversal de los electrodos, esto

normalmente implicaría que para un área suficientemente

grande la eficiencia del proceso será del 100(%), sin embargo,

el aumento del área implica costos, tanto por la cantidad de

material como por el tamaño del equipo necesario para

despresurizar la celda, lo que limita el valor de la misma.

La eficiencia del proceso varía de forma no lineal e inversa,

aproximadamente logarítmica, con respecto a la distancia de

separación de los electrodos, esto se debe a que la disminución

de la distancia entre electrodos, reduce el valor de la resistencia

del electrolito.

La eficiencia del proceso con respecto a la amplitud de la señal

de excitación del mismo es máxima para valores de entre 0,8

y 2,66 (V), para valores mayores se pierde energía debido a la

ley de Joule.

El caudal de hidrógeno producido varía de forma lineal con la

corriente del mismo, esta variación está dada por la ley de

Faraday.

Finalmente, la obtención matemática de la función de

transferencia o a su vez del modelo en espacio de estados de

una celda electrolítica simple, podría ser un tema de mucho

interés para futuras investigaciones a fin de dar paso a la

aplicación de una teoría de control clásica o moderna que

permita elevar aún más, de forma eficiente, la producción de

hidrógeno a partir de vapor sobrecalentado de agua. Además,

la construcción de una celda electrolítica simple que trabaje a

alta temperatura sería conveniente a fin de analizar y comparar

los resultados experimentales con los resultados aproximados

obtenidos en el presente trabajo.

REFERENCIAS

Alvarado, J. (2013). Materiales para ánodos, cátodos y electrolitos utilizados

en celdas de combustible de óxido sólido (SOFC). Revista Mexicana de

Física, 59, 66-87.

Atkins, P. D. (2008). Química Física. Madrid, España: Médica Panamericana.

Brisse, A., Schefold, J., & Zahid, M. (2008). High temperature water

electrolysis in solid oxide cells. International journal of hydrogen

energy, 33, 5375-5382. doi:10.1016/j.ijhydene.2008.07.120

Cachadiña, I., Solier, D., & Domínguez, A. (1998). Comportamiento eléctrico de monocristales 4.7 mol % Y-PSZ. Boletín de la Sociedad de Cerámica

y Vidrio, 37(3), 238-242. doi:Obtenido de la base de datos idus, Depósito

de Investigación Universidad de Sevilla. Cachadiña, I., Solier, D., Fatuarte, I., Sánchez, F., & Domínguez, A. (1995).

Espectroscopía de impedancia de monocristales 4.7 mol % Y-PSZ.

Boletín de la Sociedad Española de Cerámica y Vidrio, 34(4), 391-394. doi:Obtenido de la base de datos idus, Depósito de Investigación

Universidad de Sevilla.

Cachadiña, I., Solier, D., Fatuarte, I., Sánchez, F., & Domínguez, A. (2001). Electrical properties of ZrO2-12 mol % Y2O3 single crystals. Third

Euro Ceramics Properties of Ceramics, 2, 259-264

Escobedo, E., & Zamora, L. (Agosto de 2006). Modelo dinámico de celdas de combustible. Tesis para obtener el grado de maestro en ciencias en

ingeniería mecatrónica, cenidet. Obtenido de

http://www.cenidet.edu.mx/subplan/biblio/seleccion/Tesis/MK%20Enr

ique%20Escobedo%20Hernandez%202006.pdf.

Fernández, D., Goodwin, A., Lemmon, E., Levelt, J., & Williams, R. (1997).

A Formulation for the Static Permittivity of Water and Steam at Temperatures from 238 K to 873 K at Pressures up to 1200 MPa,

Including Derivatives and Debye–Hückel Coefficients. Journal of

Physical and Chemical Reference Data, 26(4), 1125-1166. doi:10.1063/1.555997

Jorcin, J. B., Orazem, M., Pébère, N., & Tribollet, B. (2005). CPE analysis by local electrochemical impedance spectroscopy. Electrochimica Acta,

51(8-9), 1473-1479. doi:10.1016/j.electacta.2005.02.128

López, J. (2010). Modelo dinámico de un electrolizador alcalino. Universidad de Sevilla, 1-110. Obtenido de

http://bibing.us.es/proyectos/abreproy/4703 (Junio,2018)

Macdonald, J. (1984). Note on the parameterization of the constant-phase admittance element. Solid State Ionics, 13(2), 147-149.

doi:10.1016/0167-2738(84)90049-3

Macdonald, J. (1992). Impedance Spectroscopy. Annals of Biomedical Engineering, 20, 289-305. Obtenido de

http://jrossmacdonald.com/jrm/wp-

content/uploads/187ImpSpectroscopy.pdf (Junio,2018) Mazloomi, K., Sulaiman, N., & Moayedi, H. (2012). An Investigation into the

Electrical Impedance of Water Electrolysis Cells – With a View to

Saving Energy. International Journal of Electrochemical Science, 7, 3466-3481. doi:Obtenido de la base de datos ResearchGate

Mazloomi, K., Sulaiman, N., & Moayedi, H. (2012). Electrical Efficiency of

Electrolytic Hydrogen Production. International Journal of Electrochemical Science, 7, 3314-3326. doi:Obtenido de la base de datos

ResearchGate

Nechache, A., Cassir, M., & Ringuedé, A. (2014). Solid oxide electrolysis cell analysis by means of electrochemical impedance spectroscopy: A

review. Journal of Power Sources, 258, 164-181.

doi:10.1016/j.jpowsour.2014.01.110 Orlik, M. (2012). Application of Impedance Spectroscopy to Electrochemical

Instabilities. Self-Organization in Electrochemical Systems, 1, 111-195.

doi:10.1007/978-3-642-27673-6_3

Pierre, J., & Capitaine, A. (2006). Hydrogen Production by High Temperature

Electrolisis of Water Vapour and Nuclear Reactors. 16 World Hydrogen

Energy Conference, France, 38, 1-16. doi:Obtenido de la base de datos IAEA-INIS

Savova, B., & Stoynov, Z. (1987). Analysis of the inductance influence on the

measured electrochemical impedance. Journal of Applied Electrochemistry, 17(6), 1150-1158. doi:10.1007/BF01023598

Yong, B., & Moon, S. (2010). Electrochemical Impedance Spectroscopy.

Annual Review of Analytical Chemistry, 3(1), 207-229. doi:10.1146/annurev.anchem.012809.102211 • Source: PubMed

Zoulias, E., Varkaraki, E., Lymberopoulos, N., Christodoulou, C., &

Karagiorgis, G. (2004). A review on water electrolysis. Centre for Renewable Energy Sources (CRES), 1-18. Obtenido de

https://www.researchgate.net/publication/284618929_A_Review_on_

Water_Electrolysis (Junio,2018)

50



BIOGRAFÍAS

Edwin René Aguilar Jaramillo.

Ingeniero en Electrónica Automatización y

Control por la Universidad de las Fuerzas

Armadas ESPE, Ecuador, Magíster en

Diseño, Producción y Automatización

Industrial, en la Escuela Politécnica

Nacional, Ecuador. Ha participado en

proyectos de ingeniería en las áreas de

Control e Instrumentación y posee

investigaciones realizadas en los campos

de modelamiento de sistemas y diseño de

sistemas de control empleando técnicas de control tales como

espacio de estados, LMIs, teoría de Lyapunov y normas de

sistemas H_2 y H_∞. Actualmente es docente del

Departamento de Eléctrica y Electrónica en la Universidad de

las Fuerzas Armadas ESPE y está participando en el proyecto

de investigación “Construcción de un simulador de

desorientación espacial de 4GDL” para la Escuela de Aviación

de la Armada del Ecuador.

Álvaro Gonzalo Xavier Aguinaga

Barragán. Ingeniero mecánico por la

Escuela Politécnica Nacional,

Ecuador, Máster en Tecnologías de la

Información en la Fabricación en la

Universidad Politécnica de Madrid,

España, Doctor Ph.D. en Ingeniería

Mecánica en Politécnica Bialostocka,

República de Polonia. Ha sido

participe de varios proyectos de

ingeniería y consultorías. Posee varias investigaciones

realizadas en ingeniería de software, lenguajes de

programación, sistemas operativos, Smart Grid, y en otros

campos de la ingeniería electrónica, automatización y control

e ingeniería mecánica. Actualmente labora como docente,

investigador y Director de la Escuela de Formación de

Tecnólogos de la Escuela Politécnica Nacional.

51

Modelación Numérica del Flujo Rasante en una Rápida Escalonada Aplicando la Dinámica de Fluidos Computacional (CFD) Mediante el Uso de

Flow-3D





11. INTRODUCCIÓN

El estudio de una rápida permite determinar su eficiencia como

estructura de disipación de energía del flujo en proyectos como

el caso de una central hidroeléctrica o la descarga de agua

lluvia hacia un cuerpo receptor, en el cual se requiere vencer

grandes desniveles topográficos. Las investigaciones en la

actualidad se centran en discusiones del flujo aireado en la

región donde se presenta el flujo rasante. La rápida escalonada

presenta un mayor porcentaje de disipación de energía en

comparación con la rápida lisa lo cual es ventajoso ya que se

reduce la profundidad y tamaño del cuenco disipador ubicado

al pie de dicha estructura (Castro M., 2015).

El análisis de la rápida escalonada resulta complejo debido a

que se tiene diferentes regímenes de flujo y regiones a lo largo

de toda la estructura. Las características del flujo en la rápida

escalonada actualmente son abordadas mediante el uso de

relaciones empíricas y modelos físicos. Una herramienta

complementaria de la modelación física hidráulica es la

modelación numérica, la cual permite reducir tiempo y costo

en la fase de modificaciones del diseño original de las

estructuras hidráulicas, dado que un modelo numérico

calibrado con el modelo físico permite desarrollar las posibles

modificaciones al modelo hidráulico en forma eficiente y

oportuna. En una rápida escalonada existen tres tipos de

vertido que son: flujo de escalón en escalón, flujo en transición

y flujo rasante (Khatsuria, R.M., 2005).

En una rápida escalonada, el flujo rasante se produce cuando

el flujo de agua llena todas las cavidades de los escalones

formándose un fondo virtual. El flujo que circula por encima

Modelación Numérica del Flujo Rasante en una Rápida Escalonada

Aplicando la Dinámica de Fluidos Computacional (CFD) Mediante

el Uso de Flow-3D

Casa E. 1; Hidalgo X. 1; Castro M. 1; Ortega P1; Vera P1

1CIERHI- Escuela Politécnica Nacional, Facultad de Ingeniería Civil y Ambiental, Quito, Ecuador

Resumen: El presente estudio tiene como objetivo principal desarrollar la modelación numérica del flujo rasante en

una rápida escalonada aplicando el paquete comercial FLOW-3D. Actualmente el diseño de este tipo de estructuras

se realiza con el uso de expresiones empíricas obtenidas con base en la modelación física y estudios complementarios

en la modelación numérica del flujo sobre la rápida escalonada con apoyo de un código CFD. Con el modelo numérico

se busca estimar la velocidad del flujo en la región uniforme, y el coeficiente de fricción para cuatro caudales de

operación de la rápida escalonada (ϴ=45º, Hd=4.61m). La representación de un flujo autoaireado es complejo, por lo

que el programa aproxima la solución con ciertas limitaciones, para ello emplea el modelo air entrainment, el modelo

drift flux y un modelo de turbulencia k-ԑ RNG. Los resultados obtenidos con la modelación numérica y la modelación

física acerca del inicio de autoaireación natural del flujo y la profundidad del flujo bifásico en la región uniforme

presentan desviaciones de hasta el 10% debido a que el flujo es altamente turbulento.

Palabras clave: Rápida Escalonada, Flujo Rasante, Introducción de Aire, Modelación Numérica, FLOW-3D.

Numerical Modeling of Flush Flow in a Rapid Step Applying

Computational Fluid Dynamics (CFD) Using Flow-3D

Abstract: The main objective of this project is to develop the numerical modeling of the skimming flow in a stepped

spillway using FLOW-3D. The design of these structures is based on the use of empirical expressions obtained from

physical modeling and complementary studies in the numerical modeling of flow over the stepped spillway with

support of CFD code. The numerical model is used to estimate the flow velocity in the uniform region and the friction

coefficient of the stepped spillway (ϴ = 45º, Hd=4.61m). The representation of auto aeration a flow is complex, so

the program approximates the solution with certain limitations, using an air entrainment model; drift flux model and

turbulence model k-ԑ RNG. The results obtained with numerical modeling and physical modeling at the beginning of

natural auto aeration of flow and depth of the biphasic flow in the uniform region presents deviations above to 10%

perhaps the flow is highly turbulent.

Keywords: Stepped spillway, skimming flow, air entrainment, drift flux, numerical modeling, FLOW-3D.

53



del fondo virtual se conoce como flujo rasante completamente

desarrollado y se caracteriza por un flujo aireado en la

superficie libre, formación de burbujas y vórtices estables en

los escalones (ver figura 1). (Khatsuria, R.M., 2005)

Figura 1. Flujo Rasante sobre una rápida escalonada.

En la figura 2, se observan las principales regiones existentes

en un flujo rasante y estas son: (1) Flujo no aireado, se presenta

el crecimiento de la capa límite hasta la superficie libre del

flujo; (2) Punto de inicio de la autoaireación, se desarrolla un

flujo parcialmente aireado; (3) Flujo gradualmente variado y

(4) Flujo en estado de equilibrio, uniforme y completamente

aireado.

Las consideraciones más relevantes para el diseño de una

rápida escalonada con flujo rasante son: Propiedades del flujo

aireado, la estimación de la resistencia al movimiento del agua,

determinación del coeficiente de fricción, introducción y

arrastre de aire, disipación de energía, consideración de las

fluctuaciones de presión y cavitación. (Khatsuria, R.M., 2005)

Figura 2. Principales regiones existentes en un flujo rasante.

2. MODELO FÍSICO

El modelo físico-hidráulico del flujo bifásico agua-aire,

altamente turbulento cumple la similitud restringida de

Froude. En el caso de flujos a superficie libre se considera

como fuerza predominante la Gravedad. Estrictamente se debe

cumplir que el número de Froude en prototipo sea igual al del

modelo. A la vez el modelo físico debe satisfacer los

parámetros adimensionales del número de Reynolds (Re >105)

y el número de Weber (√We > 110) para evitar los efectos de

escala producidos por la viscosidad y la tensión superficial

para un flujo bifásico (Pfister M., Chanson H., 2013).

El número de Reynolds relaciona las fuerzas de inercia con las

fuerzas viscosas y el número de Weber relaciona las fuerzas de

inercia con las fuerzas de tensión superficial. En la tabla 1, se

presentan los valores de estos parámetros adimensionales que

cumplen con los rangos establecidos para evitar los efectos de

escala en el modelo físico y así obtener resultados confiables.

Tabla 1. Parámetros adimensionales. Similitud restringida de Froude

CAUDAL

MODELO

NÚMERO DE

FROUDE

NÚMERO DE

REYNOLDS

NÚMERO

DE WEBER

Qm Frp=Frm Rem (Wem)^0.5

l/s - - -

13.97 7.3 3.9E+05 119

22.36 6.8 6.2E+05 161

27.95 6.9 7.5E+05 180

41.92 6.2 1.3E+06 266

El modelo numérico se elabora con las dimensiones del

modelo físico a escala 1:20 construido en el Laboratorio del

Centro de Investigaciones y Estudios en Recursos Hídricos

(CIERHI) de la Escuela Politécnica Nacional del Ecuador. En

la tabla 2, se describen las dimensiones más relevantes de la

rápida escalonada en prototipo y en modelo físico. (CIERHI,

EPN TECH, 2016)

Tabla 2. Dimensiones en prototipo y en modelo físico

RÁPIDA ESCALONADA

Estructur

a Descripción

Unida

d

Prototip

o Modelo

Rápida

escalonad

a (Una

cámara)

Longitud L m 98.85 4.94

Altura H

d m 92.26 4.61

Ancho b m 5.00 0.25

Huella l m 1.00 0.05

Contrahuella h m 1.00 0.05

# escalones e u 87 87

Inclinación ° 45 45

Caudal

diseño Q mᶟ/s 75.00 0.0419

En la figura 3 se observa un esquema del perfil longitudinal de

la rápida escalonada con el cuenco de disipación de energía al

pie de la estructura y en la figura 4 la representación digital del

modelo físico en tres dimensiones considerando las zonas

representativas de ingreso, estructura de disipación y salida.

54


Flow-3D


Figura 3. Dimensiones de la rápida escalonada el Batán.

Figura 4. Modelo físico en 3D de la rápida escalonada el Batán.

3. MODELO NUMÉRICO

El software FLOW-3D resuelve en tres dimensiones la

ecuación de la continuidad, la ecuación de la cantidad de

movimiento y para la definición de la superficie libre del agua

el modelo VOF, estas ecuaciones son las siguientes

(Fernández Oro J. M., 2012).

Vf

ρ

∂ρ

∂t+

1

ρ∇(ρu⃗ Af) = −

∂Vf

∂t (1)

∂u⃗⃗

∂t+

1

Vf(u⃗ Af. ∇u⃗ ) = −

1

ρ[∇P + ∇(τ. Af)] + G⃗⃗ (2)

∂F

∂t+

1

Vt∇. (F. u⃗ . Af) = −

F

Vf

∂Vf

∂t (3)

Donde, ρ es la densidad del fluido, u⃗ es el vector de velocidad

del flujo, Vf es la fracción de volumen, Af es la fracción de

área, P es la presión, τ es el tensor de esfuerzos viscosos, G es

la gravedad y F es la fracción del fluido.

Las variables fluido-dinámicas en un punto en el espacio están

conformadas por una serie de fluctuaciones de distintas

escalas, por tal razón, el análisis de la turbulencia se realiza

desde el punto de vista estadístico o sea con valor medio de la

velocidad e intensidad promediada de las fluctuaciones

“promedios de Reynolds.” Esta consideración ayuda en la

resolución del problema de cierre del sistema de ecuaciones

algebraicas. Para la solución numérica de los flujos turbulentos

se aplican las ecuaciones de Navier-Stokes promediadas por

Reynolds (RANS), en el cuál se calculan la energía cinética

turbulenta (k) y la tasa de disipación de energía (ε) para

obtener el tensor de esfuerzos viscosos o de Reynolds (τij =

−uiuj̅̅ ̅̅ ̅) y la viscosidad turbulenta (υt).

−uiuj̅̅ ̅̅ ̅ = υt (∂ui

∂xj+

∂uj

∂xi) −

2

3δij. k (4)

υt = cuk2

ε (5)

El tiempo establecido para que la simulación numérica del

flujo sobre la rápida escalonada alcance la convergencia, se

estima considerando el tiempo que se demora el flujo en entrar

y salir del sistema hasta que la solución numérica llegue a un

equilibrio o se estabilice tanto en la masa como en la energía

cinética. En la modelación numérica la variación del

parámetro de la energía cinética turbulenta presenta

frecuencias periódicas debido a que el flujo rasante resulta ser

intermitente y desciende en forma de series por la rápida

escalonada. Para el presente estudio se observa que desde el

tiempo simulado de 20 segundos el modelo numérico alcanza

la estabilidad y convergencia con un paso de tiempo promedio,

∆t=2.5 x 10-4s. La duración en tiempo real de las simulaciones

numéricas hasta los 20 segundos establecidos es de alrededor

de 5 horas. Las características del ordenador utilizado son las

siguientes: Memoria RAM de 16 GB y Procesador Intel ®

Core ™ i7 de 3.69 GHz (8 procesadores).

3.1 Discretización espacial y condición de frontera.

En este análisis el flujo rasante en una rápida escalonada se

desarrolla en 2 dimensiones “2D”, debido a la capacidad del

ordenador disponible y bajo la hipótesis de que las

características hidrodinámicas del flujo no varían

considerablemente en el ancho (B), es decir, se cumple la

condición B > 4y (y es la profundidad del flujo) (Bombardelli

et al., 2010). Sin embargo, esta consideración simplifica el

fenómeno físico ya que puede existir presencia de ondas en el

sentido transversal del flujo rasante (Valero, D., Bung, D.,

2015).

El solucionador emplea el método de diferencias finitas para

la discretización espacial de las ecuaciones del flujo. El

mallado se ingresa en un solo bloque de malla y divide el

espacio físico con celdas rectangulares. Las celdas utilizan sus

nodos y caras para almacenar valores de las incógnitas como

la presión y la velocidad. En este estudio se analiza la

influencia del tamaño de malla para celdas uniformes de 2, 3,

4 y 5 milímetros.

Las condiciones de frontera son datos conocidos para la

solución del campo hidrodinámico, se colocan en los límites

del bloque de malla que discretiza la región de interés como se

muestra en la figura 5. Se debe tener en cuenta que una correcta

elección de las condiciones de frontera es indispensable para

la estabilidad y convergencia numérica de la solución.

55



Figura 5. Condiciones de frontera para la modelación.

En la frontera Xmin se coloca la condición de caudal de entrada

(Q). Se ingresa el valor del caudal y del nivel de agua.

En la frontera Ymin, Ymáx se colocan la condición de simetría

(S). Es decir en la superficie del sólido el fluido tiene libertad

de deslizamiento.

En la frontera Xmáx se coloca la condición de salida del flujo

(O).

En la frontera Zmáx se coloca la condición de abierto a la

presión atmosférica (P).

La rugosidad absoluta del material en el modelo físico

(acrílico), tanto para las paredes y el fondo es de 0.0015mm.

3.2 Sub-modelo: Incorporación de aire.

El sub-modelo de ingreso de aire identifica cuando el aire se

incorpora al flujo. El flujo turbulento genera fácilmente la

introducción y transporte de aire. El sub-modelo de

incorporación y arrastre de aire en FLOW 3D evalúa un

equilibrio entre la perturbación Pt (Energía cinética turbulenta

por unidad de volumen) y las fuerzas estabilizadoras Pd

(Fuerza de gravedad y tensión superficial) (Khatsuria, R.M.,

2005). Cuando Pt > Pd, el aire rompe la superficie libre del

flujo y es transportado en la masa de agua, formando de esta

manera un flujo autoaireado. (Ver figura 6).

Figura 6. Mecanismo de la introducción de aire. (A) modelo de la introducción de aire. (B) Fuerzas actuantes en el líquido elemental.

3.3 Sub-modelo: Esponjamiento del flujo, modelo drift-flux.

El Sub-modelo de esponjamiento de flujo identifica un flujo

compuesto de agua con burbujas de aire en una rápida

escalonada, produciendo variaciones de la velocidad debido a

las diferentes densidades. Se considera que el aire es

transportado en forma de burbujas con un diámetro

característico así las burbujas producen una fuerza de arrastre

que generan una oposición al movimiento del agua. El sub-

modelo drift-flux calcula una velocidad relativa para la mezcla

agua-aire, fase continua y la fase dispersa respectivamente. La

velocidad relativa se obtiene de la siguiente ecuación: (Flow

Science, Inc., 2012)

δur

δt+ u2.∇u2 − u1.∇u1 = (

1

ρ1−

1

ρ2) ∇P − (

1

fρ1−

1

(1−f)ρ2) kd. ur (6)

Donde, f es la fracción de volumen de la fase continua, kd es

un coeficiente de fricción que da cuenta de la interacción de

las dos fases, ρ1, ρ2 es la densidad de la fase continua y

dispersa respectivamente y ur es la velocidad relativa.

Debido a la complejidad del fenómeno los resultados de las

investigaciones numéricas reportan resultados no muy

precisos en el transporte de aire y en el aumento de volumen

del flujo mixto en la zona uniforme. La limitación radica en el

transporte de aire ya que el modelo numérico representa en el

contorno de la superficie libre una concentración de aire

C=90%, es decir, aireado al 90%. En algunos casos, el flujo

sobre una rápida escalonada tiene una concentración de aire

menor al 90% y el programa asume como frontera la superficie

libre con concentración de aire al 90% y es así como ocurren

errores en la superficie libre generando un mayor aumento de

volumen del flujo autoaireado (Sarfaraz, M. and Attari, J.,

2011).

3.4 Parámetros iniciales, físicos y numéricos de ingreso en el

programa FLOW-3D.

56


Flow-3D


A continuación, se presentan los parámetros iniciales, físicos

y numéricos que se activan para el desarrollo de la simulación

numérica del flujo rasante sobre una rápida escalonada

considerando un flujo bifásico “agua-aire”.

PARÁMETROS INICIALES

Tiempo de simulación 20 s.

Fluido incompresible.

Flujo a superficie libre.

Fluido: Agua a 20°C

Unidades: Sistema Internacional.

Temperatura: Celsius.

PARÁMETROS NUMÉRICOS

Método VOF.

Resolver la ecuación de momento

Segundo grado de aproximación.

Algoritmo de acoplamiento

presión-velocidad (GMRES)

PARÁMETROS FÍSICOS

Gravedad.

Viscosidad y turbulencia.

Modelo k- RNG.

TLEN= 7% profundidad del flujo.

Incorporación y transporte de aire.

Coeficiente de arrastre de aire= 0.5

Coeficiente de tensión superficial= 0

Densidad del aire= 1.20 kg/m3

Densidad del flujo.

Segundo grado de aproximación.

Esponjamiento del flujo.

Radio promedio de la burbuja de aire= 0.025mm a 1mm

Viscosidad dinámica del flujo mixto= 0.001 kg/(m.s)

Viscosidad dinámica del agua= 0.001 kg/(m.s)

Viscosidad dinámica del aire= 1.8x10-5 kg/(m.s)

Densidad del agua= 1000 kg/m3

3.5 Calibración de la Modelación Numérica.

Se analiza la sensibilidad de los resultados para cuatro

diferentes tamaños de celda y se procede con la selección del

modelo de turbulencia más apropiado para representar el

fenómeno físico con el fin de alcanzar la calibración de la

modelación numérica.

3.5.1 Influencia del modelo de turbulencia.

Los modelos RANS constan de dos ecuaciones que son la

energía cinética turbulenta (K) y la disipación de energía

turbulenta (). En el estudio se analiza el modelo k- (RNG) y

el modelo K-. El parámetro de la viscosidad turbulenta

interviene en el proceso de la disipación de energía a lo largo

de la rápida escalonada con flujo rasante (Lucio I. et al., 2015),

por lo tanto, se debe verificar los valores de la longitud de

mezcla turbulenta (TLEN), para no tener valores erróneos de

la disipación de energía (Flow Science, Inc., 2012). El valor

de TLEN depende del tamaño de celda por lo que existe la

necesidad de calibrar este valor si se observan valores erróneos

de la viscosidad turbulenta (T).

El programa FLOW-3D calcula el valor de TLEN

automáticamente. Este valor lo determina como el 7% del

tamaño de celda más pequeño, lo cual no es aconsejable. Lo

más habitual y recomendable es adoptar el valor de TLEN

como el 7% a 15% del tirante de agua en un canal (y), según

se ajuste con la calibración del modelo numérico (ARAGUA,

2013). En la tabla 3 se presentan las simulaciones numéricas

planteadas para el análisis de la influencia del modelo de

turbulencia en un flujo rasante para un tamaño de malla de 4

mm.

Tabla 3. Plan de simulaciones. Influencia del modelo de turbulencia

Nombre

Simulación

Caudal

modelo

Modelo de

Turbulencia TLEN

- l/s - m

S1

27.95

k- (RNG) Automático

S4 0.0028

S5 K-

Automático

S6 0.0028

En la figura 7 se presenta las profundidades del flujo bifásico

a lo largo de la rápida escalonada, los resultados muestran que

aguas abajo del punto de inicio de la autoaireación se tiene

inestabilidad de la superficie libre. El modelo de turbulencia

K- sobreestima la profundidad del flujo bifásico y no

representa adecuadamente el flujo real. Con este modelo se

reporta una disminución de la velocidad del flujo debido a que

los valores de la viscosidad turbulenta son diferentes a los

valores reales lo cual genera mayores pérdidas de energía.

Se evidencia en un modelo RANS, el efecto que tiene el valor

de TLEN, sin realizar la debida calibración correspondiente

pues los resultados de la viscosidad turbulenta, físicamente

representan una mayor resistencia de las partículas del fluido

al movimiento. Por tal motivo se generan velocidades de flujo

menores y sobrelevaciones del nivel de flujo bifásico como se

puede observar en la figura 7.

El modelo k- Renormalizado con un valor de TLEN

correspondiente al 7% de la profundidad del flujo rasante e

igual a 0.0028, presenta resultados semejantes con los

obtenidos en el modelo físico con una desviación menor al 5%

debido a la fluctuación del flujo. La ventaja del modelo k-

RNG es que representa adecuadamente flujos altamente

turbulentos y es una característica del flujo rasante, mientras

que el modelo K- permite representar o modelar de mejor

manera el comportamiento del flujo en las paredes o

contornos.

57



Figura 7. Influencia del modelo de turbulencia. Qmodelo=27.95l/s.

3.5.2 Influencia del tamaño de malla.

Las simulaciones se desarrollan con el modelo de turbulencia

k- Renormalizado con un valor de TLEN=0.0028. Para este

análisis los tamaños de celda se seleccionan considerando los

siguientes criterios: Un tamaño de celda muy grande causa una

discontinuidad del flujo, lo cual afecta la estabilidad numérica

de la solución, un tamaño de celda muy pequeño genera un

costo computacional alto y además considerando que la

profundidad del flujo es aproximadamente 3 cm para el caudal

operación más bajo. Se consideran cuatro diferentes tamaños

de celda uniforme que son: 5, 4, 3 y 2mm correspondientes a

las simulaciones: S0, S1, S2 y S3 respectivamente. Los

resultados obtenidos acerca de las profundidades del flujo a lo

largo de la rápida escalonada indican claramente que el tamaño

de malla no influye en la profundidad del flujo en la zona no

aireada, es decir en los primeros 14 escalones. En la zona del

flujo completamente autoaireado se observa niveles de agua

inferiores para el tamaño de celda más grande que para el

tamaño de celda más pequeño. Si bien con un tamaño de celda

más pequeño se tiene una mayor precisión en la resolución del

flujo, se observa en la figura 8 que existe una sobreelevación

de las profundidades del flujo en la zona completamente

autoaireada para el tamaño de celda más pequeño, esto se debe

a la influencia del modelo de esponjamiento del flujo (drift

flux). En la realidad, durante el transporte de aire, el tamaño

de las burbujas de aire es variable en el tiempo y en espacio,

sin embargo, este algoritmo resuelve con un tamaño constante

de las burbujas de aire. Para una celda de mayor tamaño la

profundidad del flujo autoaireado es menor y estos resultados

son acordes con los registrados en el modelo físico. Para el

caso particular de la rápida escalonada el Batán, según el

análisis de la influencia del tamaño de celda se obtuvo buenos

resultados con un mallado uniforme de 5mm.

Figura 8. Influencia del tamaño de malla. Qmodelo=27.95l/s,

QPrototipo=50m3/s.

Del análisis de calibración del modelo numérico se

recomienda que para futuros estudios de investigación se

debería analizar con más detalle la relación entre el modelo de

turbulencia con el valor de TLEN considerando el transporte

de las burbujas de aire con tamaños variables durante el

descenso del flujo sobre la rápida escalonada.

4. ANÁLISIS DE LOS RESULTADOS NUMÉRICOS

Definidos los parámetros numéricos para la simulación del

flujo rasante en la rápida escalonada usando el paquete

comercial FLOW-3D y calibrados conforme el modelo

experimental, se plantea un rango de caudales de operación

en una cámara. En la tabla 4 se muestra el plan establecido

para las simulaciones numéricas.

Tabla 4. Plan de simulaciones numéricas

Caudal

prototipo

Caudal

modelo Observaciones

Tipo de

flujo

Qp (m3/s) Qm (l/s)

25 13.97 33.4% del Qdiseño

Rasante 40 22.36 53.3% del Qdiseño

50 27.95 66.7% del Qdiseño

75 41.92 Qdiseño

4.1 Comparación entre los resultados obtenidos con la

modelación numérica y la modelación física.

En el presente ítem se compara los resultados de las

mediciones experimentales y registros fotográficos obtenidos

del modelo físico con los resultados de las simulaciones

numéricas usando el programa FLOW-3D.

4.1.1 Comparación del flujo de aproximación hacia la rápida

escalonada.

58


Flow-3D


Aguas arriba de la rápida escalonada se tiene un cuenco

disipador de energía, esta estructura disipa energía y entrega

un flujo subcrítico y uniformizado hacia la rápida escalonada.

El parámetro hidráulico del nivel de agua en el cuenco se

utiliza como condición de frontera para la simulación

numérica. El flujo de aproximación hacia la rápida escalonada

para el caudal de 27.95 l/s, presenta una profundidad de agua

en la corona de la rápida de 11cm de acuerdo con el registro

experimental. En la figura 9 se presenta el resultado de la

simulación numérica y el valor de la profundad de agua en la

sección de control (corona de la rápida) es de 10.7cm. La

comparación de los resultados de la profundidad de agua en la

sección de control (profundidad crítica) muestran una

desviación del 2.8% (Mohammad S. et al., 2012).

Figura 9. Flujo de aproximación hacia la rápida escalonada. Modelación

numérica en FLOW-3D. Qmodelo=27.95 l/s.

4.1.2 Ubicación del punto de inicio de la autoaireación.

La localización del punto de inicio de la autoaireación (Li) se

determina siguiendo el desarrollo de la capa límite. Cuando la

capa límite llega hasta la superficie libre del flujo se encuentra

el punto de inicio de la autoaireación (Flow Science, Inc.,

2012). Inmediatamente hacia aguas abajo de este punto la

turbulencia generada introduce y transporta el aire durante el

descenso del flujo por la rápida escalonada.

Figura 10. Simulación numérica en FLOW-3D. Qmodelo=27.95l/s,

QPrototipo=50m3/s.

En la figura 10, se reportan los resultados para el caudal en

modelo de 27.95 l/s equivalente a 50 m3/s en prototipo. La

aparición constante de vórtices en los escalones se produce a

partir del escalón número 11. Desde esta posición hacia aguas

abajo se presenta un flujo bifásico, con concentraciones de aire

que van en aumento conforme el flujo desciende por la rápida

escalonada. En la figura 11 se muestra los resultados teóricos,

experimentales y numéricos de la localización del punto de

inicio de la autoaireación para diferentes caudales de

operación.

Figura 11. Comparación entre los resultados experimentales y los

resultados de la simulación numérica.

4.1.3 Profundidad del flujo mixto con concentración del 90%

de aire.

La región uniforme se alcanza a una distancia alejada de la

cresta de la rápida, donde se establece el equilibrio entre las

fuerzas de gravedad y fricción. En esta zona se mantienen

constantes las características del flujo como concentración

media de aire, velocidad y profundidad del flujo.

En la simulación numérica del flujo rasante los resultados de

la concentración aire se basan en el análisis de los valores de

la fracción volumétrica de aire en las celdas. En las figuras

12A y 12B se presenta el flujo observado en el modelo físico

con el caudal simulado de 27.95 l/s. En la región de equilibrio

o flujo uniforme, se puede visualizar aproximadamente una

concentración de aire alrededor del 90%. La profundidad del

flujo rasante perpendicular al fondo virtual es

aproximadamente igual a 5.5 cm.

Figura 12. (A) Modelo físico, vista lateral (B) Modelo físico, vista frontal

(C) Modelación numérico en FLOW-3D, vista lateral. Qmodelo=27.95l/s,

QPrototipo=50m3/s.

59



Los resultados de la modelación numérica según la figura 12C

presentan una profundidad media del flujo aireado al 90%

igual a 5.3 cm en la región uniforme. Comparando el valor

promedio experimental con el valor obtenido de la modelación

numérica se tiene una diferencia porcentual del 3.8%. Aguas

abajo de la sección donde se forma el flujo rasante, se observa

que existe una alta aireación del flujo. La coloración del flujo

es blanquecina como resultado de la intensa autoaireación,

característica del flujo rasante que llega hasta el pie de la

rápida como se puede observar en la figura 12B.

Figura 13. Comparación entre los resultados experimentales y los

resultados de la simulación numérica. dc/h vs y90.

En la figura 13 se reporta en el eje de las abscisas el parámetro

adimensional “dc/h” y en el eje de las ordenadas la

profundidad del flujo con concentración de aire al 90% “y90”,

se concluye que los resultados experimentales y numéricos

tienen una similar tendencia. Se puede visualizar que a medida

aumenta dc/h (medida del caudal de la descarga en la rápida

escalonada) la profundidad del flujo mixto también aumenta.

Los resultados obtenidos mediante del estudio experimental y

modelación numérica sobre la profundidad de flujo aireado en

la región de equilibrio para la serie de caudales muestran una

semejanza entre sí. Estos resultados presentan una desviación

inferior al 10 %.

4.2 Características hidrodinámicas del flujo rasante.

Se analizan los resultados numéricos de las variables

hidrodinámicas más relevantes del flujo rasante como: i)

Distribución de presiones en el escalón; ii) Distribución de

velocidades en el escalón; iii) Variación de la velocidad del

flujo en toda la longitud de la rápida escalonada. Ubicación de

la región del flujo uniforme; iv) Profundidades de flujo a lo

largo de la rápida escalonada; v) Concentración de aire,

cálculo del factor de fricción de Darcy- Weisbach y vi)

Disipación de energía al pie de la rápida.

4.2.1 Distribución de presiones en el escalón.

Los resultados acerca de la distribución de presiones en el

escalón como se observa en la figura 14, indican un mayor

valor de presión en el tercio final de la superficie horizontal a

la salida del escalón. Mientras en la esquina del escalón, las

presiones relativas son bajas, debido a la separación del flujo.

En la esquina superior del escalón para el caudal de 50 m3/s,

en la región del flujo uniforme (Ver figura 14), el valor de la

presión relativa o manométrica es de aproximadamente -685

Pa.

En la zona de las cavidades de los escalones donde se genera

la recirculación del flujo y en la capa superior del flujo rasante

las presiones son bajas debido a la alta aireación que adquiere

el flujo debido al descenso hasta el cuenco disipador de

energía.

Figura 14. Distribución de presiones en el escalón. Qmodelo=27.95l/s,

QPrototipo=50m3/s.

4.2.2 Distribución de velocidades en el escalón.

La magnitud del vector velocidad crece desde el fondo virtual

hacia la superficie libre como se muestra en la figura 15. La

distribución de la velocidad en la zona del flujo rasante

presenta una tendencia exponencial según la ecuación 7.

u

u90= (

y

y90)

1

n (7)

Dónde: u90 es la velocidad máxima del flujo libre, y90 es la

profundidad de flujo aireado al 90%, medida perpendicular al

pseudo-fondo, (u, y) son las coordenadas de la velocidad y

profundidad del flujo y n es el valor del exponente que según

investigaciones experimentales varía entre 3.5 y 6 (Flow

Science, Inc., 2012).

60


Flow-3D


Figura 15. Esquema de la distribución de velocidades.

Para la serie de caudales presentados en el plan de las

simulaciones numéricas la velocidad del flujo en la zona

uniforme oscila entre 4.5 m/s y 6.4 m/s. En la figura 16 se

compara la distribución de velocidades del flujo obtenido para

el caudal en modelo de 41.93 l/s, con la curva que resulta de la

ecuación 7. De acuerdo con los resultados se obtuvo un mejor

ajuste de la curva exponencial teórica de la distribución de

velocidades con el valor de n=6.

Figura 16. Distribución de velocidades en el escalón N° 70.

Qmodelo=41.93l/s, QPrototipo=75m3/s.

4.2.3 Variación de la velocidad a lo largo de la rápida

escalonada. Localización de la región uniforme.

Por el efecto de la gravedad, la velocidad del flujo aumenta

conforme desciende por la rápida escalonada hasta una cierta

región donde la velocidad tiende a un valor constante. Esta

zona se le conoce como la región de flujo uniforme. El aire

transportado reduce la fricción del flujo en las paredes a lo

largo de la rápida. A partir de la figura 17, se ubica

aproximadamente la región de flujo uniforme. La distancia a

lo largo de la rápida escalonada desde la cresta de la misma

hasta la localización de la sección donde se presentan

velocidades de flujo cuasi-uniformes lo nombraremos

distancia (Lr).

En la tabla 5 se muestra los resultados del modelo numérico

sobre la distancia hasta la formación de la zona uniforme, la

velocidad media del flujo y el parámetro adimensional (dc/h)

relaciona la profundidad crítica (dc) del flujo en la rápida

respecto a la altura del escalón (h).

Figura 17. Variación de la velocidad del flujo en la rápida escalonada.

Qmodelo=27.95l/s, QPrototipo=50m3/s.

Tabla 5. Características hidráulicas del flujo rasante

Caudal

modelo

Parámetro

adimensional

Distancia hasta la

formación de la zona

uniforme

Velocidad

media

Qm dc/h Lr Vm

l/s - m m/s

13.97 1.37 3.70 4.54

22.36 1.87 4.34 5.3

27.95 2.17 4.62 5.5

41.92 2.84 5.12 6.4

En la figura 18 se muestra la relación entre el parámetro

adimensional (dc/h) vs la velocidad media del flujo y la

distancia hasta la formación de la región uniforme. La

tendencia indica que a medida que se incrementa la relación

dc/h la velocidad de flujo también aumenta y la localización

de la región uniforme se desplaza hacia agua abajo de la

rápida, es decir, se incrementa.

61



Figura 18. Resultados de la modelación numérica en FLOW-3D.

4.2.4 Profundidades del flujo a lo largo de la rápida

escalonada.

Aguas abajo del punto de inicio de la autoaireación se observa

inestabilidad de la superficie libre de agua debido a la gran

turbulencia generada y al alto porcentaje de concentración de

aire en el flujo, teniendo así flujo de mezcla agua-aire. La

figura 19 reporta para la serie de caudales simulados los

resultados numéricos obtenidos sobre la profundidad de agua

en toda la longitud de la rápida y en particular en la región del

flujo uniforme (flujo aireado al 90%).

Figura 19. Profundidades de flujo a lo largo de la rápida escalonada.

4.2.5 Determinación del factor de fricción de Darcy-

Weisbach, en la región uniforme del flujo rasante.

Los resultados obtenidos del modelo numérico indican una

concentración de aire del 90 % en la región uniforme del flujo

rasante. Con fines de calcular el factor de fricción (fe) entre el

fondo virtual y el flujo que pasa por encima del pseudo fondo,

se determina previamente la profundidad de agua clara en la

región uniforme y la concentración media de aire. En la figura

20 se reportan los resultados del factor de fricción (fe)

determinados según la ecuación 8, con base en los resultados

numéricos obtenidos (Chanson H., 1993).

fe =8∗g∗sen θ∗de

3

qw2

(8)

Donde, fe es el factor de fricción de Darcy-Weisbach, de es la

profundidad de agua clara en la región uniforme, y90 es la

profundidad de agua con el 90% de concentración de aire, qw

es el caudal unitario y es el ángulo de inclinación de la rápida

escalonada.

Figura 20. Factor de fricción para la rápida escalonada (ϴ=45º) para el

rango de caudales analizados.

Según la literatura técnica el valor de fe varía en un rango de

0.11, 0.17 y 0.30 (resultados determinados experimentalmente

bajo las siguientes condiciones: θ > 10°, h > 0.02 m y Re > 1

× 105) (Chanson H. et al., 2015). Los resultados obtenidos en

este estudio indican que el factor de fracción varía entre 0.11

a 0.148 para el rango de caudales analizados.

4.2.6 Disipación de energía al pie de la rápida escalonada.

En este análisis se ubican dos secciones de control. La sección

1 se ubica aguas arriba de la rápida y la sección 2 se ubica en

el cuenco disipador de energía ubicado al pie de dicha

estructura (Ver figura 3). Los porcentajes de disipación de

energía al pie de la rápida escalonada para la serie de caudales

simulados se reportan en la tabla 6.

Tabla 6. Disipación de energía al pie de la rápida escalonada.

- Caudal

modelo

Energía

z1+y1+V2/2g ∆H

N Qm E1 E2 (E1-E2)/E1

- l/s m m %

S8 13.97 4.768 1.094 77.0

S9 22.36 4.809 1.499 68.8

S10 27.95 4.835 1.612 66.7

S11 41.92 4.888 2.164 55.7

Los resultados muestran que la disipación de la energía al pie

de la rápida es baja para el rango de caudales altos, debido a la

62


Flow-3D


gran cantidad de energía cinética que se tiene al pie de la

rápida. El rango de los porcentajes de disipación de energía

varía entre 55.7 % a 77.0%.

5. CONCLUSIONES

Para este estudio del flujo rasante en una rápida escalonada con

inclinación a 45º y dimensiones de los escalones de 5 x 5cm se

obtuvo que la combinación del modelo de turbulencia K-

RNG con el algoritmo TruVOF más el sub-modelo de

aireación y el sub-modelo de esponjamiento de aire, representa

el flujo rasante incluyendo la introducción espontánea de aire

así como el movimiento del flujo mixto sobre un fondo virtual.

Sin embargo, se debe tener en cuenta las simplificaciones de

los algoritmos que posee el programa y que condicionan el real

comportamiento del fenómeno de autoaireación.

El incremento del caudal de descarga en la rápida produce un

aumento de la velocidad del flujo rasante, este

comportamiento genera un desplazamiento hacia aguas abajo

la ubicación del punto de inicio de la autoaireación y la

ubicación de la región del flujo uniforme.

Para los cuatro caudales de operación la velocidad del flujo en

la región uniforme varía entre 4.5 m/s a 6.4 m/s y el factor de

fricción varía entre 0.105 a 0.148 para un rango de caudales

que va desde 13.97 l/s hasta los 41.92 l/s.

Además, con base en los resultados obtenidos de la

modelación numérica sobre la profundidad del flujo bifásico y

la velocidad al pie de la rápida escalonada se ha determinado

que la energía residual al pie de la misma es del 44% para el

caudal de diseño, Qprototipo= 75m3/s.

En la región de uniforme del flujo rasante las características

hidrodinámicas fueron calculadas con una buena precisión en

el modelo numérico utilizando el FLOW 3D, ya que los

resultados obtenidos se asemejan a los registrados en el

modelo físico experimental. Las comparaciones de estos

resultados presentan una desviación inferior 10% debido a que

el flujo es altamente turbulento.

AGRADECIMIENTOS

Los autores agradecen al Centro de Investigaciones y Estudios

en Recursos Hídricos CIERHI-EPN por las facilidades y

financiamiento brindado para le ejecución de este proyecto de

investigación.

REFERENCIAS

ARAGUA. (2013). “Modelación numérica y experimental de flujos aire-agua

en caídas en colectores.”, Laboratório Nacional de Engenharia Civil, I. P. Av do Brasil 101 • 1700-066 Lisboa.

Bombardelli, F.A., Meireles, I. and Matos, J., (2010), “Laboratory

measurement and multi-block numerical simulations of the mean flow and turbulence in the non-aerated skimming flow region of steep stepped

spillways”, Environ Fluid Mechanics.

Castro M. (2015) “Análisis Dimensional y Modelación física en Hidráulica”.

Escuela Politécnica Nacional. Quito Ecuador. 50 p.

Chanson H., D. B. Bung., J. Matos (2015). “Stepped spillways and cascades”.

IAHR Monograph. School of Civil Engineering, University of Queensland, Brisbane, Australia.

Chanson H. (1993). "Stepped Spillway Flows and Air Entrainment." Can. Jl

of Civil Eng., Vol. 20, No. 3, June, pp. 422-435 (ISSN 0315-1468).

CIERHI, EPN TECH, (2016). “Estudio experimental en modelo físico de las

rápidas con perfil escalonado y liso de la quebrada el Batán Fase I y Fase II”, Escuela Politécnica Nacional, Quito Ecuador.

Fernández Oro J. M. (2012)., “Técnicas Numéricas en Ingeniería de Fluidos:

Introducción a la Dinámica de Fluidos Computacional (CFD) por el Método de Volúmenes Finitos”. Barcelona: Reverté.

Flow Science, Inc. (2012). “FLOW 3D 10.1.0 Documentation Release.

Manual de Usuario”, Los Alamos National Laboratory. Santa Fe, New México

Khatsuria, R.M., (2005)., “Hydraulics of Spillways and Energy Dissipators”.

Department of Civil and Environmental Engineering Georgia Institute of Technology Atlanta.

Lucio I., Matos J., Meireles I. (2015). “Stepped spillway flow over small

embankment dams: some computational experiments”. 15th FLOW-3D European users conference.

Mohammad S., Jalal A. and Michael P., (2012). “Numerical Computation of

Inception Point Location for Steeply Sloping Stepped Spillways” 9th

International Congress on Civil Engineering. Isfahan University of

Technology (IUT), Isfahan, Iran

Pfister M., Chanson H., (2013), “Scale Effects in Modelling Two-phase Air-water Flows”, Proceedings of 2013 IAHR World Congress.

Sarfaraz, M. and Attari, J. (2011), “Numerical Simulation of Uniform Flow

Region over a Steeply Sloping Stepped Spillway”, 6th National Congress on Civil Engineering, Semnan University, Semnan, Iran.

Valero, D., Bung, D., (2015), “Hybrid investigation of air transport processes in moderately sloped stepped spillway flows”, E-proceedings of the 36th

IAHR World Congress 28 June – 3 July, 2015, The Hague, the

Netherlands.

BIOGRAFÍAS

Edwin Patricio Casa Tipán (1988). Ingeniero Civil – Hidráulico de la

Universidad Central del Ecuador.

Magister en Recursos Hídricos mención

“Diseño de Proyectos Hidráulicos”,

Escuela Politécnica Nacional, EPN,

Quito. Profesor Ocasional a tiempo

completo en la Facultad de Ingeniería

Civil y Ambiental de la EPN. Realiza estudios de investigación

en el Laboratorio del Centro de Investigaciones y Estudios en

Recursos Hídricos, CIERHI de la EPN, sobre temas de diseño,

evaluación y optimización de estructuras hidráulicas con base

en la modelación física-hidráulica y modelación numérica.

Ximena Hidalgo Bustamante (1960). Ingeniera Civil especialidad Hidráulica de

la Escuela Politécnica Nacional, Master of

Science en Diseño de Estructuras

Hidráulicas y Modelos Hidráulicos en el

Instituto de Hidromecánica de la

Universidad de Karlsruhe, Alemania.

Actual Decana y Profesora titular a tiempo

completo de la Facultad de Ingeniería Civil y Ambiental de la

EPN desde 1988. Directora del Centro de Investigaciones y

Estudios en Recursos Hídricos CIERHI-EPN. Investigadora

reconocida en el área de modelación física de estructuras

hidráulicas, diseño y optimización de proyectos desarrollados

en Ecuador.

63



Marco Antonio Castro Delgado (1948). Doctor Ingeniero en Ingeniería Civil en

Karlsruhe Institute of Technology - KIT -

Karlsruhe, Alemania. Ingeniero Civil

especialización Hidráulica de la EPN.

Profesor titular a tiempo completo desde

el año 1981 en la Facultad de Ingeniería

Civil y Ambiental de la EPN.

Investigador, Director técnico y Jefe de diseño en Estudios y

Diseños en más de 100 proyectos hidráulicos, hidrológicos,

hidroeléctricos, hidro-sanitarios y de infraestructura a nivel de

factibilidad y diseño definitivo en todo el país.

Patricio Rubén Ortega Lara (1986). Ingeniero Civil especialidad Hidráulica

graduado en la Universidad Central de

Ecuador. Magister en Recursos Hídricos

mención Diseño de Proyectos Hidráulicos

en la EPN. Profesor Titular Tiempo

Completo desde 2015 en la Facultad de

Ingeniería Civil y Ambiental de la EPN.

Coordinador del Centro de Investigaciones y Estudios en

Recursos Hídricos CIERHI-EPN. Presidente de la Red de

jóvenes Investigadores de la IAHR. Editor Asociado de la

Revista Hidrolatinoamericana de Jóvenes Investigadores y

profesionales. Investigador acreditado que cuenta con 3

patentes en estructuras separadoras de caudal.

Pablo Alberto Vera Romero (1981). Ingeniero Civil y Magister en Recursos

Hídricos mención en Diseño de proyectos

hidráulicos de la Escuela Politécnica

Nacional. Tiene experiencia docente en el

Departamento de formación básica y el

laboratorio de hidráulica de la Facultad de

Ingeniería Civil y Ambiental de la EPN.

Actualmente es integrante del equipo técnico del Centro de

Investigaciones y Estudios en Recursos Hídricos CIERHI-

EPN, desarrollando proyectos de investigación básica y

aplicada. Miembro la red de jóvenes investigadores IAHR y

editor asociado en la revista Hidrolatinoamericana de jóvenes

investigadores y profesionales.

64

Preparación de Artículos para la Revista Politécnica

_________________________________________________________________________________________________________________________

Revista Politécnica- XXXX 201X, Vol. XX, No. X

11. SECCIÓN I

Este documento es una plantilla para versiones Microsoft

Word 2013 o posteriores. Si está leyendo una versión impresa

de este documento, por favor descargue el archivo electrónico,

revistapolitécnicaformato2016.docx. En caso de que el autor

desee enviar el artículo en formato LaTex por favor

comunicarse con la coordinación de edición

([email protected]). Por favor, no coloque numeración

ni pie de página en el documento presentado.

No cambie los tamaños de fuente o espaciado de renglones

para ajustar el texto a un número limitado de páginas.

Utilice cursiva o negrita para dar énfasis a un texto, no

subrayado.

2. SECCIÓN II

Para las pautas de presentación, siga las instrucciones emitidas

por el sistema del sitio web de la revista de la EPN.

Colocar el correo electrónico del autor de correspondencia.

La presentación inicial debe tomar en cuenta todas las

indicaciones que se presentan en la plantilla, para de esta

manera tener una buena estimación de la longitud del artículo

a publicarse. Además, de esta manera el esfuerzo necesario

para la presentación final del manuscrito será mínimo.

Como sugerencia, es importante tomar en cuenta que, el primer

autor es el investigador que hizo la mayor parte del trabajo,

mientras que el último autor suele ser el profesor quien es el

líder intelectual y, a menudo edita y presenta el borrador final

del documento.

La Revista Politécnica pondrá en marcha un sistema de

transferencia electrónica de derechos de autor en su momento.

Por favor, "no" enviar formularios de derecho de autor por

correo o fax. A continuación se detallan las consideraciones

que se deben tener en cuenta para la presentación final del

artículo.

3. SECCIÓN III

Preparación de Artículos para la Revista Politécnica Utilizar

Mayúsculas en cada Palabra en el Caso del Título

Apellido, Nombre1; Apellido, Nombre2; Apellido, Nombre3

1Escuela Politécnica Nacional, Facultad de Ingeniería Mecatrónica, Quito, Ecuador

2Escuela Politécnica del Litoral, Facultad de Ingeniería Industrial, Guayaquil, Ecuador 3Universidad de Cuenca, Facultad de Ciencias Exactas, Cuenca, Ecuador

Resumen: Las siguientes instrucciones establecen las pautas para la preparación de artículos para la Revista

Politécnica. Los artículos pueden ser escritos en español o en inglés, pero tendrán un resumen de máximo 250 palabras

en los dos idiomas. Los autores pueden hacer uso de este documento como una plantilla para componer su artículo si

están utilizando Microsoft Word 2013 o superior. Caso contrario, este documento puede ser utilizado como una guía

de instrucciones. El número mínimo de páginas será 6 y el máximo 15, Para el envío de los artículos, los autores

deben seguir las instrucciones colocadas en el sistema de recepción de artículos del sitio web de la Revista Politécnica

(www.revistapolitecnica.epn.edu.ec). En caso de que su artículo sea en inglés colocar el título y el resumen en los

dos idiomas.

Palabras clave: Incluir una lista de 3 a 6 palabras.

Title of Manuscript

Abstract: These instructions give you guidelines for preparing papers for EPN Journal. Papers can be written in

Spanish or English; however, an abstract of maximum 250 words and written in both languages is required. Use this

document as a template to compose your paper if you are using Microsoft Word2013 or later. Otherwise, use this

document as an instruction set. The minimum number of pages will be 6 and the maximum will be 15. For submission

guidelines, follow instructions on paper submission system from the EPN Journal

website(www.revistapolitecnica.epn.edu.ec).

Keywords:Include a list of 3 to 6 words.


http://www.revistapolitecnica.epn.edu.ec/

http://www.revistapolitecnica.epn.edu.ec/

Apellido Nombre; Apellido Nombre; Apellido Nombre; Apellido Nombre

_______________________________________________________________________________________________________________________________


3.1 Figuras, tablas y márgenes

Todas las figuras deben ser incorporadas en el documento. Al

incluir la imagen, asegúrese de insertar la actual en lugar de un

enlace a su equipo local. Los archivos de: figuras, dibujos,

fotografías, etc., deberán enviarse en formato bmp o jpg, con

al menos 1200 puntos (resolución) en uno de sus ejes, con

leyendas legibles y de tamaño adecuado. El artículo debe

contener entre tablas y figuras un máximo de 10.

Las etiquetas de los ejes de las figuras son a menudo una fuente

de confusión. Utilice las palabras en lugar de símbolos. Por

ejemplo, escriba la cantidad "Magnetización," o

"Magnetización M" no sólo "M".

Las figuras y tablas deben estar en la parte superior e inferior

de las columnas. Evite colocarlas en medio de ellas. Las

figuras y tablas grandes pueden extenderse a lo largo de ambas

columnas. Las leyendas de las figuras deben estar centradas

debajo de las figuras, los títulos de las tablas deben estar

centrados sobre ellas. Evite colocar figuras y tablas antes de su

primera mención en el texto. Para la mención de figuras, tablas

o ecuaciones utilice las palabras completas con la primera letra

en mayúscula, por ejemplo "Figura 1".

Coloque las unidades entre paréntesis. No etiquete los ejes sólo

con unidades. Por ejemplo, escriba "Magnetización (A/m)" o

"Magnetización (Am-1)", no sólo "Magnetización A/m." No

etiquete los ejes con una relación de cantidades y unidades. Por

ejemplo, escriba "Temperatura (K)", no "Temperatura K".

Los multiplicadores pueden ser especialmente confusos.

Escriba "Magnetización (kA/m)" o "Magnetización

(103A/m)". No escriba "Magnetización (A/m) x 1000" porque

el lector no sabrá si la etiqueta del eje de arriba significa 16000

A/m o 0,016 A/m. Las etiquetas de las figuras deben ser

legibles, con un valor de 8 y sin espacio de separación con la

figura.

Figura 1. Distribución Weibull de 60 Hz voltajes de ruptura11 cables α =

45,9 kV picoβ = 5,08.Intervalo de Confidencia 95%

Los autores deben trabajar activamente con los márgenes

solicitados. Los documentos de la revista serán marcados con

los datos del registro de la revista y paginados para su inclusión

en la edición final. Si la sangría de los márgenes en su

manuscrito no es correcta, se le pedirá que lo vuelva a

presentar y esto, podría retrasar la preparación final durante el

proceso de edición.

Por favor, no modificar los márgenes de esta plantilla. Si está

creando un documento por su cuenta, considere los márgenes

que se enumeran en la Tabla 1. Todas las medidas están en

centímetros.

Tabla 1.Márgenes de página

Página Superior Inferior Izquierda/

Derecha Primera 2,0 2,5 1,5 Resto 2,0 2,5 1,5

3.2 Ecuaciones

Si está usando MSWord, sugerimos utilizar el Editor de

ecuaciones de Microsoft o el MathTypeadd-on para las

ecuaciones en su documento (Insertar/Objeto/Crear

Nuevo/Microsoft Ecuación o Ecuación MathType). La opción

"flotar sobre el texto" no se debe elegir.’

Enumere las ecuaciones consecutivamente con los números de

la ecuación en paréntesis contra el margen derecho, como en

(1). Utilice el editor de ecuaciones para crear la ecuación y esta

debe estar localizada en el margen derecho, como se muestra

en el ejemplo siguiente:

)]2(/[),( 020

2

rddrrFr

(1)

Asegúrese de que los símbolos en su ecuación han sido

definidos antes de que aparezcan en la ecuación o

inmediatamente después. Ponga en cursiva los símbolos (T

podría referirse a la temperatura, pero T es la unidad tesla).

Para referirse a la ecuación se escribe por ejemplo “Ecuación

(1) "

3.3 Unidades

Utilice el SI como unidades primarias. Otras unidades pueden

ser utilizadas como unidades secundarias (en paréntesis). Por

ejemplo, escriba "15 Gb/cm2 (100 Gb/in2)". Evite combinar las

unidades del SI y CGS, como la corriente en amperios y el

campo magnético en oerstedios. Esto a menudo lleva a

confusión porque las ecuaciones no cuadran

dimensionalmente. Si tiene que usar unidades mixtas, aclare

las unidades para cada cantidad en una ecuación.

Por ejemplo, en el SI la unidad de fuerza de campo magnético

Hes A/m. Sin embargo, si desea utilizar unidades de T, o bien

se refiere a la densidad de flujo magnético B o la fuerza del

campo magnético simbolizadas como µ0H. Use un punto en el

centro para separar las unidades compuestas, por ejemplo,

“A·m2.”

3.4 Abreviaturas y Siglas

Defina las abreviaciones y acrónimos la primera vez que se

utilizan en el texto, incluso después de que ya han sido

Breakdown Voltage (kV)

100 101 102

0.2

0.1

2

20

70

90

98

99.9

50

Wei

bull

Bre

akdo

wn

Pro

babi

lity

(%)

30

10

5

1

0.5

Preparación de Artículos para la Revista Politécnica

_________________________________________________________________________________________________________________________


definidos en el resumen. No utilice abreviaturas en el título a

menos que sea inevitable.

3.5 Otras recomendaciones

Para expresar valores decimales se usarán comas, por

ejemplo 3,45. Use un cero antes del decimal.

Se incluirá un espacio entre números para indicar los

valores de miles, por ejemplo 463 690.

Utilice notación científica para expresar números con

más de 3 cifras hacia la derecha o izquierda, es decir,

mayores a 2,50E+05 o menores a 4,8E-03.

Finalmente, de ser necesario y de manera opcional, se

pueden incluir conclusiones, recomendaciones y

agradecimiento.

REFERENCIAS

La lista de referencias debe estar en Formato APA

ordenada alfabéticamente de acuerdo con el apellido del

primer autor del artículo. El agregado et al no debe ir en

cursiva. Por favor nótese que todas las referencias listadas aquí

deben estar directamente citadas en el cuerpo del texto usando

(Apellido, año). Las notas al pie deben evitarse en la medida

de lo posible.

El artículo debe contener un mínimo de 6 referencias.

Seguir el formato indicado a continuación de acuerdo al tipo

de referencia a:

Formato básico para referenciar libros:

Apellido, Inicial Nombre. (Año). Título del libro. Ciudad,

País: Editorial.

Libros con un autor:

En las referencias: King, M. (2000). Wrestling with the angel: A life of Janet Frame. Auckland,

New Zealand: Viking.

Cita en el texto:

(King, 2000) o King (2000) argumenta que ...

Libros con dos autores:

En las referencias: Treviño, L. K., y Nelson, K. A. (2007). Managing business ethics: Straight

talk about how to do it right. Hoboken, NJ: Wiley

Cita en el texto:

(Treviño y Nelson, 2007) oTreviño y Nelson (2007)

ilustran…

Libros con dos o más autores:

En las referencias: Krause, K.-L., Bochner, S., y Duchesne, S. (2006). Educational psychology

for learning and teaching (2nd ed.). South Melbourne, VIC., Australia:

Thomson.

Cita en el texto:

De acuerdo con Mezey et al. (2002) o ... (Mezey et al.,

2002).

Formato básico para referenciar artículos científicos

Apellido, Inicial Nombre. (Año). Título del Artículo.

Título/Iniciales de la Revista. Número de Volumen (Tomo),

páginas

Artículos en revistas:

En las referencias: Sainaghi, R. (2008). Strategic position and performance of winter

destinations. TourismReview, 63(4), 40-57.

Cita en el texto:

(Sainaghi, 2008) oSainaghi (2008) sugiere ...

Artículos con DOI

En lasreferencias: Shepherd, R., Barnett, J., Cooper, H., Coyle, A., Moran-Ellis, J., Senior, V.,

& Walton, C. (2007). Towards an understanding of British public attitudes concerning human cloning. Social Science& Medicine, 65(2), 377-392.

http://dx.doi.org/10.1016/j.socscimed.2007.03.018

Cita en el texto:

Shepherd et al. (2007) o Shepherd et al. (2007) resaltan la...

Artículos sin DOI

En las referencias Harrison, B., & Papa, R. (2005). The development of an indigenous

knowledge program in a New Zealand Maori-language immersion

school. Anthropology and EducationQuarterly, 36(1), 57-72. Obtenido de la base de datos AcademicResearch Library

Cita en el texto:

(Harrison y Papa, 2005) o En su investigación, Harrison y

Papa (2005) establecieron...

Artículos en línea

En lasreferencias: Snell, D., & Hodgetts, D. (n.d.). The psychology of heavy metal communities

and white supremacy. Te KuraKeteAronui, 1. Obtenido de: http://www.waikato.ac.nz/wfass/tkka. (Mayo, 2015).

Cita en el texto:

(Snell y Hodgetts, n.d.) oSnell y Hodgetts (n.d.) identificaron

"..."

http://dx.doi.org/10.1016/j.socscimed.2007.03.018

http://www.waikato.ac.nz/wfass/tkka

Apellido Nombre; Apellido Nombre; Apellido Nombre; Apellido Nombre

_______________________________________________________________________________________________________________________________


INFORMACIÓN ADICIONAL

Sistema de Arbitraje:

Todos los artículos cumplen con una revisión por pares, la cual

consiste en:

Selección de dos o tres árbitros, actualmente la Revista

Politécnica cuenta con revisores internos, externos e

internacionales, quienes envían al editor

su evaluación del artículo y sus sugerencias acerca de

cómo mejorarlo.

El editor reúne los comentarios y los envía al autor

Con base en los comentarios de los árbitros, el editor

decide si se publica el manuscrito.

Cuando un artículo recibe al mismo tiempo evaluaciones

tanto muy positivas como muy negativas, para romper

un empate, el editor puede solicitar evaluaciones

adicionales, obviamente a otros árbitros.

Toda la evaluación se realiza en un proceso doble ciego,

es decir los autores no conocen quienes son sus

revisores, ni los revisores conocen los autores del

artículo.

Instructivo para publicar un Artículo

1. Crear un usuario y contraseña para acceder al portal

web de la Revista Politécnica, para mayor

información está el correo [email protected]

2. Ingresar al portal web e iniciar el proceso de envío

3. Comenzar el envío

4. Colocar requisitos de envío

Lista de comprobación de preparación de envíos

Como parte del proceso de envío, se les requiere a los

autores que indiquen que su envío cumpla con todos

los siguientes elementos, y que acepten que envíos

que no cumplan con estas indicaciones pueden ser

devueltos al autor.

- La petición no ha sido publicada previamente, ni

se ha presentado a otra revista (o se ha

proporcionado una explicación en Comentarios

al Editor).

- El fichero enviado está en formato OpenOffice,

Microsoft Word, RTF, o WordPerfect.

- Se han añadido direcciones web para las

referencias donde ha sido posible.

- El texto tiene interlineado simple; el tamaño de

fuente es 10 puntos; se usa cursiva en vez de

subrayado (exceptuando las direcciones URL); y

todas las ilustraciones, figuras y tablas están

dentro del texto en el sitio que les corresponde y

no al final del todo.

- El texto cumple con los requisitos bibliográficos

y de estilo indicados en las Normas para

autoras/es, que se pueden encontrar en "Acerca

de la Revista".

Nota de copyright

Los autores que publican en esta revista están de

acuerdo con los siguientes términos:

- Los autores conservan los derechos de autor y

garantizan a la revista el derecho de ser la

primera publicación del trabajo al igual que

licenciado bajo una Creative Commons

Attribution License que permite a otros

compartir el trabajo con un reconocimiento de la

autoría del trabajo y la publicación inicial en esta

revista.

- Los autores pueden establecer por separado

acuerdos adicionales para la distribución no

exclusiva de la versión de la obra publicada en la

revista (por ejemplo, situarlo en un repositorio

institucional o publicarlo en un libro), con un

reconocimiento de su publicación inicial en esta

revista.

- Se permite y se anima a los autores a difundir sus

trabajos electrónicamente (por ejemplo, en

repositorios institucionales o en su propio sitio

web) antes y durante el proceso de envío, ya que

puede dar lugar a intercambios productivos, así

como a una citación más temprana y mayor de

los trabajos publicados (Véase The Effect of

Open Access) (en inglés).

Declaración de privacidad

- Los nombres y direcciones de correo-e

introducidos en esta revista se usarán

exclusivamente para los fines declarados por esta

revista y no estarán disponibles para ningún otro

propósito u otra persona.

5. Subir el envío

6. Introducir metadatos

7. Subir ficheros adicionales

8. Confirmar el envío

http://es.wikipedia.org/wiki/Evaluaci%C3%B3n


TEMÁTICA - revistapolitecnica.epn.edu.ec · Se autoriza la reproducción total o parcial de su...

Documents

Transcript of TEMÁTICA - revistapolitecnica.epn.edu.ec · Se autoriza la reproducción total o parcial de su...