Introduction

The Haiti alluvial plains system has been recognized for a long time as the country’s greatest potential for directly exploitable groundwater resources (Adamson et al. 2016; Margat 1990; Martin and Moussu 1977). Stretching over the departments of North and North-East, the alluvial system of the Plaine du Nord/Massacre transboundary alluvial aquifer, hereby called Plaine du Nord aquifer for brevity, is considered a major aquifer (UNDP 1991), extending from the plain of the Limbé River in the west to the Massacre River, whose transboundary watershed extends into the Dominican Republic (ISARM and UNESCO 2009). The alluvial aquifer system is used to supply water to smaller towns such as Fort Liberté and Terrier Rouge, as well as larger communities such as Cap-Haitien, the second largest city of the country where an extension and rehabilitation program of the water supply network is ongoing. Hundreds of hand-pumps and thousands of domestic wells supply all other smaller communities from the upper part of the aquifer (Gutierrez 2011). The population of Cap-Haitien (Arrondissement) is estimated to be 357,000 inhabitants (IHSI 2015) and could reach nearly 600,000 by 2025 (Archambault and Biscan 2018). The Plaine du Nord aquifer is also the main water source for industrial sites, power plants, and agricultural irrigation. Indeed, the alluvial aquifer plays a major role in the supply of water for the whole north and north east of Haiti.

In addition to growth in demand for water, climate change in the northern Caribbean is expected to result in intensified droughts due to rising temperatures and decreasing annual rainfall (Hall et al. 2013). A rise in sea level is also expected in the region (Le Cozannet et al. 2020), which will exacerbate saltwater intrusion in coastal aquifers (IPCC 2007). Adamson et al. (2016) and Gutierrez et al. (2012) have already documented saltwater intrusion in the Plaine du Nord aquifer. Collectively, these threaten water security for all users.

The Appui à la Valorisation du Potentiel Agricole du Nord, pour la Sécurité Économique et Environnementale (AVANSE) [Support for the Valuation of the Agricultural Potential of the North, for Economic and Environmental Security] project began in 2013, funded by USAID and implemented by Development Alternatives Inc. (DAI). The principal aim of the project is to promote and develop a sustainable approach to agriculture within the region (DAI Global, LLC 2020). One of the objectives of the project is to increase irrigated area from 1,500 to 4,000 ha in a sustainable manner (Molnar et al. 2015). This expansion of irrigated land will assist in securing the production of the value chains of banana, rice, and cacao, which are vital to the regional and national economy. To assist in the sustainable management of water resources, a groundwater model of the shallow water in the Plaine du Nord was developed to support the decision-making process. A thorough understanding of the water balance in the area, as well as the many different interactions between water sources and extractions was necessary for the conceptual framework of the model. To develop a robust decision tool, it was necessary to search for and compile available data, acquire relevant key additional data, develop a conceptual and numerical model, calibrate the model in the first instance at a steady-state, and then develop some prospective scenario simulations related to the extension of irrigated areas for banana crops.

The lack of functional climatic and water monitoring networks (Adamson et al. 2016) was the most challenging hurdle to overcome. A local team of researchers mobilised the various water stakeholders in the study area to make the best use of historical data and current status of knowledge. Critical to the success of this study, and indeed to the furthering of hydrogeological knowledge in northern Haiti was the close collaboration between DAI, BRGM (French Geological Survey), and the Universities of Port-au-Prince and Limonade. Each milestone of the modelling process was supported by a training session in Port-au-Prince, bringing together BRGM hydrogeologists, current graduate students, professors, and DAI technicians and engineers.

Materials and methods

Study area

The study area stretches across the North and North-East departments of Haiti. The northern alluvial plain concerns five surface watersheds, from west to east: Haut du Cap, Grande Rivière du Nord, Trou du Nord, Marion, and Jassa (Fig. 1). The watershed of the Jassa River is a tributary of the Massacre River which has a watershed that extends beyond the border with the Dominican Republic. The average annual rainfall varies between 1,000 and 1,300 mm, with a winter dry season and a short dry season in July–August. The average annual temperature varies between 22 and 27 °C.

Fig. 1
figure 1

Location map of the study area (SRTM90m produced by NASA)

Full size image

Geological and hydrogeological setting

