Application of Time Domain Electromagnetics (TDEM) for Refining the Geometry of a Coastal Aquifer ()
1. Introduction
The vulnerability of aquifers in coastal cities is now well established. Groundwater, which exists in a very fragile balance between land and sea, is a highly sought-after freshwater resource in Africa for drinking water supply [1]. Among the major problems common to these cities in West Africa is the decline in the piezometry of the most exploited aquifers, which is closely related to rainfall and climate change [2] [3], and the intrusion of saline water lenses beneath freshwater aquifers [4]. In terms of quality, these waters are polluted by nitrates, and geochemical studies show that in some regions, redox reactions and dilution affect the quality of the water in these aquifers, with facies that are mostly sodium chloride [2]. These areas, which are particularly influenced by their proximity to the sea, attract populations for various economic activities [5], leading to intensive pumping and the depletion of aquifers. According to [6], approximately 3 billion people, or about half of the world’s population, live within 200 kilometers of a coastline, and these figures are believed to have doubled today. In Benin, the coastal sedimentary basin, which covers only 10% of the country’s total area, is home to 60% of the population [6], whose water needs must be met. Although the captured aquifer, the Miocene-Pliocene Continental Terminal (CT), is very productive, accounting for nearly 35% of the country’s groundwater reserves [4] [7], its monitoring remains essential for optimal and sustainable management of the resource. Given the importance of this aquifer, several studies in recent years have focused on improving knowledge of its geometry and piezometry [4] [7] [8], the physicochemical quality and dating of the water [7]-[9]; and the possible interactions between surface water and groundwater [10]. In 2022, the project to build a model for managing the Miocene-Pliocene aquifer on the scale of the Allada plateau (a plateau neighboring that of Sakété with almost identical hydrogeological and hydrodynamic conditions) with the AGIRES/INE projects [11], NOEVA projects, has highlighted some rather worrying projections regarding the influence of pumping cones on catchment fields between now and 2030 and 2050. Water flows induced by pumping will create water demand from Lake Nokoué, which is already subject to seasonal seawater intrusion. Our research aims to design a fairly accurate geometry of the CT aquifer at the scale of the Sakété Plateau in order to reproduce an exquisite hydrogeological operating model for the sustainable management of groundwater resources.
2. Study Area
2.1. Location of the Study Area
The Sakété Plateau is one of the three southern plateaus of Benin’s Coastal Sedimentary Basin (BSC). It is located in the southeast of the basin on the Nigerian border, bordering the coastal plain, and covers an area of 1,500 km2. According to Benin’s administrative divisions, it covers 11 municipalities and 57 districts in the departments of Ouémé and Plateau, with an estimated population of 1,096,850 in Ouémé and 624,146 inhabitants for the department of Plateau, representing a total of 17.20% of Benin’s population according to General Population and Housing Census (RGPH4) [12]. It is bordered to the north by a large depression called the Lama Depression, to the south by the coastal plain, to the west by a large valley of the Ouémé-So River complex, and to the east by Nigeria (Figure 1).
Figure 1. Map of the study area.
2.2. Hydrogeology and Hydrodynamics
Four aquifers can be identified in the coastal sedimentary basin. The oldest is the Upper Cretaceous aquifer, followed by the Eo-Paleocene aquifer, then the Miocene-Pliocene aquifer, and finally the Quaternary surface aquifer Figure 2. These different aquifer levels are interspersed with more or less continuous layers of marl and clay [4] [13]. The hydrogeological column of the basin identifies two powerful aquifers. The first and most productive is the Continental Terminal, which corresponds to hydrogeological units V to VII and whose water table is free on the Allada and Sakété plateaus with static levels of around 50 m, followed by the confined water table of the Upper Cretaceous (unit I), supplemented by units VIII and IIb (Holocene, and Lower Paleocene) [14]. The hydrogeological formations of the Continental Terminal, highlighted through various hydrogeological sections on the Sakété plateau, consist of sand, clay, clayey sand, clayey sand, and gravel with average permeability in the vicinity of the catchment area of approximately 8 × 10−4 m2·s−1 as horizontal transmissivity is between 4 × 10−3 and 14 × 10−3 m2·s−1 [8] [13]. Piezometric maps of the CT carried out on the Sakété plateau in 2004 and 2023 show a drop in level of approximately 10 m. The groundwater discharge outlets are notably the Porto-Novo Lagoon and the Ouémé River, with a main dome to the north of the Plateau [8] [10].
![]()
Figure 2. Hydrogeological cross-section of the BSC [14].
3. Materials and Methods
The materials and methodology used to construct the 3D model of the aquifer are presented in this section. First, a geophysical method (Time Domain Electromagnetics) was used with the TEM-FAST 48 HPC system (AEMR technology) for surveys following three profiles defined on the basis of the absence of borehole logs in certain areas. Coincident and central devices were set up at 32 sites and the raw data collected was inverted and modeled using the TEM-Res program. Then, Feflow software was used to design the 3D model of the Miocene-Pliocene aquifer.
3.1. Time Domain Electromagnetics (TDEM)
Three profiles were chosen for a series of TDEM surveys: a west-east profile on the northern edge of the plateau; a north-south profile on the western edge of the plateau; and a west-east profile to the south of the plateau. A total of 32 surveys were carried out across the three profiles, four of which were too noisy to be interpreted; the SNR (Signal to Noise Ratio) value is low. Three criteria guided the choice of these profiles.
a-Poor distribution and scarcity of borehole logs at the edges of the plateau; the spatial distribution of the few borehole logs that exist in this area is not well distributed in order to establish hydrogeological sections highlighting the variation in lithological facies.
b-Availability of previous data on the geometry of the aquifer; hydrogeological sections were carried out along sections A, B, C, and D by [8] [13] (Figure 3), and therefore fairly accurate information on the geometry of the aquifer exists along these sections.
c-The southern borehole profile was carried out in the vicinity of section C (Figure 3) in order to verify the consistency of the boreholes log data and the TDEM data [7].
Figure 3. Distribution profiles of TDEM surveys and existing cross-sections.
Our target for this study is the Miocene-Pliocene Continental Terminal aquifer. Previous studies at the BSC scale and at the Sakété and Allada plateau scales have shown that this aquifer has been extensively studied and has very strong potential. It is also the most important in terms of thickness, reaching up to 120 m [15] [16]. The Time Domain Electromagnetic method was preferred for its many advantages, which are well suited to our study context. Numerous studies list the advantages of TDEM, which include its high sensitivity to variations in the resistivity of conductive soils such as clays; the controlled nature of the EM source; the ability to vary the size of the transmission loop, allowing exploration of deep layers; and the very short deployment and measurement time, allowing several surveys to be conducted per day [17] [18]. In short, it provides good resolution of conductive terrain to depths of up to several hundred meters, depending on the loop size [18]-[20]. TDEM surveys are sensitive to electrically conductive layers up to several hundred meters deep. It has been tested by several authors in a sedimentary context with conclusive result for accurately determining the depth of the conductive clay aquiclude that generally marks the aquifer wall [21], then in conjunction with other geophysical methods [22] [23]. A coincident loop system was used with the “TEMFAST 48” system (AEMR technology), which allows the subsoil to be explored in the range of 0 m to 150 m with a low current intensity varying between 1 and 4 amps at a frequency of 50 Hz.
In the field, conducting a TDEM survey involves sending an electromagnetic pulse into the ground using a transmitter loop and studying the electromagnetic response of the ground to this pulse. This response, measured using a receiver loop, depends on the electrical conductivity of the ground at the survey location. These loops consist of one (or more) insulated electrical cable(s). During the survey, they are laid on the ground, but can also be airborne in certain applications. Each is connected to the TDEM measurement unit, which controls the entire data acquisition process. At each site, two surveys using square cable loops measuring 25 and 150 m on each side, with the same dimensions for the transmitter and receiver (“coincident” devices), and one survey using a 150 m transmitter and a 25 m receiver (“central” device) were carried out. The configuration of a survey is shown in Figure 4.
According to the TEM-FAST 48 HPC User Manual, Version 8a [24] it is a “simple, fast, and robust” system. Like its predecessors, the new instrument is designed for electromagnetic surveys of rock masses within the first 300 - 500 m based on the technology. Plates 1-3 show the device and some accessories used for a TDEM survey.
The parameter sought by a TDEM prospector, as with most geophysical methods, is resistivity. The presence or absence of water, and its nature, in a formation is reflected in the resistivity of that formation [17]. For TDEM, this parameter is not measured directly in the field. The initial parameter for TDEM prospecting is the TDEM signal (V/A) as a function of time (µs) of the vertical component of the secondary field. The iterative TEM-Res inversion software is used for the inversion and modeling of survey data. The raw data are inverted into resistivity as a function of time, and modeling is performed by varying the parameters (resistivities and thicknesses) as shown in Figure 5. For each survey at each site, the inversion of the TDEM data consisted of adjusting the experimental data with a theoretical model by minimizing the difference between the model and the data. These curves are then adjusted until they are as close as possible or superimposed. The quality of the inversion is assessed using this curve adjustment, which is closely linked to the Root Mean Squared (RMS) criterion. All inversions are validated at Root Mean Squared values below 4%, which is an acceptable level of accuracy in inversion models [23].
![]()
Figure 4. Configuration of a TDEM loop for a coincident survey with a 150 m loop.
Plate 1. TDEM equipment and accessories.
Plate 2. Cables.
Plate 3. Compass.
Figure 5. Diagram showing the natural boundaries of the Miocene-Pliocene aquifer on the Sakété plateau.
Archi’s law [24] on saturation allows us to take into account the influence of certain parameters such as porosity, clay content, the percentage of voids filled with water and its ion concentration, temperature, and the state of the water [17]. Although the ranges of resistivity values for formations have been broadly defined by several authors [25]-[27], it is more practical to contextualize them using a formula that takes into account some of these parameters that influence the resistivity of facies. Archie’s law is an experimental relationship used in geophysics that links the electrical resistivity of a saturated rock to its porosity and the resistivity of the water that impregnates it. The hydrogeology of Miocene-Pliocene aquifer as consisting of continental deposits composed of sand, sandy clay, clay, and more detrital and coarse sandstone [8] [15] [28] [29]. For these unconsolidated formations, the reformulated Archi Equations (1)-(3) are used to determine the resistivity of the formations.
(1)
(2)
(3)
is the resistivity of the rock,
is resistivity of the inhibiting water contained in the pores,
a cementation factor with a value of 1.37, and
is the porosity of the rock, estimated between 18 and 35% for moderately cemented sedimentary rocks and between 40 and 50% for clays.
is the saturation coefficient and its value is 0.88 according to Keller’s 1988 work cited by [17]. These data were used in various formulas with the aquifer’s water conductivity to calibrate the resistivity of each facies crossed by the boreholes. The electrical conductivity (EC) of the water in the Mio-Pliocene aquifer varies from 24 µS/cm to 920 µS/cm during the rainy season, with an average conductivity of 333.25 µS/cm. During the dry season, conductivity is slightly higher (27 µS/cm to 1205 µS/cm), with an average of 320.46 µS/cm [10] [30]. Thus, the deduced resistivity ranges are as follows: between 100 and 600 for bar soil; between 2 and 50 Ohm.m for clays; and between 700 and 2000 Ohm.m for layers of sand, fine to coarse sand.
3.2. 3D Digital, Model of the Study Area
The geographical boundaries of the study area presented above are not the boundaries of the Miocene-Pliocene aquifer under study. This aquifer is shared by Togo on the western side of Benin, as well as by Nigeria to the east. In order to better define the boundary conditions of the aquifer for flow simulations, the natural boundaries are shown in Figure 5. The simulation part is not presented in this article; it only presents the constitution of the aquifer geometry and therefore no flow conditions (upper recharge, lateral flows, river loads) are yet imposed as well as the sensitivity study of the inversion of surveys on the aquifer walls, which is more closely related to the modeling software used.
FEFLOW (Finite Element subsurface FLOW and transport system) is an interactive modeling system with a finite element mesh. Numerical models in hydrogeological modeling are often associated with a discretization of an aquifer system into geometric elements of various shapes described by differential equations [18]. One of the major problems encountered by hydrogeological modelers on fairly complex aquifers is the identification of the aquifer wall [20]. The reliability of the numerical model and the quality of the outputs will depend on a good reconstruction of its system geometry. Drilling logs have always been one of the data from which it is possible to reconstruct the chronology of subsoil formations and to extend to a hydrogeological section to see the geometric structure of a given area. For more or less complex contexts, it would be more suitable for finite element discretization [31]. The finite element discretization equation in Feflow for a coastal aquifer can be written in the form of Equations (4)-(7) in a saturated porous medium:
[32] [33](4)
where
Ss is the mass matrix;
K is the stiffness matrix;
H is the hydraulic load vector;
Q is the vector of sinks and sources.
Galerkin method is a numerical technique used to solve partial differential equations in finite elements. In the context of Feflow this flow equation can be written as:
[32](5)
Finite element discretization consists of dividing the domain ω into finite elements, denoted ωe. On each element, we define local basis functions, denoted
.
The approximate solution h is then expressed as a linear combination of the local basis functions:
[34](6)
Galerkin formula with source term can be written in the form:
[35] [36] (7)
Q represents the source term at the boundary
of the elements
.
The discretized super mesh is characterized by 272,892 prism-triangle elements, 193,864 and 6 nodes per element. To account for the heterogeneities of the layers in the model across the different interfaces defined by the drilling sections, the depths of the interfaces are related to the altitudes of each point. For the mesh, at pumping wells; especially the boreholes of National Water Company of Benin (SONEB), the hydraulic gradients are expected at the center of the well cone. To better represent them, a locally fine discretization is applied in Figure 6.
Figure 6. Numerical finite element model of the Sakété Plateau.
For this study, 527 boreholes logs were used and three layers were reconstructed from the lithological descriptions of these logs based on TDEM survey results and previous hydrogeological knowledge of the aquifer geometry at the BSC Plateau scale [7] [8] [30] [37] [38]. The interfaces of these three layers are defined as follows: a first interface Z1 corresponding to the Digital Terrain Model (the topography of the natural terrain); an interface Z2 corresponding to the Upper Mio-Pliocene bar-land interface; an interface Z3, which is the Upper Mio-Pliocene-Lower Mio-Pliocene interface; and Z4, which corresponds to the depth of the substratum or base of the Lower Mio-Pliocene. Interfaces Z2, Z3, and Z4 are marked by clayey facies and are therefore impermeable layers. The regionalization of these elevation data with the geographical coordinates of the drilling logs allows the raw geometry of the model to take its actual form. Z1 is extracted from a CGIAR-CSI v4.1 Digital Terrain Model with a spatial resolution of 90 m × 90 m downloaded from https://srtm.csi.cgiar.org/. These elevation values can be corrected using Equation (5),
[38](5)
where a and b are factors in the regression line equation Figure 7 between the elevation values of the CGIAR-CSI v4.1 DTM (retained) and those of the geodetic points of the National Geographic Institute of Benin.
[37](6)
4. Results and Discussions
The results will be presented in two sections. The first section presents the results obtained from the interpretation of TDEM surveys, and the second section presents the results of the construction of the 3D geometry of the aquifer in the Feflow numerical code.
4.1. Results
TDEM results
Figure 8 shows the resistivity curves as a function of time. Measurement points with large error margins were removed from the curves in order to better model the surveys.
Figure 7. Regression line of the elevation values of the CGIAR-CSI v4.1 DTM and those of the geodetic points collected at the IGN in the study area.
Figure 8. Geotechnical models of TDEM boreholes on the North profile.
On the North profile, the survey curves (Resistivity vs. depth) in Figure 9 show an initial resistivity formation varying from approximately 90 to 400 Ohm.m and with a thickness between 10 and 20 Ohm.m, with the exception of survey S1, where all the formations crossed are conductive (Ro < 10 Ohm.m). The second layer is more resistant than the first (1000 and 2000 Ohm.m) with a thickness between 20 and 50 m. The formations below, between 40 and 100 m deep, are more conductive (Ro < 10 > 1 Ohm). The conductive wall of the Miocene-Pliocene aquifer is therefore estimated to be between 50 and 100 m deep on the northern edge of the Plateau near the Lama depression.
Figure 10 shows the distribution of resistivity survey curves as a function of depth along the profile. These curves on this profile, as well as on the other profiles, are presented with their actual field spacing. The northern profile covers a distance of 48,000 m. The spatial distribution of resistivity values on this vertical section Figure 11 provides a visual representation of the variation in resistivity along the profile. A conductive layer with a resistivity between (Ro < 10 > 1 Ohm)
Figure 9. Resistivity curves of the terrain with a depth.
Figure 10. Spatialization of boreholes on the northern profile at DTM elevation for each site.
can be observed at (S1) in Igbodogui in the municipality of Pobè. This thick impermeable layer is also very close, at a depth of about 10 m in the municipality of Adja-Ouèrè (S8). A deeper aquifer formation is detected by boreholes S6 and S8, which could be an Upper Cretaceous aquifer.
On the west profile, the geotechnical models of the boreholes are shown in Figure 12. A first layer of resistivities varying between 200 Ohm.m and 500 Ohm.m can be observed, as shown in Figure 13. The second layers have resistivities of 1000 and 2000 Ohm.m. Conductive layers (Ro < 10 > 1 Ohm.m) are observed between depths of 90 and 120 m depending on the site; this is the esti mated depth of the aquifer wall on the western edge of the plateau. Figure 14 shows the spatial distribution of the boreholes according to the actual topography and borehole depths along the entire profile.
Spatialization of the layers with their resistivity values Figure 15 highlights the presence of conductive layer lenses (Ro <10 >1 Ohm.m) in the center of the profile towards the southern part of the western profile.
Figure 16(a) and Figure 16(a) show the geotechnical models of the boreholes and the soil resistivity curves as a function of depth on the western profile, respectively.
Figure 11. Spatialization of soil resistivity on the northern profile at DTM elevation for each site.
Figure 12. Geotechnical models of TDEM surveys on the western profile.
Figure 13. Resistivity curves of the terrain as a function of depth on the western profile.
Figure 14. Spatial distribution of surveys on the western profile at DTM elevation for each site.
Figure 15. Spatial distribution of soil resistivity on the West profile at DTM elevation for each.
The spatialization of resistivities in Figure 17 and Figure 18 provides an overall view of the contrast in resistivity between the layers. The conductive layer (Ro < 30 m) begins at an average depth of 80 m, except in the S26 zone common to Ifangni and Daagbé, where this impermeable layer is much deeper (110 m). On this profile, the depth of the aquifer wall is estimated at 110 m to the east of the profile and remains almost constant at 80 m towards the east.
Figure 19 shows in a, b, c, and d the interfaces reconstructed by regionalization of elevations with the coordinates of the borehole logs and surveys. The natural heterogeneities of the aquifer are evident at each of the interfaces.
Figure 16. (a) Geotechnical models of TDEM boreholes on the South profile; (b) Soil resistivity curves as a function of depth on the West profile.
Figure 17. Spatial distribution of boreholes on the southern profile at DTM elevation for each site.
Figure 18. Spatialization of soil resistivity on the southern profile at DTM elevation for each site.
Figure 19. Schematic diagrams of slices a, b, c, and d.
The superimposition of the four slices (3D geometry) shown in Figure 20 reveals that the different layers of the aquifer are thicker in the southeast of the plateau and become thinner towards the Ouémé River, the Lama Depression and Nigeria.
4.2. Discussions
In general, TDEM surveys highlight a contrast in the resistivity of conductive terrain between depths of 50 m and 100 m on the northern edge of the Plateau near the Lama depression; these are the estimated depths of the Miocene-Pliocene aquifer wall in this area. It is much shallower in places (around 20 m) and outcrops towards Nigeria on this profile. On the western profile, it is located between 90 and 120 m depending on the site. On the southern profile, the depth of the aquifer wall is estimated at 110 m to the east of the profile and remains almost constant at 80 m towards the east. The southern profile is defined in the same transect as the cross section DD in Figure 21 [28]. These surveys highlight the local heterogeneity of the aquifer. Analysis of the resistivity of the boreholes drilled on this profile allows us to deduce, using Archie’s law [24], that the formations below the barren soil consist of sandy clay of medium to coarse sand and sandstone sutured with water in places and limited by a dis-continuous clay layer of variable thickness. The depth of the conductive layer detected by TDEM is in the same range as the depth of the clay layer highlighted on this section, confirming the reliability of the TDEM data in this hydrogeological context of the study area. According to the work of [7] [13], the first layer corresponds to more or less permeable bar soil allowing the recharge of the underlying aquifer consisting of sand with clay lenses in places. According to [15] [28], on the southern plateaus of the BSC, the deposits are continental and consist of sand, sandy clay, clay, and more detrital and coarse sandstone. The layers are discontinuous and irregular, with lenticular clays. TDEM surveys confirm this discontinuity of the impermeable clay layer, the thickness of which has not been estimated.
![]()
Figure 20. 3D model of the CT aquifer on the Sakété Plateau.
Figure 21. Hydrogeological cross-section showing the hydrogeological formations and their thicknesses [7].
In addition, two surveys (S27 and S28) were also carried out near the boreholes to verify the consistency of the surveys and logs Figure 22. At survey S27, a layer of clay approximately 5 m thick, highlighted by the borehole log, was not detected by the TDEM. This layer is therefore not laterally continuous. For borehole S28, located 530 m from the log, the clay layer (Ro = 3 Ohm.m) with a thickness of 3 m, detected at 65 m by the borehole, is visible from 62 m by the drilling log. A joint analysis of these two data sets suggests that the first level corresponding to the bar soil layer has a resistivity ranging from 300 to 600 Ohm.m and an average thickness of 35 m, consisting of fine to coarse sand 35 m thick at this borehole. The thickness of the layers appears to be underestimated by the TDEM surveys.
![]()
Figure 22. Resistivity curves of the boreholes as a function of depth and lithology section of two boreholes (S27 and S28) close to these boreholes. (a) Survey S28 located 530 meters from the borehole; (b) Survey S27 located 769 meters from the borehole.
5. Conclusion
Hydrogeological modeling is based not only on quantitative and qualitative input data, but also on a fairly detailed and accurate reconstruction of the reservoir geometry. In order to model the underground flow of the Miocene-Pliocene aquifer on the Sakété Plateau, the first step, which is the construction of the 3D geometry of the aquifer, was developed in this study. Four interfaces and three layers were identified, revealing the natural heterogeneity of the aquifer. Lateral thinning was also identified at the edges of the Plateau, which constitute the discharge zones of underground flows. TDEM surveys conducted at the edges of the plateau proved to be an effective tool for refining the geometry. These TDEM surveys highlight a contrast in the resistivity of the conductive terrain between 50 and 100 m deep at the northern edge of the plateau, near the Lama depression. These are the estimated depths of the Miocene-Pliocene aquifer wall. It is much shallower in some places (around 20 m) and outcrops towards Nigeria on this profile. At the western end, the aquifer wall is between 90 and 120 m deep, depending on the site. On the southern profile, the depth of the wall is estimated at 110 m east of the plateau and remains almost constant at 80 m towards the east. These surveys highlight the local heterogeneity of the aquifer. The 3D model highlights low-lying areas with Quaternary formations in the lower southern part of the plateau and in low-lying areas, while the central and upper parts feature Miocene-Pliocene formations. For the rest of our work, boundary conditions will be imposed on the 3D model and simulations can be carried out to understand the behavior of the aquifer in order to plan sustainable management of groundwater resources on the Plateau.
Acknowledgements
The authors thank the African Center of Excellence for Water and Sanitation in Benin.