Caracterización resistiva 3D del área geotérmica de Paipa, Colombia

Carlos Eduardo González-Idárraga

Caracterización resistiva 3D del área geotérmica de Paipa, Colombia

Boletín de Geología, vol. 42, núm. 3, 2020

Universidad Industrial de Santander

Carlos Eduardo González-Idárraga

Servicio Geológico Colombiano, Colombia


Recibido: 10 Febrero 2020

Aceptado: 19 Agosto 2020

Resumen: Con el fin de realizar la caracterización resistiva del área geotérmica de Paipa, se llevaron a cabo diferentes campañas de adquisición de sondeos audiomagnetotelúricos (AMT) y magnetotelúricos (MT). La estimación de las curvas de resistividad aparente se encuentra en un rango de periodo entre 0,0001 y 10 s. El análisis de dimensionalidad geoeléctrica realizado a los tensores de impedancia estableció un comportamiento 1D en bajos periodos y un comportamiento 2D-3D dominante en periodos altos. Dado lo anterior, se realizó una inversión 3D mediante el código ModEM y se interpretó bajo el modelo conceptual del sistema geotérmico. En inmediaciones de los domos Alto Los Volcanes y Alto Los Godos, donde se sugiere la ubicación del reservorio, se identificó la presencia de un cuerpo conductor que se asocia a la posible capa sello del sistema. Adicionalmente, se observaron una serie de conductores al centro y NO del área, bajo aparente control estructural, asociados al tránsito de fluidos hasta las zonas de descarga en el cruce de fallas; así como la identificación de un cuerpo resistivo asociado al basamento al sureste del área.

Palabras clave: Método magnetotelúrico, Inversión 3D, Exploración geotérmica, Paipa, Colombia.

Abstract: In order to performe the electrical characterization of the geothermal area of Paipa, different campaigns were carried out to acquire audiomagnetotellurics (AMT) and magnetotelluric soundings (MT). The estimation of the apparent resistivity curves is in a period range between 0.0001 and 10 s. The geoelectric dimensionality analysis performed to the impedance tensors established a 1D behavior in low periods and a dominant 2D-3D behavior in high periods. In the light of this, a 3D inversion of the area’s electrical structure was performed using the ModEM code and interpreted under the conceptual model of the geothermal system. The presence of a conductive body was identified in the vicinity of the Alto Los Volcanes and Alto Los Godos domes, where the reservoir location has been previously suggested, and is associated with the possible seal layer of the system. Moreover, a series of conductive bodies in apparent structural control, were observed in the center and NW of the area. These bodies are associated with the transit of fluids to the discharge zones at fault intersections. Finally, a resistive body associated with the basement to the southeast of the area was also identified.

Keywords: Magnetotelluric method, 3D inversion, Geothermal exploration, Paipa, Colombia.

Introducción

La fase exploratoria de un proyecto geotérmico se encuentra caracterizada por los estudios de reconocimiento, prefactibilidad y factibilidad (OLADE y BID, 1994); donde se conjugan diversas técnicas geológicas, geoquímicas y geofísicas en la identificación de las estructuras características de un sistema geotérmico (fuente de calor, reservorio permeable, capa sello, zonas de recarga, etc.). Los métodos geofísicos se implementan para la identificación indirecta de las condiciones estructurales y la geometría del sistema. Particularmente, el método magnetotelúrico (MT) ha sido ampliamente utilizado en la caracterización resistiva de diversos ambientes geotermales (Volpi et al., 2003; Abdul-Azeez y Harinarayana, 2007; Arango et al., 2009; Cumming y Mackie 2010; Peacock et al., 2012); ya que la resistividad eléctrica del subsuelo está en gran medida controlada por la porosidad, permeabilidad, saturación del medio, salinidad, presión y temperatura del medio (Chave y Jones, 2012; Simpson y Bahr, 2005).

Trabajos como los de Spichak y Manzella (2009), Muñoz (2014) y Van Leeuwen (2016) muestran una recopilación de estudios magnetotelúricos realizados en diferentes ambientes geotérmicos, así como las posibles firmas resistivas asociadas a estos. Por ejemplo, los sistemas hidrotermales de alta temperatura (Cumming y Mackie 2010; Volpi et al., 2003), presentan una capa sello conductora dada por la presencia de arcillas sobre un reservorio geotermal más resistivo, reflejando un aumento de la temperatura con la profundidad (Ussher et al., 2000). Sin embargo, en algunos sistemas geotérmicos asociados a ambientes volcánicos la capa sello puede estar ausente (Muñoz, 2014). Por otro lado, en sistemas geotérmicos nomagmáticos (Abdul-Azeez y Harinarayana, 2007; Arango et al., 2009; Peacock et al., 2012) existe una amplia gama de configuraciones geométricas, lo que hace difícil la identificación de patrones resistivos a seguir (Muñoz, 2014). Cada sistema debe ser analizado rigurosamente con el apoyo de las diferentes técnicas de estudio.

El área geotérmica de Paipa, ubicada al noreste de Colombia (Figura 1), fue catalogada con una prioridad media-alta a partir de estudios de geología de superficie y de la geoquímica de la fase acuosa de las fuentes termales (OLADE, GEOTÉRMICA ITALIANA, ICEL y CONTECOL, 1982). Por otro lado, un estudio de prefactibilidad sugirió la existencia de un reservorio geotermal con potencial para la generación de energía eléctrica (Japan Consulting Institute, 1983). Dado lo anterior, el Servicio Geológico Colombiano ha realizado diversos estudios geológicos, geoquímicos y geofísicos, encaminados a la caracterización del sistema geotérmico de Paipa.