The northern plain corresponds to a filling of the collapsed northeast slope of the northern massif by detrital materials. The geology of the region is affected by the proximity of the northern edge of the Caribbean plate in collision-subduction with the North American plate. The great northern fault (sinister strike-slip in the direction of WNW–ESE) which affects the north of the Dominican Republic passes offshore of the northern plain which happens to be in a tectonic strip, limited in the south by the fault of la Española (Monthel et al. 2013). The nature of the deposited alluvium is very heterogeneous as evidenced by the high variability of the available geological logs (Tomasi 1990a). Alluvium is made up of alternating sands, gravels and clays. Silty facies over at least 6 m are found on the surface (Monthel et al. 2013). The coarsest facies were found in the old beds of Grande Rivière du Nord and Trou du Nord. In contrast, clay facies seem to predominate in Phaëton, Fort-Liberté, Terrier Rouge and Carrefour Oranger. The alluvium that make up the northern plain is thick, increasing from south to north, up to more than 300 m at the edge of the ocean in the central zone, a thickness estimated by geophysics (Radstake and Chery 1992; Miner et al, Northwater International (USA), ‘Geophysical reconnaissance of the western basin of the Plaine du Nord aquifer, Haiti’ (Hydrogeology Journal, this issue), 2022). Boreholes reaching the basement are rare. Among the UNDP drilling campaign results in 1990, only two drillings reached the basement (Phaeton and Chabert, respectively at depths of 63 and 37 m). The local presence of the basement at shallower depth reflects the existence of thresholds between the Trou du Nord, Terrier Rouge and Fort-Liberté basins, which could represent an obstacle to underground flows (Tomasi 1990b).

In the area of Cap-Haitien, the Balan municipal boreholes that serve Cap-Haitien intersected over a 100 m of alternations of sand and clay which correlate quite poorly despite the proximity of the drillings (a few hundreds of meters) but in which coarser layers (3–5 m thick) can be identified around the depths of 10, 30, 45 m and sometimes at deeper levels but at varying depths depending on the boreholes. The inadequacy of the position of the well screen and the filter beds led to the silting up of old boreholes and the construction of two new ones (F10 and F11) in 2005 (BRGM-BRL-PC Consortium 2011). In 2020, the rehabilitation of two of the four boreholes and the drilling of two new ones in the same area was planned (Archambault and Biscan 2018; Adamson and Miner 2018). These sandy-gravel levels are framed by silty and clay formations which confine them. Some of the deeper boreholes, located in the downstream part, are even artesian (La Rue, Carrefour Oranger, Quartier Morin); however, the water body is never far from the topographic surface. Between Cap Haitien and La Grande Rivière du Nord, an inventory undertaken by BRGM in 2011 shows a first sandy-silty level tapped by hundreds of shallow domestic wells, 4 m deep on average (Gutierrez 2011).

Borehole production as documented by Tomasi (1990b) are on the order of 35–70 m3/h, reaching higher values in the old river beds (Trou du Nord, Grande Rivière du Nord). Pumping rates of the Balan boreholes in 2018 were on the order of 100–180 m3/h—F8: 100–140 m3/h, F9: 180 m3/h, F10: 112 m3/h, F1: 126 m3/h (Northwater International, unpublished data, 2018), slightly lower than in 2010 (BRGM-BRL-PC 2011). Borehole F8 is now abandoned due to a break in the casing (J. Adamson, Northwater International, personal communication, 2020). Well testing performed on the Balan boreholes by Inter-American Development Bank in 2018 recorded discharge at 284 m3/h with 43 m of drawdown at F9, and 205 m3/h with 36.5 m of drawdown at F8 (Adamson and Miner 2018).

To the south of the plain, the Massif du Nord is made up of intrusive volcanic formations (rhyolites, dacites, andesites, diorites, gabbros, etc.), rocks that are by nature not very permeable (secondary porosity). Many springs are exploited from these rocks, with average flow rates between 0.18 and 3.6 m3/h (Tomasi 1990c). Large-diameter wells are also dug in the weathered shallow horizon. The productivity of boreholes intercepting the underlying fissured horizon remains limited (maximum 3.6 m3/h). The hypothesis of a largely dominant runoff on the northern slope of the massif appears realistic. This runoff contributes to the recharge of the groundwater which flows from the south to the north with an average gradient of 0.45%, with a stronger gradient upstream an\d near the border of the massifs (up to 2%) and weaker gradient downstream in the flat areas. The reported piezometry is shown on the hydrogeological map in Fig. 2 (MARNDR and UNDP 1990).

Fig. 2
figure 2

Extract from the hydrogeological map of Haiti at the scale of 1:250,000 (MARNDR and UNDP 1990)

Full size image

Hydrochemical characterization of groundwater carried in the Cap Haitien region (Gutierrez et al. 2012) confirms that the alluvial plain is a multilayer aquifer, the most superficial level of which is in oxidizing conditions (unconfined water table). On the other hand, the semi-deep hand-pump wells (15–20 m deep in this sector) tap deeper levels that are confined and protected by clay or silty layers. The semi-deep (15–20 m) and deep (up to 100 m deep) wells benefit from good natural protection and denitrification under reducing conditions.

The chemical analyses and piezometric measurements carried out in October 2010 by BRGM on the sector of the Grande Rivière du Nord that crosses the plain on a north–south axis show without ambiguity that the river strongly drains the aquifer upstream of Pont Parois (Gutierrez et al. 2012). The piezometric map in Fig. 2 also shows the same flow trend for all other rivers of the study area. On the other hand, on the upstream end of the plain, the valley bottoms are topographically higher than the water table; the flow of watercourses would thus participate in recharging the groundwater by infiltration (Tomasi 1990c).