El área de estudio se encuentra en un ambiente sedimentario con presencia de depósitos piroclásticos y cuerpos efusivos (domos riolíticos) (Cepeda y Pardo, 2004; Rueda-Gutiérrez, 2017). Los manantiales termales presentes en la zona han registrado temperaturas máximas entre 70-76°C (Alfaro, 2002a; Ortiz y Alfaro, 2010). De acuerdo a la localización de las manifestaciones termales (Figura 1), estas estarían controladas estructuralmente por cruces de fallas en el área de influencia de formaciones sedimentarias (Velandia, 2003; Alfaro et al., 2012), permitiendo el ascenso de fluidos. Adicionalmente, Alfaro (2002a; 2002b) establece la presencia de sulfato sódico en el análisis de las muestras de los manantiales, proveniente de evaporitas al oeste de la secuencia sedimentaria, el cual se mezclaría con los fluidos geotérmicos, enmascarando la química del fluido hidrotermal y dificultando la correcta estimación de la temperatura en el reservorio. Sin embargo, a partir de geotermómetros de sílice, se ha estimado un reservorio con temperaturas entre 120°C y 230°C (Alfaro, 2002b).

Llanos et al. (2015), a partir de modelación geológicogeofísica con información de métodos potenciales (Vásquez, 2012), identificaron anomalías de alta densidad correspondientes a basamento profundo e intrusiones en la zona central de trabajo (domos Alto Los Volcanes y Alto Los Godos, en la Figura 1), así como anomalías de baja densidad asociadas a la secuencia sedimentaria al centro-oeste del área.

El objetivo de este trabajo es realizar la caracterización resistiva del área geotérmica de Paipa, a partir de la inversión 3D de datos magnetotelúricos, utilizando el código ModEM (Egbert y Kelbert, 2012; Kelbert et al., 2014), e identificar las principales estructuras resistivas del sistema geotérmico.

Mapa geológico y estructural del área geotérmica de Paipa. Unidades del Cuaternario: Depósitos piroclásticos (DP), Depósitos del Cuaternario (DC). Unidades del Neógeno: Brecha hidrotermal (BH), Formación Tibatá (FT). Unidades del Paleógeno: Formación Bogotá (FB). Unidades del Cretácico: Formación Labor y Tierna (FLyT), Formación Guaduas (FG), Formación Los Pinos (FLP), Formación Plaeners (FP), Formación Conejo (FC), Formación Churuvita (FCh), Formación Une (FU), Formación Tibasosa (FTi). ATT: Anticlinal Tibasosa-Toledo. Modificado de Velandia y Cepeda (2004), Pardo (2004) y Rueda-Gutiérrez (2017).
Figura 1
Mapa geológico y estructural del área geotérmica de Paipa. Unidades del Cuaternario: Depósitos piroclásticos (DP), Depósitos del Cuaternario (DC). Unidades del Neógeno: Brecha hidrotermal (BH), Formación Tibatá (FT). Unidades del Paleógeno: Formación Bogotá (FB). Unidades del Cretácico: Formación Labor y Tierna (FLyT), Formación Guaduas (FG), Formación Los Pinos (FLP), Formación Plaeners (FP), Formación Conejo (FC), Formación Churuvita (FCh), Formación Une (FU), Formación Tibasosa (FTi). ATT: Anticlinal Tibasosa-Toledo. Modificado de Velandia y Cepeda (2004), Pardo (2004) y Rueda-Gutiérrez (2017).


Marco geológico

El área de estudio se encuentra principalmente en un ambiente sedimentario al sur del municipio de Paipa (Boyacá, Colombia). De acuerdo a la caracterización geológica y estructural realizada por Velandia (2003), al NO y SO del área de estudio (Figura 1) se encuentran sedimentos del Cuaternario y del Neógeno (Formación Tilatá), representados por secuencias de arenas, limos arcillas, conglomerados correspondientes a la actividad más reciente de tipo aluvial, lacustre y fluvio-lacustre. En el sector central y oriental, se encuentran rocas sedimentarias del Paleógeno (Formación Bogotá) y Cretácico (Formaciones Guaduas, Labor y Tierna, Los Pinos, Plaeners, Conejo, Churuvita Une y Tibasosa), constituidas por intercalaciones de areniscas cuarzosas, limolitas, arcillolitas, niveles esporádicos de calizas e impregnaciones de sales y lutitas. El espesor de esta secuencia sedimentaria se estima entre 1 – 1,5 km aproximadamente al NO del área (Velandia, 2003). Bajo la secuencia sedimentaria se encuentra el basamento formado por rocas sedimentarias (niveles arenosos fracturados, intercalados con arcillolitas y limolitas) y metamórficas (filitas, esquistos y gnéises) del Paleozoico. Dicho basamento correspondería a la continuación en profundidad del Macizo de Floresta, el cual aflora fuera del área de estudio al NE, formando el núcleo del Anticlinal Tibasosa-Toledo (Figura 1). Adicionalmente, en la Tabla 1 se puede observar una breve descripción litológica de cada una de las formaciones geológicas, así como su importancia hidrogeológica.

Tabla 1
Descripción litológica de las formaciones sedimentarias (Velandia, 2003) y su clasificación hidrogeológica (SGC, 2016).
Descripción litológica de las formaciones sedimentarias (Velandia, 2003) y su clasificación hidrogeológica (SGC, 2016).


Sobre la cubierta sedimentaria, Cepeda y Pardo (2004) observaron una secuencia de unidades piroclásticas correspondientes a dos etapas eruptivas posteriores al Eoceno tardío, con intercalaciones de sedimentitas de ambientes fluvial, lacustre y coluvial (Figura 1). Por otro lado, al oeste del área, Velandia (2003) describe la presencia de una brecha hidrotermal (Figura 1) caracterizada por fragmentos muy alterados de posible origen ígneo, posterior a los eventos volcánicos, con un alto grado de alteración argílica avanzada (González et al., 2008; Rodríguez y Alfaro, 2015), sugiriendo un posible aporte de calor radiogénico al sistema hidrotermal (Alfaro et al., 2017), producto de una intrusión granítica en la zona.

Rueda-Gutiérrez (2017) caracterizó y cartografió tres sectores con cuerpos dómicos de origen de arco volcánico, denominados Alto Los Volcanes, Alto Los Godos y Quebrada Honda; con edades estimadas de 1,76 y 1,8 Ma para el Alto Los Volcanes y Quebrada Honda, respectivamente; y 2,8, 2,7 y 2,6 Ma para el sector del Alto Los Godos. Sin embargo, refiere la edad más reciente de estas intrusiones en 1 millón de años a partir de dataciones complementarias por trazas de fisión.

El área se encuentra deformada por fallas regionales con dirección NE-SO, asociadas a la apertura de una cuenca marina en el Cretácico y la inversión tectónica consecuencia de un régimen compresivo que dio lugar al levantamiento de la Cordillera Oriental en el Neógeno (Ortiz, 2014). Velandia (2003) observó una serie de estructuras longitudinales y transversales a los lineamientos regionales (Figura 1). Para las fallas longitudinales con orientación NNE-SSO (El Bizcocho, El Batán, Rancho Grande, Buena Vista, Agua Tibia, Lanceros) describe sus trazos como sinuosos y localmente discontinuos por estar cubiertos por depósitos neógenos y cuaternarios, tanto sedimentarios como vulcano-clásticos, dificultando la interpretación de relaciones estructurales entre las unidades preneógenas (Velandia, 2003). Las fallas transversales con orientación NO-SE (Cerro Plateado, Paipa-Iza, Las Peñas) controlan morfología y drenaje, posiblemente facilitando el paso de fluidos hidrotermales (Alfaro et al., 2012), y permitiendo el ascenso de los magmas que dieron origen al volcanismo reconocido en la zona. Entre tanto, las fallas transversales con orientación NESO (El Hornito, Canocas, Santa Rita, El Tuno) afectan sólo la secuencia sedimentaria (Velandia, 2003).

Metodología

El modelo resistivo del área geotérmica se estimó por medio de la inversión de datos magnetotelúricos. El método magnetotelúrico es un método geofísico de fuente pasiva, el cual registra en superficie las variaciones temporales del campo electromagnético natural, por medio del cual se puede estimar la distribución de la resistividad del subsuelo (Vozoff, 1991) desde decenas de metros hasta centenas de kilómetros. A partir de la relación entre las componentes horizontales del campo magnético y eléctrico, en el dominio de las frecuencias, se halla un tensor de impedancia complejo (Z) de 2x2, cuya expresión matemática está dada por:

Cada una de las componentes de . es transformada en resistividad aparente (ρa) y fase (φ), que posteriormente son modeladas.

donde μ0 es la permeabilidad magnética del vacío y . es la frecuencia angular.

Adicionalmente, la medida de la componente vertical del campo magnético y su relación con las componentes horizontales del mismo campo, genera la función de transferencia geomagnética (VTF en inglés). La cual permite estimar las variaciones laterales de la resistividad en el medio (Vozoff, 1972).

En el área de estudio se adquirieron 78 sondeos audiomagnetotelúricos y magnetotelúricos (AMT+MT) (Figura 2) a lo largo de diferentes campañas entre los años 2012 y 2016 (Moyano y Vallejo, 2015; González-Idárraga y Rodríguez-Rodríguez, 2017). La información se adquirió con equipos receptores de la marca Phoenix Geophysics y sensores magnéticos con capacidad de registro hasta 10000 s. Los dipolos eléctricos y bobinas fueron orientadas al norte para el eje x y al este para el eje y. Los tiempos de medición fueron de una hora para las mediciones AMT y no menos de 12 horas para las mediciones MT. A partir de un procesamiento robusto de las señales temporales adquiridas, bajo un esquema de referencia remota local (Gamble et al., 1979), se obtuvieron las respectivas curvas de resistividad aparente, fase y VTF.

En la Figura 3 se presentan las curvas de resistividad aparente, fase y VTF de 6 sondeos representativos de las adquisiciones realizadas en el área de estudio. Estos muestran buenas estimaciones hasta periodos de 10 s para todas las componentes del tensor de impedancia, mientras que en 17 sondeos, la componente Zyx alcanza estimaciones iguales o mayores a 100 s (PA0064, PA0077, PA0130, PA0077, PA0145, PA0150, PA0160, PA0162, entre otros. Ubicaciones en la Figura 2). Entretanto, las curvas de VTF presentan buenas estimaciones hasta periodos entre 1 y 10 s. También se identificó en algunas mediciones una baja estimación en el rango de banda muerta del AMT, como se observa en las componentes de VTF de los sondeos PA0020, PA0073 y PA0098. Para periodos mayores a 10 s las curvas presentaron mayores errores y dispersión en los datos, por lo que no se tuvieron en cuenta al momento de realizar la inversión.

Finalmente, sondeos como PA0020, PA0027 y PA0077 presentaron corrimiento estático (static shift, en inglés) en las curvas de resistividad aparente. Sin embargo, la inversión 3D puede corregir la distorsión galvánica que genera dicho corrimiento, así como los efectos inductivos debido a estructuras cercanas a la superficie (Siripunvaraporn, 2012).

Ubicación de 78 sondeos magnetotelúricos (círculos negros) en el área de interés geotérmica de Paipa. Sondeos resaltados con circunferencia verde son representados en la Figura 3.
Figura 2
Ubicación de 78 sondeos magnetotelúricos (círculos negros) en el área de interés geotérmica de Paipa. Sondeos resaltados con circunferencia verde son representados en la Figura 3.


Ejemplos de curvas de resistividad aparente y de fase de las componentes Zxy (curva roja) y Zyx (curva azul), y VTF de la zona de estudio. De forma general, se observan buenas estimaciones de las curvas hasta un periodo aproximado de 10 s. La ubicación de los sondeos mostrados se encuentra resaltada en la Figura 2.
Figura 3
Ejemplos de curvas de resistividad aparente y de fase de las componentes Zxy (curva roja) y Zyx (curva azul), y VTF de la zona de estudio. De forma general, se observan buenas estimaciones de las curvas hasta un periodo aproximado de 10 s. La ubicación de los sondeos mostrados se encuentra resaltada en la Figura 2.