Development of the groundwater flow model

Computer code

The finite-difference based MARTHE software was used for the model (Arnaud and Amraoui 2018; Kloppmann et al. 2012; Thiéry 1993, 2010a, b, c; Thiéry et al. 2018, Thiéry and Picot-Colbeaux 2020), thus enabling three-dimensional (3D) simulations of groundwater flow, mass and energy transport in porous media. The software allows one to model hydrosystems by integrating groundwater flow, hydro-climatic budgets and surface-water flows. The underlying schemes can be simple or complex (unsaturated flow, multiphase flow, density-dependent fluid flow, consideration of vegetation, aquifer–river interactions, etc.). It is also possible to make a reactive transport calculation using a coupling with the PHREEQC code (Thiéry 2015).

The software uses a finite-difference numerical scheme using irregular parallelepiped or rectangular stencils. This can represent complex multi-layer aquifer systems, with pinchouts and the occurrence of free surfaces in any layer. MARTHE addresses modelling of regional multi-layer systems covering several thousand square kilometres, as well as local-scale modelling of porous media systems down to a few cubic centimetres, with millimetre-size cells. It is especially designed to evaluate and manage groundwater resources. Finally, the MARTHE software presents all the functions necessary for modelling of the alluvial groundwater of the North Plain.

In the MARTHE code, there is no concept of “confined layer”, “unconfined layer” or “convertible layer”. At any time, any cell can become confined, unconfined or dewatered, depending upon the current hydraulic head. The transmissivity is calculated from the hydraulic conductivity assigned to each mesh and the water height. In the transient mode, MARTHE uses two storage coefficients: the storage coefficient in an unconfined state (equivalent to the effective porosity) when the cell is not fully saturated with water, and the confined storage coefficient when the cell is saturated with water.

The MARTHE software includes a package that fully couples the groundwater flow calculation with the river network, as well as lakes and wetlands. Depending upon the relative position of the water level in the river and of the aquifer’s groundwater level, water flows can be simulated from the river reach to the aquifer, or, inversely, from the aquifer to the river network. Such flow exchanges are controlled, in a nonlinear manner, by linking riverbed and banks proportionally to their surface area, thickness and hydraulic conductivity. Different configurations of exchanges between an aquifer and a river can be simulated by MARTHE as briefly presented in Fig. 3.

Fig. 3
figure 3

Different configurations of exchanges between an aquifer and a river with MARTHE and related Eqs. (1), (2) and (3) (see text)

Full size image

Grid design and geometry

The model domain overlays a regular grid of 136 rows and 343 columns with 200-m2 cells. In total, the model is made up 20,431 active cells, covering an area of 817 km2.

Given the very high heterogeneity of the facies involved (clay, sand, gravel) and the rapid lateral variations of facies (lenticular geometries), it was impossible to produce a 3D geological model of the alluvium due to the lack of localised geological information. Therefore, after performing a sensitivity analysis, it was decided that the alluvial aquifer would be equally well represented by a single-layer geological model. This analysis consisted of a preliminary simulation of different aquifer structures—one layer, two layers, and three layers with two aquifers and one aquitard.

To the north, the conceptual model is limited by the natural limits of the alluvial outcrops and the coastline, while towards the south, the model extends to the foot of the slopes of the Massif du Nord, and towards the east, the model was extended over the administrative border with the Dominican Republic in order to include the entire watershed of the Jassa and Massacre rivers.

A 5-m resolution digital elevation model derived from LiDAR was used to represent the topographic surface; it was provided by CNIGS (Centre National de l’Information Géo-Spatiale). Within the modelled area, elevations vary between 0.15 and 72 m, with a mean value of 19 m above sea level. The geophysical surveys carried out by BRGM in 1980 (588 vertical electrical soundings), coupled with the drilling data, allowed the estimation of a depth of the alluvial bedrock (Radstake and Chery 1992). The thickness of the alluvium varies between 36 and 300 m, with a mean value of 141 m (Fig. 4).

Fig. 4
figure 4

Map of estimated thickness of the alluvial aquifer of the Plaine du Nord

Full size image

Hydraulic boundary conditions

The boundary conditions were defined according to the hydrogeological configuration of the Plaine du Nord as follows:

  • The groundwater level was fixed to zero along the coastline to represent the sea boundary condition (389 cells concerned).

  • The groundwater level was fixed along the natural limits of alluvial deposits and limestone or along the limits of volcanic formations to represent the inflow from the latter (599 cells concerned). The imposed hydraulic heads are first reconstructed from the linear relationship between the groundwater levels and the topography, then adjusted during the calibration process. The linear regression results from observations made by DAI on 792 wells (piezometric level = 0.919 × elevation – 1.615, R2 = 0.97). In the transient state, the imposed hydraulic heads are replaced by imposed flows.

  • No flux (Neumann-type boundary conditions) was fixed in the lateral boundaries to represent groundwater divide lines.

  • For the five main rivers, the exchange between groundwater and surface water was represented by river limits according to Eqs. (1), (2) and (3) (Fig. 3).