Análisis de dimensionalidad

Con ánimo de estimar la dimensionalidad geoeléctrica del medio, se realizó un análisis a partir del Tensor de Fase (Caldwell et al., 2004). El tensor de fase está definido a partir de las relaciones de fase entre las componentes horizontales de los campos eléctrico y magnético. La representación gráfica de las cantidades escalares valor máximo del tensor (φmax), valor mínimo del tensor (φmin) y el ángulo de inclinación skew β (medida de la asimetría del tensor) reflejan la dimensionalidad del medio, así como la dirección de strike del contraste de conductividad (Caldwell et al., 2004). En la Figura 4 se puede observar la representación gráfica del tensor de fase, donde las circunferencias definen un medio 1D con un ángulo skew de cero, un medio 2D es representado por una elipse con valores de skew de cero o cercanos (±3°) y valores diferentes para el skew reflejan un medio 3D (Caldwell et al., 2004; Bibby et al., 2005; Booker, 2014).

Se representaron los tensores de fase para los periodos 0,01, 0,1, 1 y 10 s (Figura 4). En el periodo 0,01 s se observa, de forma general, circunferencias con valores de skew cercanos a cero indicando una dimensionalidad 1D en las zonas más someras del área de estudio. Este comportamiento es esperado, al encontrarse la zona de estudio en un ambiente sedimentario. A periodos más altos (0,1, 1 y 10 s) se observan dimensionalidades 2D y 3D. Por otra parte, las flechas negras representan los vectores de inducción (Parkinson, 1959) que proporcionan información de contrastes en la conductividad del medio. Este parámetro refleja un gradiente conductivo hacia el NO del área de estudio a un periodo de 1 s, apoyado por la orientación del eje mayor de las elipses del tensor de fase en la misma dirección.

Análisis de dimensionalidad a partir del tensor de fase (Caldwell et al., 2014) a diferentes periodos, las elipses con código de color representan la distribución geoeléctrica del medio. Las flechas negras representan los vectores de inducción, de acuerdo a la convención de Parkinson (1959) apuntan hacia el conductor.
Figura 4
Análisis de dimensionalidad a partir del tensor de fase (Caldwell et al., 2014) a diferentes periodos, las elipses con código de color representan la distribución geoeléctrica del medio. Las flechas negras representan los vectores de inducción, de acuerdo a la convención de Parkinson (1959) apuntan hacia el conductor.


Inversión 3D

El algoritmo de inversión magnetotelúrica 3D ModEM (Egbert y Kelbert, 2012; Kelbert et al., 2014) fue empleado para generar el modelo de resistividad 3D del área de estudio. Se realizó la inversión conjunta del tensor de impedancia y de la función de transferencia geomagnética de 78 sondeos, previamente depurados en los periodos con una baja calidad. Se estableció un valor de 0,3 de covarianza. Se implementó un modelo inicial con una resistividad de 100 ohm.m, un mallado horizontal de 33 celdas en el eje N-S, 30 celdas en el eje E-O, dimensiones de 400 m x 400 m y 10 celdas de atenuación a cada lado con un incremento logarítmico de 1,3. El eje vertical se discretizó en 79 celdas. Entre el sondeo realizado a mayor altura y el de menor altura, las celdas poseen una resolución de 50 m y posteriormente se realizó un incremento logarítmico de 1,2. Adicionalmente, los sondeos se proyectaron a la superficie acorde a la topografía y se centraron en cada celda.

La inversión fue realizada con 18 periodos del tensor de impedancia (Zxx, Zxy, Zyx, Zyy) y el VTF, en un rango desde 0,001 s hasta 10 s. Se realizaron pruebas de inversión de sólo Z, sólo VTF y Z+VTF. La secuencia se inició invirtiendo el tensor de impedancia con errores de piso de hasta 50%, posteriormente se utilizó el modelo resultante como modelo inicial de una nueva inversión con errores de piso menores, esto con el ánimo de ir evaluando los ajustes de las anomalías en el modelo (Romero et al., 2019). Finalmente se definieron errores de piso en 10% para la antidiagonal de Z, 20% para la diagonal de Z y un valor de 3% para el VTF, obteniendo un modelo resistivo del área con un valor RMS de 1,68 después de 118 iteraciones.

La Figura 5 muestra curvas de ajuste del modelo resistivo calculado. De manera general, se detalla un ajuste coherente para las resistividades aparentes de las componentes Zxy y Zyx del tensor de impedancia y para las curvas de VTF, en el rango de periodo modelado. Las componentes de la diagonal principal del tensor de impedancia (Zxx y Zyy) presentan desajustes superiores. Estas componentes están más sujetas a errores y se les ha asignado un error de piso superior. Sin embargo, se recomienda la inversión de todas las componentes del tensor de impedancia, ya que las estructuras pueden ser mejor cartografiadas (Kiyan et al., 2014).

Resultados

Los resultados de la modelación 3D realizada se muestran en las Figuras 6 y 7. Al oeste del área se observa la presencia de una serie de anomalías conductoras someras (C1, C2 y C3) que no superan los 500 m de profundidad. Dichas anomalías son reflejo del análisis de los vectores de inducción que mostraron un gradiente de conductividad al NO de la zona (Figura 4).

Al centro-este y sur del área modelada (Figuras 6 y 7) se encuentra una respuesta resistiva (R1) de grandes dimensiones, por debajo de 500 metros de profundidad para las secciones al sur. Sobre dicha anomalía resistiva se encuentra una capa de resistividad intermedia (entre 10 y 100 ohm.m), aparentemente con una mayor profundidad al oeste del área. Sin embargo, como se observa a partir del análisis de sensibilidad presentado en la siguiente sección, dado el rango de periodos modelado, al norte y noroeste del área no es posible asegurar la caracterización de las estructuras por debajo de 1470 msnm.

Análisis de sensibilidad

Al modificar los valores de resistividad a diferentes profundidades, hasta la base del modelo, y correr modelos directos para observar el desajuste en los periodos altos respecto al de las curvas calculadas en el modelo final, se puede estimar la profundidad de estudio alcanzada por el modelo.

El área modelada tiene una altura media de 2780 msnm, tomando este valor como cero, se modificó la resistividad en un valor de 1000 ohm.m a profundidades de -3150 msnm, -2230 msnm, -1470 msnm, -100 msnm y +1010 msnm. Adicionalmente, se realizaron análisis de sensibilidad con modelos escalonados (Figura 8A), donde el resistivo se encontraba más cercano a la superficie al norte y oeste del área modelada, dado que en dicha zona se encuentran cuerpos conductores superficiales y, al noroeste, sondeos MT con la menor longitud en términos de periodo, como se describió anteriormente. Al observar los diferentes desajustes en los sondeos se concluyó que, de manera general los ajustes a las diferentes profundidades, el modelo posee una profundidad de estudio de hasta la cota de -1470 msnm al centro-este y toda la zona sur del área (sondeos PA0020, PA0077 y PA0098, en la Figura 8B), mientras que al norte y noroeste solo tiene buena precisión hasta la cota +1010 msnm (sondeo PA0073 en Figura 8B).

Para el análisis de sensibilidad de los cuerpos conductores C1, C2 y C3, se modificó su resistividad a valores de 30 ohm.m (Figura 8C) y se realizó una inversión utilizando este como modelo inicial. En la Figura 8C, en la sección iii, se observa como los conductores se recuperan, confirmando la presencia de estos dentro del modelo.

Ejemplo del ajuste de las curvas de resistividad aparente y de fase de las componentes Zxy (curva roja) y Zyx (curva azul), y VTF de la zona de estudio. Datos observados (puntos), datos calculados (líneas).
Figura 5
Ejemplo del ajuste de las curvas de resistividad aparente y de fase de las componentes Zxy (curva roja) y Zyx (curva azul), y VTF de la zona de estudio. Datos observados (puntos), datos calculados (líneas).


Vistas horizontales del modelo resistivo 3D a diferentes profundidades. Se observa el presunto control estructural sobre las anomalías conductoras al centro-oeste del área (vistas A, B y C), así como la presencia de un cuerpo resistivo al sureste del área. Las vistas son referidas sobre el nivel del mar. Los puntos negros hacen referencia a la ubicación de los sondeos MT. Las líneas continuas y punteadas corresponden a las fallas en la Figura 1. La zona en blanco de la vista A hace referencia a celdas sin información (aire), dada por la topografía en el modelo. El área opacada al norte, en las vistas E y F, corresponde a una zona no sensible del modelo.
Figura 6
Vistas horizontales del modelo resistivo 3D a diferentes profundidades. Se observa el presunto control estructural sobre las anomalías conductoras al centro-oeste del área (vistas A, B y C), así como la presencia de un cuerpo resistivo al sureste del área. Las vistas son referidas sobre el nivel del mar. Los puntos negros hacen referencia a la ubicación de los sondeos MT. Las líneas continuas y punteadas corresponden a las fallas en la Figura 1. La zona en blanco de la vista A hace referencia a celdas sin información (aire), dada por la topografía en el modelo. El área opacada al norte, en las vistas E y F, corresponde a una zona no sensible del modelo.


Discusión

Las anomalías conductoras al oeste del área, como se observa en la Figura 6, parecen tener un control estructural dado por las fallas de bajo ángulo de buzamiento. El cuerpo conductor C1, formaría un aparente conducto de fluidos a través de las formaciones Guaduas, Los Pinos y particularmente Plaeners (al ser acuífero), con un posible control a lo largo de las fallas El Bizcocho y El Batán, terminando al norte con la descarga en los manantiales del sector de La Playa, los cuales han registrado temperaturas de hasta 76°C (Alfaro, 2002a; Ortiz y Alfaro, 2010). Por su parte, el cuerpo conductivo C2 representaría tránsito de fluidos termales a través de los depósitos Cuaternarios, a lo largo de las fallas Canocas y El Hornito, hasta su descarga en los manantiales del sector de ITP-Lanceros (temperatura máxima de descarga de 71°C). Dicha descarga se da en el cruce de fallas de El Hornito y El Bizcocho al norte del área de estudio (Figura 1). Estas dos anomalías (C1 y C2) parecieran estar conectadas, como se observa en las vistas en planta a 2484 y 2258 msnm (Figura 6). Por otro lado, de acuerdo a la extensión de la anomalía C2, se esperaría que el manantial de SALPA tuviera aportación de fluidos hidrotermales como en el sector de ITP-Lanceros. Sin embargo, este manantial posee altas concentraciones de sulfato de sodio y una temperatura baja (Alfaro, 2002a; Ortiz y Alfaro, 2010). Alfaro et al. (2017) proponen que este sector seria la descarga de fluidos provenientes de niveles inferiores de la secuencia sedimentaria, los cuales originarían disolución de evaporitas al oeste del área. Adicionalmente, dados los cortes superficiales del modelo (Figura 1 y Figura 6), no se observa una diferenciación resistiva entre las formaciones geológicas, ya que estas presentan grandes similitudes litológicas e hidrogeológicas.

La anomalía resistiva R1 hace referencia a la respuesta del basamento metamórfico al sureste del área, donde muestra una mayor presencia. Por otro lado, las respuestas resistivas R2 y R3 (Figura 6) se asocian a las respuestas superficiales de los sectores dómicos Alto Los Volcanes y Alto Los Godos, respectivamente.

La capa resistividad intermedia (entre 10 y 100 ohm.m), sobre la anomalía resistiva del basamento, se ha interpretado como el bloque de formaciones sedimentarias (González-Idárraga y Rodríguez- Rodríguez, 2017). Sin embargo, dado el análisis de sensibilidad, al noroeste del área no permitiría caracterizar las estructuras por debajo de las anomalías conductoras. Por lo anterior, esto no permite verificar la respuesta resistiva de la fuente salada del manantial de SALPA.