$$ {Q}_{\mathrm{Percol}}={\mathrm{AREA}}_{\mathrm{Exch}}\times {K}_{\mathrm{R}}\times \frac{\left({H}_{\mathrm{R}}-{H}_{\mathrm{b}}\right)}{\mathrm{Thickn}} $$
(1)
$$ {Q}_{\mathrm{Exch}}={\mathrm{AREA}}_{\mathrm{Exch}}\times {K}_{\mathrm{R}}\times \frac{\left({H}_{\mathrm{R}}-{H}_{\mathrm{G}}\right)}{\mathrm{Thickn}} $$
(2)
$$ {Q}_{\mathrm{Overfl}}=\left({\mathrm{AREA}}_{\mathrm{cell}}-L\times W\right)\times {K}_{\mathrm{R}}\times \frac{\left({H}_{\mathrm{G}}-{H}_{\mathrm{Topo}}\right)}{\mathrm{Thickn}} $$
(3)

where the terms are defined as follows:

L:

Length of the river reach

W:

Width of the river reach

KR:

Hydraulic conductivity of the river bed and river banks

Thickn:

Thickness of the river bed (and the river banks)

HR:

Absolute elevation of the free water surface in the river

HG:

Hydraulic head in the underlying groundwater

Hb:

Elevation of the bottom of the river bed

HTopo:

Elevation of the topographic surface of the cell containing the reach

AREACell:

Horizontal area of the cell in which the reach is located

AREAExch:

Area of the groundwater-river exchange =L x [W + 2 x (HR − Hb)]

QExch:

Flow of the river exchange towards groundwater

QPercol:

Flow of the river percolation towards disconnected groundwater

QOverfl:

Complementary overflow of groundwater to the river

Groundwater recharge

Due to the lack of recent and homogenous climatic data, it was not possible to include the calculation of hydro-climatic balances within the model to establish the groundwater recharge potential within the study area. It was decided to rely on the national recharge map established by Miner and Adamson (2017). This geographic information system (GIS) model is achieved with a water balance calculation, and applies a coefficient based on country-wide GIS layers that include soil and geologic permeability, vegetation cover, relief, and karstification.

Based on the work of Adamson et al. (2016), the modelled domain was finally divided into five homogeneous recharge zones (Fig. 5), in which each mesh cell receives an identical recharge flow. The average value, calculated from the grid at 30-m intervals, was initially assigned to each of the zones, i.e. 100 mm/year for zones 1 and 2, 170 mm/year for zone 3, 200 mm/year for zone 4, and 250 mm/year for zone 5. These values were subsequently adjusted during the calibration of the model.

Fig. 5
figure 5

Homogeneous recharge zones defined in the Plaine du Nord model and initial values for the steady sate conditions

Full size image

Pumping wells

In terms of current use of the wells for groundwater abstraction (Fig. 6), groundwater uses were grouped in three categories: drinking water supply, irrigation use, and industrial uses (such as Caracol industrial park, Larue distillery, ice factories, water purifiers, etc.). Drinking water is supplied through a few deep boreholes and a large number of semi-deep boreholes (hand pumps). For the deep boreholes supplying the municipalities of Cap Haitien and Ouanaminthe, the abstracted volumes were collected from DINEPA (DIrection Nationale de l’Eau Potable et de l’Assainissement).

Fig. 6
figure 6

Map of average annual groundwater withdrawals in the Plaine du Nord

Full size image

For hand-pump wells, the estimation of the quantities required one to cross check three main sources of information to arrive at a reliable average estimation of pumped water volumes. These three sources are: the current distribution of population by municipality, the average water consumption per capita, and the pumped water volumes extracted from the wells identified under the ‘Haiti Outreach’ database (Haiti Outreach 2018).

For the average water consumption per capita, a distinction was made between rural areas (20 L/day per inhabitant) and urban areas (70 L/day per inhabitant); these estimations were extracted from a DINEPA technical leaflet (DINEPA 2013). Each well listed as functional in the Haiti Outreach database was assigned an average annual withdrawal volume—an example of the calculation is provided for the municipality of Fort-Liberté in Table 1.

Table 1 Method used to calculate average annual withdrawals for drinking water supply—example of Fort-Liberté
Full size table

The same approach to estimate the withdrawals from irrigation wells proved to be difficult due to the lack of regular monitoring and reliable data. The decision was made to ignore their influence for the time being, with the exception of two working wells listed in a BRGM study performed in 2010 (Gutierrez et al. 2012). These withdrawals for irrigation are unknown but represent small volumes (family farms, very shallow wells).