Secciones O-E del modelo resistivo 3D. Las profundidades se refieren a msnm. Los cortes A y B muestran las anomalías conductoras superficiales asociadas a los sitios de descarga de los manantiales termales. En el corte C se observa la ubicación del posible reservorio geotermal y de la capa sello. A lo largo de todos los cortes, se observa el basamento (R1) con una mayor presencia al SE del área. Los números dentro de las secciones, al lado de las líneas, corresponden a valores de isoresistividad. El área opacada al norte, en las secciones A y B, corresponde a una zona no sensible del modelo.
Figura 7
Secciones O-E del modelo resistivo 3D. Las profundidades se refieren a msnm. Los cortes A y B muestran las anomalías conductoras superficiales asociadas a los sitios de descarga de los manantiales termales. En el corte C se observa la ubicación del posible reservorio geotermal y de la capa sello. A lo largo de todos los cortes, se observa el basamento (R1) con una mayor presencia al SE del área. Los números dentro de las secciones, al lado de las líneas, corresponden a valores de isoresistividad. El área opacada al norte, en las secciones A y B, corresponde a una zona no sensible del modelo.


Alfaro et al. (2017) proponen un modelo conceptual para el área geotérmica a partir de estudios de geología, geoquímica y geofísica. La Figura 9 muestra una representación de dicho modelo conceptual sobre el modelo resistivo 3D obtenido en este trabajo. La recarga del sistema (flecha a en la Figura 9) sucedería en formaciones sedimentarias permeables (Tibasosa y Une) que afloran en el anticlinal de Tibasosa-Toledo al sureste del área, donde el terreno basculado al NNO permite la circulación de los fluidos en la mismadirección (Alfaro et al., 2017). La fuente de calor estaría asociada al calor remanente de las intrusiones del Altos los Volcanes y Alto Los Godos (Alfaro et al., 2017), cuya edad de evento más reciente es del orden de 1 millón de años (1,0±0,25 y 1,0±0,35 Ma) (Rueda- Gutiérrez, 2017). Los fluidos, dada la permeabilidad de las capas sedimentarias, así como la presencia de fallas de basamento (Paipa-Iza, Las Peñas), infiltrarían la zona intermedia entre los sectores dómicos (fuentes de calor) (Alfaro et al., 2017). Al encontrarse el sistema geotermal en un ambiente sedimentario es difícil identificar la respuesta resistiva de la capa sello y del reservorio geotermal (Van Leeuwen, 2016). En el modelo presentado (Figura 9), y bajo la suposición de la ubicación del reservorio entre los sectores dómicos (Alfaro et al., 2017), es posible observar de una forma más clara la presencia de una capa conductora somera (C3 en Figura 6 y 7) que se asociaría al sello del sistema y bajo esta se ubicaría el reservorio geotermal (Figura 7 y 8) con un aumento en la resistividad, como lo menciona Muñoz (2014).

Posteriormente, los fluidos termales continuarían su tránsito al norte por medio de la permeabilidad dispuesta por las fallas de Cerro Plateado y El Bizcocho hasta su descarga en los sectores de La Playa e ITP- Lanceros (flechas b y c en la Figura 9). Alfaro et al. (2017) establece que a lo largo de este tránsito, los fluidos recibirían aportaciones de aguas salobres más frías. Adicionalmente, proponen un segundo sistema de circulación, donde la fuente correspondería al calor radiogénico de la intrusión El Durazno (Figura 1) y que tendría como conducto el corredor entre las fallas Canocas y El Hornito (flecha d en Figura 9).

Análisis de sensibilidad del modelo resistivo obtenido. A. Ejemplos de modelos modificados para análisis de sensibilidad en profundidad i) modelo con resistividad de 1000 ohm.m a partir de la cota -1470 msnm, ii) modelo con resistividad de 1000 ohm.m a partir de la cota +1010 msnm al norte-oeste del modelo y -1470 msnm al centro-este y sur de la zona de estudio. B. Ejemplos de los ajustes de curvas del modelado directo del ejemplo A.ii donde los sondeos PA0020, PA0077 y PA0098 fueron sensibles a la cota -1470 msnm y el sondeo PA0073 fue sensible a la cota +1010 msnm. Las líneas punteadas corresponden a respuestas del modelo final obtenido, mientras que las líneas continuas corresponden a las respuestas de modelado directo realizado a los modelos modificados en profundidad. C. Secciones NS del modelo resistivo donde se observa i) la presencia de las anomalías conductoras C1, C2 y C3, ii) enmascaramiento de dichas anomalías con una resistividad de 30 ohm.m y iii) resultado de inversión del modelo enmascarado donde las anomalías conductoras se recuperan.
Figura 8
Análisis de sensibilidad del modelo resistivo obtenido. A. Ejemplos de modelos modificados para análisis de sensibilidad en profundidad i) modelo con resistividad de 1000 ohm.m a partir de la cota -1470 msnm, ii) modelo con resistividad de 1000 ohm.m a partir de la cota +1010 msnm al norte-oeste del modelo y -1470 msnm al centro-este y sur de la zona de estudio. B. Ejemplos de los ajustes de curvas del modelado directo del ejemplo A.ii donde los sondeos PA0020, PA0077 y PA0098 fueron sensibles a la cota -1470 msnm y el sondeo PA0073 fue sensible a la cota +1010 msnm. Las líneas punteadas corresponden a respuestas del modelo final obtenido, mientras que las líneas continuas corresponden a las respuestas de modelado directo realizado a los modelos modificados en profundidad. C. Secciones NS del modelo resistivo donde se observa i) la presencia de las anomalías conductoras C1, C2 y C3, ii) enmascaramiento de dichas anomalías con una resistividad de 30 ohm.m y iii) resultado de inversión del modelo enmascarado donde las anomalías conductoras se recuperan.