For the estimation of water extraction for industrial uses, the survey carried out by the local team identified only eight abstraction points, related to the production of sachet drinking water. Once again, the information proved difficult to collect/access and cross check. The estimation for the Caracol industrial park is notably underestimated.

In total, 2109 abstraction wells were identified thanks to the Haiti Outreach database and introduced in the modelled area. These wells represent an average annual withdrawal estimated at 9.5 Mm3/year, including 9.3 Mm3/year dedicated to drinking water supply.

Hydrodynamic parameters

The aquifer of the Plaine du Nord has essentially been investigated during the United Nations project at the end of the 1980s (Tomasi 1990b). Pumping tests performed allowed an overall assessment of the hydraulic transmissivity in different points of the aquifer. Most of tests concern deep wells (depth > 40 m). The hydraulic transmissivity has been estimated to range between 5.2 × 10−4 and 2.2 × 10−2 m2/s, with an average of 6.6 × 10−3 m2/s. For this study however, no information was collected about the storage coefficient for the shallow aquifer in general.

Hydrographic network

The five main rivers crossing the alluvial plain—Haut du Cap, Grande Rivière du Nord, Trou du Nord, Marion and Jassa—were modelled explicitly with the software. For each tributary, the river data are defined by their section within the cell. The hydrological priorities are defined by a method of branching tree that arranges the tributaries’ connectivity and forms the modelled hydrographic network. A river section is described in the model by its geometry (length, width, and elevation of the riverbed), by its physical properties (thickness and permeability of the bed and the banks), by its source terms (pumping or injections) and by its state (height of water in the river). The geometric characteristics are taken from the field surveys carried out by DAI water technicians.

In addition to the river network within the study, an estimation of the upstream inflows at the boundaries of the modelled area was required for each main river. However, there are currently no flow-monitoring gauging stations upstream of the flat Plaine du Nord. The upstream river inflows were therefore estimated by calculating the effective precipitation in the watershed of each upstream mesh using the Thornthwaite method.

Steady-state model calibration

The three main parameters adjusted for the wedging were—the permeability, and the hydraulic heads on the fixed-head boundary and the layer of water infiltrated on the plain to a lesser extent. For each parameter, a trial-and-error approach was used.

A homogeneous initial permeability of 1e-5 m/s was chosen on the basis of the available data. The adjustment of permeabilities was then based on areas of high electrical resistivity revealed by geophysical investigations performed in 1990 (Radstake and Chery 1992). The areas identified as being more resistive correspond to coarser, and therefore more permeable, alluvium, something that has been occasionally confirmed by various reconnaissance drilling endeavors. During the calibration process, the delimitation of these homogeneous permeability zones was adjusted in order to best reproduce the observed piezometric levels. The aquifer is considered isotropic. This study focused on uniform value ranges rather than detailed reproduction of field data leading to a contrasting mosaic of values, because regional models can not consider heterogeneities observed at small scale. The distribution of hydraulic conductivity finally adjusted is shown in Fig. 7.

Fig. 7
figure 7

Spatial distribution of the calibrated hydraulic conductivity

Full size image

For the hydraulic heads on the fixed-head boundary, an adjustment on different segments was necessary in order to respect the consistency of the inflows into the alluvial plain and the piezometric levels observed on the upstream edge of the plain.

Regarding groundwater recharge, the zones and average annual values proposed by Adamson et al. (2016) were retained (see Fig. 5), with the exception of the central zone, where the value was lowered from 200 to 180 mm/year.

The field observations available and used to constrain the steady-state calibration process are—the piezometric map developed by Tomasi (1990c) and revised within the framework of this project (Fig. 2), the piezometric levels noted by DAI water technicians and BRGM (506 measurement points), the annual average flows at four historic stations, the gauge measurements carried out for the AVANSE project (DAI Global, LLC 2020) and maps of the wetlands and flood zones (Fig. 8). The initial hydraulic heads are only given on the fixed head boundary. For the rest of the domain, the initial hydraulic heads are fixed to 0 msl.

Fig. 8
figure 8

Head and hydrometric data available for the steady-state model calibration

Full size image

Simulation scenarios

As previously noted, the lack of monitoring data (i.e. piezometric levels, river flows, withdrawals, rainfall and evapotranspiration (ET)) represents a real challenge, and therefore does not allow for the calibration of the model for a transient state. However, the AVANSE project required an understanding of the seasonal fluctuations of the shallow water table, which are observed to be fundamental in the continuous supply of water to the many hand pumps that are used on a daily basis. To this purpose, a theoretical transient-state model has been developed based on the following approach:

  1. 1.

    Development of a “theoretical” transient regime over the period 2000–2016, at monthly time steps, which represents the seasonal and annual fluctuations (scenario 1)

  2. 2.

    Simulation of a theoretical scenario of irrigation sized by the AVANSE project agronomists, mainly based on the water needs of the banana crops according to the FAO CROPWAT irrigation demands 9 (scenario 2)

  3. 3.

    Comparison of the irrigation scenario 2 with the reference scenario 1 for the same period 2000–2016