Bosquejo del modelo conceptual (Alfaro et al., 2017) sobre modelo resistivo 3D del área geotérmica de Paipa. Las flechas negras aluden a tránsitos de fluidos termales (ver texto). La zona opacada corresponde a una zona no sensible del modelo.
Figura 9
Bosquejo del modelo conceptual (Alfaro et al., 2017) sobre modelo resistivo 3D del área geotérmica de Paipa. Las flechas negras aluden a tránsitos de fluidos termales (ver texto). La zona opacada corresponde a una zona no sensible del modelo.


Conclusiones

El modelo resistivo del área geotérmica de Paipa obtenido a partir de la inversión 3D de datos magnetotelúricos, mostró la presencia de cuerpos conductores al centro-oeste del área, aparentemente controlados por las fallas de bajo ángulo (Canocas, El Hornito y El Bizcocho), formando conductos para el tránsito de fluidos hidrotermales hacia las zonas de descarga de La Playa e ITP-Lanceros. El basamento del área es identificado como un cuerpo resistivo al este y sur del área de estudio, sobre el que se identificó, en inmediaciones de los domos Alto Los Volcanes y Alto Los Godos, las respuestas eléctricas asociadas a la capa sello (cuerpo conductivo) y del reservorio geotermal.

Agradecimientos

El presente trabajo se realizó dentro de los objetivos misionales del Servicio Geológico Colombiano (SGC). El autor agradece a los integrantes del Grupo de Exploración de Recursos Geotérmicos del SGC, especialmente a la Química Claudia Alfaro y al MSc. Juan Camilo Herrera, por sus aportaciones y sugerencias durante la realización de este trabajo. Así como a los revisores, que con sus observaciones se logró una mejora sustancial del artículo.

Referencias

Abdul-Azeez, K.A.; Harinarayana, T. (2007). Magnetotelluric evidence of potential geothermal resource in Puga, Ladakh, NW Himalaya. Current Science, 93(3), 323-329

Alfaro, C. (2002a). Geoquímica del sistema geotérmico de Paipa. Informe Técnico, INGEOMINAS, Bogotá, Colombia

Alfaro, C.M. (2002b). Estudio isotópico de aguas del área geotérmica de Paipa. Informe técnico, INGEOMINAS, Bogotá

Alfaro, C.M.; Monsalve, M.L.; Franco, J.V.; Ortiz, I. (2012). Contribuciones al modelo conceptual sistema geotérmico de Paipa. Informe técnico, Servicio Geológico Colombiano, Bogotá

Alfaro, C.; Matiz, J.C.; Rueda, J.; Rodríguez, G.; González, C.; Beltrán, M.; Rodríguez, G.; Malo, J. (2017). Actualización del modelo conceptual del área geotérmica de Paipa. Informe Técnico, Servicio Geológico Colombiano, Bogotá, Colombia

Arango, C.; Marcuello, A.; Ledo, J.; Queralt, P. (2009). 3D magnetotelluric characterization of the geothermal anomaly in the Llucmajor aquifer system (Majorca, Spain). Journal of Applied Geophysics, 68(4), 479-488. https://doi.org/10.1016/j.jappgeo.2008.05.006

Bibby, H.G.; Caldwell, T.; Brown, C. (2005). Determinable and non-determinable parameters of galvanic distortion in magnetotellurics. Geophysical Journal International, 163(3), 915-930. https://doi.org/10.1111/j.1365-246X.2005.02779.x

Booker, J. (2014). The magnetotelluric phase tensor: a critical review. Surveys in Geophysics, 35(1), 7-40. https://doi.org/10.1007/s10712-013-9234-2

Caldwell, T.G.; Bibby, H.M.; Brown, C. (2004). The magnetotelluric phase tensor. Geophysical Journal International, 158(2), 457-469. https://doi.org/10.1111/j.1365-246X.2004.02281.x

Cepeda, H.; Pardo, N. (2004). Vulcanismo en Paipa. Informe Técnico, INGEOMINAS, Bogotá, Colombia

Chave, A.; Jones, A. (2012). The Magnetotelluric Method, Theory and Practice. Cambridge University Press. https://doi.org/10.1017/CBO9781139020138

Cumming, W.; Mackie, R. (2010). Resistivity imaging of geothermal resources using 1D, 2D and 3D MT inversion and TDEM static shift correction illustrated by a glass mountain case history. World Geothermal Congress 2010, Bali, Indonesia

Egbert, G.; Kelbert, A. (2012). Computational recipes for electromagnetic inverse problems. Geophysical Journal International, 189(1), 251-267. https://doi.org/10.1111/j.1365-246X.2011.05347.x

Gamble, T.; Goubau, W.M.; Clarke, J. (1979). Magnetotellurics with a remote magnetic reference. Geophysics, 44(1), 53-68. https://doi.org/10.1190/1.1440923

González, L.; Vásquez, L.; Muñoz, R.; Gomes, H.; Parrado, G.; Vargas S. (2008). Exploración de recursos energéticos exploración de uranio en Paipa, Iza, Pesca, Chivata (Boyacá). Informe técnico, INGEOMINAS, Bogotá

González-Idárraga, C.E.; Rodríguez-Rodríguez, G.F. (2017). Modelo Resistivo del área geotérmica de Paipa a partir de datos magnetotelúricos. Informe técnico, Servicio Geológico Colombiano, Bogotá

Kelbert, A.; Meqbel, N.; Egbert, G.D.; Tandon, K. (2014). ModEM: a modular system for inversion of electromagnetic geophysical data. Computer & Geosciences, 66, 40-53. https://doi.org/10.1016/j.cageo.2014.01.010

Kiyan, D.; Jones, A.G.; Vozar, J. (2014). The inability of magnetotelluric off-diagonal impedance tensor elements to sense oblique conductors in threedimensional inversion. Geophysical Journal International, 196(3), 1351-1364. https://doi.org/10.1093/gji/ggt470

Japan Consulting Institute. (1983). Feasibility study report of geothermal power plant for Instituto Colombiano de Energía Eléctrica (Electrificadora de Boyacá). Bogotá. 74 p