Assumptions for the theoretical transient model (scenario 1)

The results of the steady-state calibration for the piezometry and the flow rates in the five rivers were used as initial conditions for the transient regime. The temporal variability of the meteoric recharge at the level of the plain is taken into account, based on a calculation of effective rainfall at monthly time steps over the period 2000–2016. The calculation was performed based on monthly cumulative rainfall records for Cap Haitien and Terrier Rouge (source: USAID/MARNDR), while the average monthly potential evaporation (PET) values were calculated using the FAO records. A correction coefficient was applied to get to the theoretical average monthly recharge for the period 2000–2016. The correction coefficient was applied according to the following formula:

$$ {R}_{(i)}=\left({\mathrm{PEFF}}_{(i)}\div {\mathrm{PEFF}}_{\mathrm{avg}(j)}\times {R}_{\mathrm{avg}(j)}\right) $$
(4)

where the terms are defined as follows:

PEFF(i):

Effective rainfall computed for month i with i = 1–204 (i.e. January 2000 to December 2016)

PEFFavg(j):

Mean monthly effective rainfall calculated over 2000–2016 (j = 1–12)

Ravg(j):

Mean monthly recharge used in the transient model (j = 1–12).

A recharge chronicle at monthly time steps (2000–2016) is thus assigned to each of the five recharging zones of the model. (Fig. 9). This seasonal recharge fluctuation was verified with some records for the same period collected at the AVANSE farms. The meteoric recharge introduced in the model shows a large deficit during the drought period observed between 2014 and 2016, with even zero recharge in some places.

Fig. 9
figure 9

Evolution of annual recharge in the transient model (2000–2016)

Full size image

Concerning the inflows coming into the model boundaries, the lack of climatic data on the Massif du Nord for the period 2000–2016 was overcome by adopting an even more schematic approach. A monthly distribution, identical for each year, was defined from the calculation of the effective monthly average precipitation based on old rainfall and PET data (MARNDR), representative of the Massif du Nord. For each cell concerned, the calculated inflow in steady state was assigned a distribution coefficient.

For the transient simulations, in absence of measurements, the hydraulic heads that have been fixed at the upstream limits are therefore replaced by fixed flows. This approach ensures that the inflows simulated by the modelling software do not increase unrealistically in the event that a significant drawdown affects the modelled area (under the effect of pumping near the limits, for example). All other boundary conditions remain unchanged (with loads imposed on the coast, lateral limits with zero flow, meshes with overflow over the entire domain). Additionally, groundwater withdrawals concerning water supply and industrial uses are considered constant throughout the year. Finally, in the transient mode, a storage coefficient value must be assigned to the model cells. An average value of 2% was chosen for the unconfined state.

Theoretical irrigation scenario (scenario 2)

Eleven areas targeted for potential future development of banana cultivation were studied. The water requirements were estimated from FAO guides (CROPWAT 8.0 software). In Haiti, the irrigation period lasts on average from June to October (Table 2).

Table 2 Water requirements of banana crops (FAO)
Full size table

In accordance with AVANSE agronomists, the assumption was made that irrigation would be satisfied by setting up boreholes, such that one borehole structure would satisfy the irrigation requirements for a plot of 10 ha. In the identified areas for banana crop, an additional total of 270 boreholes were thus implemented in the modelled area to meet the water needs of the 11 zones identified by AVANSE (Fig. 10).

Fig. 10
figure 10

Locations of the target areas and boreholes implemented for the theoretical irrigation scenario

Full size image

Based on these assumptions, the unit production rate of an irrigation borehole varies between 2.2 m3/h (June) and 13.6 m3/h (July). All uses combined, withdrawals are thus brought to a maximum of around 3.5 Mm3 in July, against 0.79 Mm3 for the reference model.

In the scenario, agricultural boreholes were commissioned from 2007 onwards. In order to measure the impact of this pumping scenario on the water resource of the alluvial aquifer, the reference simulation (scenario 1 without agricultural abstraction) was compared with the prospective simulation (scenario 2 with agricultural pumping) for the period 2007–2016, i.e. over a 10-year regular irrigation period.

Results and discussion

Steady-state calibration

The simulated steady-state piezometric map is shown in Fig. 11. The piezometric contours from the Haitian hydrogeological map (MARNDR and UNPD 1990) are also shown for comparison. The results obtained from the modelled piezometric levels are satisfactory with well-simulated flow directions and hydraulic gradients. The scatter diagram in Fig. 12 shows how close the simulated contours area aligned with the field measurements and the absence of an underestimation or a systematic overestimation of the model with deviations of less than 6 cm. The results for the simulation of the watersheds of Haut-du-Cap, Grande Rivière du Nord and Trou du Nord are particularly good and match the measured piezometer levels of the 1990 hydrogeology map well. The isocontour lines between 2 and 20 msl are almost superimposed. Towards the east, the simulation is less accurate around Terrier Rouge and for the Marion and Jassa watersheds. The comparison of simulated and observed piezometer levels in these areas is in any case less accurate since there were fewer field measurements undertaken in these areas, and therefore the accuracy of the 1990 potentiometric map itself is less certain. Only additional field measurements in those areas could help improve the calibration of the model in these basins.

Fig. 11
figure 11

Piezometric map calculated in steady-state conditions vs piezometric map produced by Tomasi 1990c

Full size image
Fig. 12
figure 12

Simulated heads in steady-state conditions vs measured heads

Full size image

Although the model reproduces piezometer levels well, the scarcity of field observations on the upstream end of the plain could be improved by undertaking additional baseline data and by consolidating the steady-state piezometric calibration.

From a scientifical perspective, the calibration of the steady river flow does not make sense; however, this was performed to check the consistency of the simulated surface flows with the measured values. The only available historical hydrometric data span the period 1920–1940 for the main rivers, while more recent gauging data were collected through the AVANSE during 2017 (five locations with weekly measurements from January to July) from which an average measured flow was calculated (Table 3).

Table 3 Simulated river flows in steady-state conditions vs measured flow in 2017 (GRN: Grande Rivière du Nord, HDC: Haut du Cap river, TDN: Trou du Nord river)
Full size table

Even if deviations exist for most of the measurement points, the orders of magnitude of the simulated flows are acceptable (Table 3). Deviations are inevitable given different measurement periods and uncertainties linked to the gauging operations on the one hand and the model on the other hand.

The absence of proper description of the hydrographic network does affect the calculation of groundwater and river exchanges and introduces uncertainties in the model results. The most important river parameters to survey are the stream morphology and the geometry of the riverbed. A topographic survey campaign (cross sections) of the main five rivers was launched by the end of 2017 by the AVANSE team; however, at the time of the assignment, the data were not ready. The newly acquired data could be integrated into the model at a later stage.

The simulation of surface flows can only be consolidated from a calibration of the model in the transient regime. Therefore, the modelled exchange flows between the five main rivers and the alluvial groundwater should be considered with caution.

Water balance

The diagram of Fig. 13 represents the main components of the hydrological balance of the alluvial aquifer at the scale of the model. The share of withdrawals remains moderate overall compared to the flows contributing to the groundwater supply (infiltration + contributions from the Massif du Nord). At the level of the plain, it can therefore be concluded that the availability of groundwater resources is still significant. As will be shown later, this overall result at the level of the plain does not locally exclude operational difficulties, in particular a risk of saline intrusion. Indeed, the question of readily available water volumes cannot be considered independently from spatial (location of withdrawals) and temporal (high water/low water) considerations. Most parts of the five modelled rivers are indeed draining groundwater (rivers are fed by the groundwater). The exchange flows from shallow groundwater to rivers represent 85% of the exchanged flows (Fig. 13).

Fig. 13
figure 13

Water balance in steady-state conditions (negative values represent water leaving the aquifer, whereas positive values represent water entering the aquifer). MCM million cubic meters

Full size image

The water contributions from areas outside of the modelled domain, mainly by water runoff on the Massif du Nord, are always difficult to quantify in the field. However, in order to validate the model’s calculations with regard to inflows from upstream, the continuous acquisition of hydro-climatic data is essential, at least by establishing monitoring of—(1) the flows of the main rivers at their outlet into the plain, (2) piezometric fluctuations upstream of the alluvial plain; (3) rainfall over the massif. In summary of the overall water balance in the study area, a better understanding of withdrawals, both in aquifers and rivers, is essential. Some recent borehole measurements could not be included in the inventory (Caracol industrial park, large farms, etc.).

Results of simulation scenarios

Reference scenario

In a transient regime, the model allows one to simulate and understand the dynamics of groundwater in relation to the variable infiltration of precipitation on the plain and inputs from the Massif du Nord. Annual cyclicality predominates and recharge periods are short, spanning mostly 1–3 months, and generally occur between October and December, which is then followed by a period of groundwater drainage that extends from January to September (Fig. 14).

Fig. 14
figure 14

Transient groundwater dynamics in the Plaine du Nord during the period 2000 to 2016 with a variation of cumulative aquifer storage simulated and b variation of simulated heads in Grande Rivière du Nord

Full size image

At the beginning of 2015, the two successive dry years impact the water-resource balance, with the alluvial aquifer almost continuously draining (destorage) water out of the aquifer, which results in noticeable absence of the usual rise in piezometric levels (Fig. 14).