Llanos, E.M.; Bonet, C.; Zengerer, M. (2015). 3D geological - geophysical model building and forward and inverse modeling of magnetism and gravimetry data from Paipa geothermal area, Colombia. Informe técnico, Servicio Geológico Colombiano - Intrepid Geophysics, Bogotá

Moyano, I.E.; Vallejo, E.F. (2015). Contribution to the knowledge of the Paipa geothermal system by the application of the magnetotelluric method. World Geothermal Congress 2015. Melbourne, Australia

Muñoz, G. (2014). Exploring for geothermal resources with electromagnetic methods. Surveys in Geophysics, 35(1), 101-122. https://doi.org/10.1007/s10712-013-9236-0

OLADE; BID (1994). Guia para estudios de reconocimiento y prefactibilidad geotermico. Serie: Documentos de OLADE. Quito, Ecuador

OLADE; GEOTÉRMICA ITALIANA; ICEL; CONTECOL. (1982). Estudio de reconocimiento de los recursos geotérmicos de la república de Colombia. 7 volúmenes. Pisa: Organización Latinoamericana de Energía, Geotérmica Italiana Srl, Instituto Colombiano de Energía Electrica, Contecol ltda

Ortiz, I. (2014). Informe de revisión de la geología estructural en el área geotérmica de Paipa – Boyacá. Informe técnico, Servicio Geológico Colombiano, Bogotá

Ortiz, I.; Alfaro, C. (2010). Inventario de puntos de agua y geoquímica de las áreas geotérmicas de Paipa e Iza: aguas, suelos y peloides. Informe técnico, INGEOMINAS, Bogotá, Colombia

Pardo, N. (2004). Mapa geológico de Vulcanitas de Paipa, Planchas 171 y 191. Mapa escala 1:25.000. INGEOMINAS, Bogotá

Parkinson, W.D. (1959). Directions of rapid geomagnetic fluctuations. Geophysical Journal International, 2(1), 1-14. https://doi.org/10.1111/j.1365-246X.1959.tb05776.x

Peacock. J.R.; Thiel, S.; Reid. P.; Messellier. M.; Heinson, G. (2012). Monitoring enhanced geothermal fluids with magnetotellurics, test case: Paralana, South Australia. 37th Workshop on Geothermal Reservoir Engineering, Stanford, California

Rodríguez, G.; Alfaro, C. (2015). Caracterización de núcleos de perforación en las zonas de El Durazno, Paipa y criptodomo de Iza. Informe técnico, Servicio Geológico Colombiano, Bogotá

Romero, G.; La Terre, E.F.; Paneto, L.P.; Fontes, S.L. (2019). Upper crustal structures of the southeast edge of Parnaiba basin using 3D magnetotelluric data imaging. Journal of South American Earth Sciences, 96. https://doi.org/10.1016/j.jsames.2019.102392

Rueda-Gutiérrez, J. (2017). Cartografía de los cuerpos dómicos del área geotérmica de Paipa. Informe técnico, Servicio Geológico Colombiano, Bogotá, Colombia

SGC. (2016). Modelo hidrogeológico de Boyacá Centro. Informe técnico, Grupo de Exploración de Aguas Subterráneas, Servicio Geológico Colombiano, Bogotá

Simpson, F.; Bahr, K. (2005). Practical Magnetotellurics. Cambridge University Press. https://doi.org/10.1017/CBO9780511614095

Siripunvaraporn, W. (2012). Three-dimensional magnetotelluric inversion: an introductory guide for developers and users. Surveys in Geophysics, 33(1), 5-27. https://doi.org/10.1007/s10712-011-9122-6

Spichak, V.; Manzella, A. (2009). Electromagnetic sounding of geothermal zones. Journal of Applied Geophysics, 68(4), 459-478. https://doi.org/10.1016/j.jappgeo.2008.05.007

Ussher, G.; Harvey, C.; Johnstone, R.; Anderson, E. (2000). Understanding the resistivities observed in geothermal systems. World Geothermal Congress, Kyushu, Japan

Van Leeuwen, W.A. (2016). Geothermal exploration using the magnetotelluric method. PhD Thesis, Utrecht University, Netherlands

Vásquez, L. (2012). Aplicación geofisica de métodos potenciales en el área geotérmica Paipa-Iza. Informe técnico, Servicio Geológico Colombiano, Bogotá

Velandia, F. (2003). Cartografía geológica y estructural sector sur del municipio de Paipa. Informe técnico, INGEOMINAS. Bogotá, Colombia

Velandia, F.; Cepeda, H. (2004). Geología sector del sur del municipio de Paipa (Boyacá) planchas 171 y 191. Mapa escala 1:25.000. INGEOMINAS, Bogotá

Volpi, G.; Manzella, A.; Fiordelisi, A. (2003). Investigation of geothermal structures by magnetotellurics (MT): an example from the Mt. Amiata area, Italy. Geothermics, 32(2), 131-145. https://doi.org/10.1016/S0375-6505(03)00016-6

Vozoff, K. (1972). The magnetotelluric method in the exploration of sedimentary basins. Geophysics, 37(1), 98-141. https://doi.org/10.1190/1.1440255

Vozoff, K. (1991). The magnetotelluric method. In: M. Nabighian (ed.). Electromagneticmethods in applied geophysics (pp. 641-712). Volume 2, application, parts A and B. Society of Exploration Geophysicists. https://doi.org/10.1190/1.9781560802686.ch8

Información adicional

Forma de citar: González-Idárraga, C.E. (2020). Caracterización resistiva 3D del área geotérmica de Paipa, Colombia. Boletín de Geología, 42(3), 81-97 DOI: https://doi.org/10.18273/revbol.v42n3-2020003

Secciones
Boletín de Geología
ISSN: 0120-0283
Vol. 42
Num. 3
Año. 2020

Caracterización resistiva 3D del área geotérmica de Paipa, Colombia

Carlos Eduardo González-Idárraga
Servicio Geológico Colombiano,Colombia
Contexto
Descargar
Todas