A further example is presented in Fig. 15, which presents the simulated piezometric map for the month of June 2016, 1 year in the drought period. It is interesting to note that, in a context of a recharge deficit, the piezometric depression of Cap Haitien is spreading, in particular at the mouth of the river, where simulated groundwater levels are below sea level. The model clearly shows that under very low water conditions, as observed in 2016, the risk of saline intrusion is inevitably increased in the areas where the groundwater is unsustainably exploited.

Fig. 15
figure 15

Piezometric map calculated in transient-state conditions (June 2016)

Full size image

Irrigation scenario

The impacts of agricultural pumping on groundwater levels differ from one study area to another. The drawdowns appear to be better controlled around the Cap Haitien and Limonade sectors, with simulated maximum drawdowns of around 1 and 0.7 m, respectively (Fig. 16). In contrast, the piezometric drawdown is greater in the areas of Caracol, Terrier Rouge and Ouanaminthe, with drawdowns that can exceed 4 m and reach up to 6 m in Caracol (Fig. 16).

Fig. 16
figure 16

Simulated drawdowns for the irrigation scenario (August 2016), with a variation of simulated heads in Cap Haitien area (for a 1-m diameter well), b variation of simulated heads in Caracol area (for a 1-m diameter well), and c simulated drawdown map

Full size image

It is interesting to also note that the simulated piezometric levels rise again as soon as the irrigation is stopped, to quickly get back to the levels of the reference simulation (Fig. 16). This remains true during periods of pronounced low water such as during the drought of 2014–2016, with the exception of certain sectors located at the upstream end of the plain (Terrier Rouge). For the latter, the water inflow from Massif du Nord does not appear to compensate for the additional agricultural pumping under low recharge conditions.

Finally, from an operational point of view, the two main operating constraints highlighted by the simulation of such a pumping scenario are as follows:

  1. 1.

    Along the coast, localised operating difficulties could be expected due to a significant risk of seawater intrusion. Indeed, under the effect of agricultural pumping, upconing phenomena and/or reversal of groundwater flow direction cannot be excluded. In particular, the areas of Cap Haitien and Caracol that have simulated piezometric levels close to 0 msl, would be immediately impacted by seawater intrusion when the water table drops below sea level. This conclusion is based on the simulation of drawdown by the model. To develop this point, it might be interesting in the future to introduce into the model variable-density groundwater flow (possibility provided by the MARTHE software).

  2. 2.

    The simulated drawdowns amount to several meters in most cases (without taking into account the quadratic head loss). It is therefore likely that it will not be possible to rely on shallow wells to meet crop water requirements. According to the model, the most sensitive areas are around Caracol, Terrier Rouge and Ouanaminthe. The cumulative impact of intensive agricultural and industrial drawdown, combined with the potential effect of climate change, could potentially lead to unsustainable use of groundwater that would not allow for enough time for recharge and a return to a safe water balance status. An integrated use of water resources is therefore recommended.

Conclusions

The Plaine du Nord aquifer represents a major water resource for Haiti. Agricultural development projects combined with urban development in the region require an integrated approach to the management of the alluvial aquifer and the river network that drains it. A stable hydrodynamic model was developed based on the scarce data available and new data acquired during the AVANSE project. Due to the lack of monitoring of the groundwater, watercourses or even climate, the calibration of the model could only be carried out under steady-state hydraulic conditions. The results, however, are promising and the work presented in this report is only the first stage of a much wider work package that could be complemented by the calibration of the transient state that would allow one to better take into account the dynamics of storage/destorage of the aquifer in the various locations in the Plaine du Nord. Additional survey work should be considered, and the new data acquired can easily be introduced into the model.

The steady-state model has, however, already shed new light on the groundwater resources of the alluvial plain. The development of the hydraulic balance of the groundwater is showing that significant water resources can be easily mobilized at the level of the plain, if carefully managed to rule out localised operational difficulties.

The potential for groundwater withdrawals is more limited for the Haut-du-Cap watershed due to an already higher pressure for withdrawal to supply drinking water to Cap Haitien, and at the Marion River watershed, due to a lower recharge potential, whether by precipitation or by water flows from the Massif du Nord. The simulated transient scenarios show that the alluvial aquifer conditions allow for the commissioning of agricultural irrigation in the identified areas by the AVANSE project. A rapid replenishment of the water resource was simulated after the irrigation campaign was stopped. This was still the case in times of drought, with the exception being made of the areas upstream of the Plaine du Nord.

The simulated groundwater drawdowns at the additional pumping wells are pluri-metric, reaching 6–7 m. Given that the current wells are very shallow (3–6 m depth on average), such drawdowns would result in the current wells becoming dry. The siting of deeper boreholes (>20 m) should be considered, especially in a context of deficit recharging such as the one observed during 2014 and 2016 to allow for a sustainable balance between extraction and recharge. Finally, agricultural areas that are close to the town of Cap Haitien and Caracol are particularly vulnerable to saline-water intrusion from the ocean.