Hydrogeochemical evolution and groundwater mineralization of shallow aquifers in the Bas-Saint-Laurent region, Québec, Canada

Abstract

This study presents the first regional groundwater hydrogeochemical portrait of the Bas-Saint-Laurent region (BSL), a region shaped by the Appalachians, a strong Quaternary glacial heritage, and coastal dynamics from the St. Lawrence Estuary. The proximity of BSL’s aquifers to St. Lawrence Estuary and its geological history with the past Goltwait Sea transgression create unique issues with respect to groundwater mineralization and sustainability in the region. The study is based on the distribution of major and trace elements and stable isotope signatures of water and inorganic carbon (δ¹⁸O, δ²H, δ¹³C_DIC) in 145 groundwater samples collected in private and municipal wells distributed evenly over the study area. Groundwater shows a wide range of composition as indicated by the seven facies revealed by their composition of major elements. Ca-HCO₃ and Na-HCO₃ facies mainly dominate the regional groundwater composition, representing respectively 66 and 20% of the samples. Nevertheless, no significant relation between the geology, the aquifer confinement, and the geochemical facies emerged. This suggests that factors other than the hydrological settings may control the chemical composition of the groundwater in the study area. A hierarchical cluster analysis (HCA), including major, minor and trace elements, was performed, allowing the water samples to be distributed into four distinct geochemical groups that reveal a gradient from less mineralized (C4 and C2 groups with a dominant Ca-HCO₃ facies) in the recharge areas to more mineralized (C1 with a Ca-HCO₃ facies, to C3 with a Na-HCO₃ facies) in the coastal discharge areas. Based on geochemical graphs and isotopic signatures, a conceptual model is proposed to explain this hydrogeochemical evolution at the regional scale. The most remarkable finding is that groundwater mineralization does not originate from modern seawater mixing despite the proximity of St. Lawrence seawater. Most of the hydrochemical evolution and groundwater mineralization is induced by the mixing with evaporated or remnant seawater originated from past transgressions, cation exchanges and mineral dissolution.

L’étude propose le premier portrait hydrogéochimique des eaux souterraines de la région du Bas-St-Laurent, une région physiographiquement marquée par la présence de la chaine Appalachienne, un fort héritage de l’histoire quaternaire et la dynamique côtière. La proximité actuelle avec les eaux salées de l’Estuaire du St-Laurent ainsi que les anciennes transgressions de la Mer de Goltwait sur le territoire soulèvent des questions particulières quant aux mécanismes de minéralisation des eaux souterraines à l’échelle régionale et à leur préservation. Ce portrait se base sur la distribution des ions majeurs, mineurs et en traces, ainsi que sur la signature isotopique de l’eau et du carbone inorganique dissous (δ¹⁸O, δ²H, δ¹³C_DIC) de 145 échantillons prélevés dans des puits municipaux et privés, répartis dans toute la région. Basé sur les ions majeurs, sept facies géochimiques se distinguent. Les facies Ca-HCO₃ et Na-HCO₃ sont prépondérants, représentant respectivement 66 et 20% des échantillons. Aucune relation significative n’a cependant pu être observée entre les facies, la géologie et le degré de confinement des aquifères, ce qui suggère que d’autres facteurs pourraient contrôler la composition géochimique des eaux souterraines dans l’aire d’étude. Une classification hiérarchique des échantillons basée sur les ions majeurs, mineurs et en trace révèle quatre groupes qui se distinguent par leur degré de minéralisation : les moins minéralisés dans les zones de recharge des aquifères (C4 et C2 dominés par un facies en Ca-HCO₃) aux plus minéralisés dans les zones de décharge côtière (C1 avec un facies en Ca-HCO₃ et C3 en Na-HCO₃). En se basant sur des graphiques binaires et les signatures isotopiques, un modèle conceptuel est proposé pour expliquer l’évolution hydrogéochimique et la minéralisation des eaux souterraines de la région. Malgré la proximité des eaux salées du St-Laurent, la minéralisation observée ne semble pas résulter du mélange avec des eaux salées actuelles. L’évolution hydrochimique et la minéralisation observée est principalement induite par le mélange avec de l’eau de mer évaporée ou restée piégée dans les pores des roches à la suite des transgressions passées, des échanges de cations qui en résulte et de la dissolution des minéraux les plus solubles.

Introduction

The origin and the geochemical background of groundwater are key components to an adequate assessment of water quality at a regional scale (Edmunds et al. 1987, 2003; Cloutier et al., 2006). In coastal regions, groundwater salinization occurs in many aquifers around the world and leads to a deterioration in water quality (Oude Essink et al. 2001; de Montety et al. 2007; Alcalà and Custodio 2008; Barlow and Reichard 2010). Understanding the origin of the geochemical background as well as the processes responsible for the mineralization of coastal groundwater resources is an important starting point for preventing further deterioration and anticipating future changes (Downing and Wilkinson 1991). Groundwater mineralization is often due to (1) the inflow of current seawater (Oude Essink 2001; Pulido-Leboeuf 2004) or deep saline water (fossil water; Cloutier et al. 2006) due to overexploitation of the aquifer, (2) water pollution due to extensive irrigation or well failures (Foster and Chilton 2003), or (3) evaporite dissolution (Celle-Jeanton et al. 2009). Because of the different causes of groundwater quality deterioration, it is necessary to identify and characterize the specific mechanisms involved.

In addition to these salinization mechanisms, the geological history also influences the regional geochemical signature of groundwater (Banks et al. 1998; Louvat et al. 1999; Edmuns et al. 2003). In North America, the Quaternary glaciations and the last deglaciation shaped the landscape and created large and complex sedimentary bodies that act as aquifers or aquitards (Occhietti 1989; Occhietti et al. 2011). The mineralogy of Quaternary deposits, residence time, and water flow paths within the deposits directly influence the geochemical signature and mineralization of groundwater (Couture 1997; Cloutier et al. 2006, 2008, 2010; Blanchette et al. 2010; Beaudry et al. 2011; Montcoudiol et al. 2015). In the St. Lawrence Lowlands, for example, Cloutier (2004), Cloutier et al. (2006, 2008, 2010) and Occhietti et al. (2011) highlighted the contribution of three different but concomitant mechanisms affecting groundwater signatures: (1) cation exchanges, which take place between fresh groundwater recharge and clay deposits left by the Champlain Sea between 13,000 and 11,000 years BP; (2) water–rock interaction with the dissolution of marine-derived particles such as carbonate and gypsum; and (3) brine dissolution at depth or mixing with ancient seawater, probably related to Champlain Sea water.

This study presents the first regional groundwater hydrogeochemical portrait of the Bas-Saint-Laurent region (BSL; Québec, Canada), a region shaped by the Appalachians, the Quaternary glacial heritage (Hétu 1994a, 1994b, 1998), and coastal dynamics of the St. Lawrence Estuary. The main goal was to identify and understand how the different mineralization mechanisms affect groundwater quality. Specifically, the study explored the role of the modern seawater intrusion (e.g. seawater from the St. Lawrence Estuary), the potential contribution of deep saline waters (e.g. remnants of the ancient Goldthwait Sea), water–rock interactions with sedimentary rocks (e.g. marine clays, Quaternary deposits, and/or sedimentary substrata), and dissolution of evaporite minerals (e.g. carbonate and salt). Major ions and isotope studies, including stable isotopes of water and carbon (δ¹⁸O, δ²H, and δ¹³C_DIC), were analyzed in 164 water samples collected in Quaternary deposits and sedimentary substrata of the Paleozoic Era. This work is a contribution to the Programme d’Acquisition des Connaissances sur les Eaux Souterraines (PACES), a vast programme studying groundwater in municipal regions in the province of Québec. The programme sought to develop regional hydrogeological portraits focusing on the quantity, quality and vulnerability of groundwater to promote sustainable management practices of groundwater resources.

Hydrogeological context of the study area

The BSL region is bounded by the St. Lawrence Estuary to the north and New Brunswick and Maine (USA) to the south (Figure 1). The study area covers 4000 km² of the 22,200 km² of the BSL and is mainly along the St. Lawrence shores. The average annual temperature is 3.5°C (average temperature of 17.5°C in July and −12.2°C in January), and the average annual precipitation is 933 mm, including the water equivalent of 325 cm of snow received each winter (Mont-Joli meteorological station, Environment Canada 2015).

Figure 1. A, Location of the study area in the Bas-Saint-Laurent region, Québec, Canada, with the geological and hydrogeological context of the study area. B, Schematic drawing of the principal hydrostratigraphic units present along cross section 1–2, from the coastal plain to the highlands.

The study area is part of the vast Appalachian Mountains region (Bostock 1967) and can be divided into two major physiographic regions: the coastal plain and the highlands. The coastal plain is generally flat and covers the area adjacent to the current shoreline and those that were flooded by the Goldthwait Sea at the end of the last glaciation. The highlands consist of rounded hills, Appalachian ridges that are parallel to the St. Lawrence Estuary, and high outcrop plateaus with an average altitude of 300 m (Hétu 1994a). Valleys filled with alluvial deposits at the downstream end of drainage networks are also elements of the physiographic configuration (Marchand et al. 2014).

The sedimentary substrata of the area are part of the Appalachian Orogen (Williams 1995) and can be divided in two tectono-stratigraphic units of different ages separated by the Neigette fault (Figure 1A; Lajoie 1972; Vallières 1975; Lavoie et al. 2003; MERNQ 2012). North of the Neigette fault, Cambrian–Ordovician sedimentary rocks from the Taconic orogeny are characterized by alternating ‘soft’ (argillite, shale, siltstone and slate) and ‘harder’ (conglomerates and sandstones) rocks. The numerous rocky outcrops characterizing the Appalachian regional relief are composed of harder rocks. Large fracture networks are absent from these bedrock aquifers, resulting in hydraulic conductivity ranging between 2.9 × 10⁻⁶ and 8.9 × 10⁻⁷ m/s. South of the Neigette fault, Silurian–Devonian sedimentary rocks from the Acadian orogeny are mainly composed of conglomerates, sandstones and limestones. These harder bedrock aquifers have a fracture network slightly more developed than that of the Taconic orogeny, resulting in higher hydraulic conductivity (from 2.0 × 10⁻⁵ to 2.7 × 10⁻⁶ m/s). These two units represent the bedrock aquifers of the area.

Most surficial deposits in the area occurred following the retreat of the ice fields from the last glacial period. The history following the retreat of the ice fields supports the division of the study area into the two main physiographic regions that characterize the dominant hydrogeological contexts.

The coastal plain is characterized by a succession of Appalachian ridges filled by marine, glacio-marine and alluvial deposits. This stratigraphic sequence, which is typical of Appalachian ridges, allowed the establishment of unconfined granular aquifers and confined bedrock aquifers (Figure 1B). Toward the end of the last ice age (i.e. between 18,000 and 15,000 BP), a powerful ice stream flowing to the northeast settled in the St. Lawrence Valley. The Goldthwait Sea then diachronically invaded the territory in parallel with the upstream migration of a calving bay in the St. Lawrence Valley (Dionne 1977, 1997; Hétu and Gray 2002). The marine transgression started around 15,000 BP in the northeast part of the territory (Figure 1A; Matane and Rimouski) and reached the southwest sector (Rivière-du-Loup; Figure 1A) around 14,300 BP. The marine limit decreases progressively from 165 m above mean sea level in the southwest part of the study area (e.g. Saint-Fabien; Dionne 1972, 1977, 2002; Locat 1978) to less than 105 m in the eastern part near Matane (Dionne and Coll 1995; Hétu and Gray 2002).

During the transgressive marine phase, thick sequences (10–40 m) of glacial marine and marine clays were deposited on the sedimentary substratum in the lower part of the coastal plain (< 40 m). The bedrock aquifers in the lower part of the coastal plain are thus largely confined and their recharge is low – or even zero in some cases – because of the thickness of marine clay deposits on top. In some locations, the marine clays cap sandy and gravelly deposits, producing confined granular aquifers that are highly productive and are used by municipalities as their main source of drinking water. One such location is the Neigette River aquifer for the city of Rimouski (Figure 1A).

From the maximum of the transgressive marine phase and during the following regressive phase, a large variety of sedimentary deposits were laid down in both the lower and upper parts of the coastal plain. In the lower part (altitude < 40 m), the glacio-marine and marine clays are topped with granular littoral deposits forming vast coastal terraces, especially in the areas of Trois-Pistoles, Rimouski, Sainte-Luce, Mont-Joli and Sainte-Flavie (Figure 1B). The littoral deposits are highly permeable (average hydraulic conductivity 5.1 × 10⁻⁵ m/s), but generally have low saturated values. These unconfined granular aquifers are frequently accessed by shallow wells. The invasion by the Goldthwait Sea also led to the establishment of vast deltas fed by glacial river and stream waters (Dionne 1972; Locat 1978; Hétu 1994b 1998). These unconfined granular aquifers are generally located in the upper part of the coastal plain; they have high hydraulic conductivity (average K = 1.4 × 10⁻⁵ m/s) but are very rarely saturated. Finally, river valleys near the coast are also generally covered with sand and fluvial gravels deposited on the top of marine clays, which are themselves based on rock. Significant glacio-fluvial and alluvial deposits are found in the Matane and Mitis valleys. In the Matane valley, which is 1 km wide and more than 50 km in length, alluvial deposits (2–5 m thick) rest on top of deep glacio-fluvial deposits (20–50 m thick). This is the largest unconfined granular aquifer in the study area.

In the highlands, the bedrock aquifer lies beneath thin till or alteration layers (1–5 m) and is thus unconfined where direct recharge of the aquifer can take place. The highlands include all areas in the southern BSL where altitudes are greater than 150 m above mean sea level. The highlands occupy 70% of the study region and correspond to high dissected plateaus with rounded hills and alternating Appalachian ridges. The bedrock aquifer is the main recharge area of the region.

Materials and methods

Sampling strategy

Water samples were collected in 164 wells during the summers of 2013 and 2014. Most samples were from private wells (82%) while the remainder were observation wells installed in 2013 and 2014 (13%) and wells that are included in the provincial groundwater monitoring network (5%). Most samples were collected from wells installed in bedrock aquifers at an average depth of 35 m (67%); the remainder were from wells in unconfined granular aquifers with an average depth of 19 m (29%) or from surface groundwater springs (4%).

All samples were collected using the protocol established by Blanchette et al. (2010). Depending on the depth of the groundwater level, sampling was done using a peristaltic pump or a submersible pump with a low flow rate (~ < 6 L/min). Physico-chemical parameters (temperature, pH, salinity, dissolved oxygen saturation [%]) were continuously measured using a calibrated multiparametric probe until stable readings were obtained. Samples were filtered onto 0.45-μm Whatman Polycap filters and stored in separate 125-mL polyethylene bottles for further analyses. Water samples for the analysis of cations and anions were stored with 70% Suprapur® HNO₃ reagent to prevent cation precipitation; water samples for nutrients were acidified with Suprapur® HSO₄; and Zn acetate was added for sulphide analyses. All water samples were stored at 4°C and shipped to the laboratory at the end of the day. Duplicate samples (~5% of the total) were also submitted to verify the data quality and accuracy. Additional water samples were collected without filtration for the analysis of stable isotopes (δ¹⁸O and δ²H) and ¹³C_DIC. For stable isotope analyses, 30-mL scintillation vials filled with groundwater were sealed and stored at room temperature; 123 water samples were taken. In 2013, 34 water samples were collected for δ¹³C_DIC: 250-mL glass amber bottles were filled with no headspace, 0.5 mL of saturated HgCl₂ was added for preservation, and bottles were stored at 4°C.

Four water collectors adapted from models developed by the International Atomic Agency (2008) and Gröning et al. (2012) were installed to collect precipitation (liquid and solid). Precipitation was sampled monthly from May to November 2014, and was preserved according to the protocol used for groundwater isotope analyses (δ¹⁸O and δ²H).

Analytical methods

Analyses of alkalinity, sulphides, and 39 ions including major ions (Ca, Mg, K, Na, Br, Cl, NO₃, SO₄) were conducted by a private company (Maxxam Analytics Inc., Montréal) using standard methods and in accordance with PACES guidelines (GRIES 2010). The concentrations of cations were determined using inductively coupled plasma mass spectrometry (ICP-MS). The analytical methods were continuously quality checked by parallel analysis of international certified reference water. The precision was generally better than 20%. Anions were analyzed using ion chromatography (IC). The precision was better than 10%. Total dissolved solids (TDS) were calculated from the content of major ions, and HCO₃ content was calculated based on pH and total alkalinity measurements. Analytical results were validated by testing the ionic balance of each of the 164 samples. According to Hounslow (1995), an ionic balance of ± 10% is considered acceptable when using the results for further statistical analysis. The balance of cations and anions suggested that 19 of the 164 samples did not meet Hounslow’s criterion; thus, a subgroup of 145 samples was used for the remaining geochemical and multivariate analyses.

Stable isotopes of water (δ¹⁸O, δ²H) were analyzed by elemental analysis–isotope ratio mass spectrometry (EA-IRMS). Precisions were ± 0.05‰ and ± 1‰ (at the 1σ level) for δ¹⁸O and δ²H, respectively. Isotopic analyses are reported compared to the international Vienna Standard Mean Ocean Water (VSMOW) as defined by Gonfiantini (1978). Reference materials were used throughout the isotopic water analyses to ensure high-quality data. The isotope ratio of total dissolved inorganic carbon, δ¹³C_TDIC, was analyzed using gas chromatography (GC)-IRMS. An internal calibration was performed using two internal standards to obtain a two-point calibration. Results are reported in parts per thousand in the δ notation relative to Vienna Pee Dee Belemnite (standardized to TS-Limestone [NEBS19] and lithium carbonate [LSVEC] AIEA referenced materials) with a precision of ± 0.1‰.

Saturation index (SI)

The solubilities of carbonate species and other mineral phases as calcite, dolomite and gypsum, respectively, were calculated using the geochemical code PHREEQC (Parkhurst and Appelo 1999). Saturation index (SI) values around zero indicate that a water sample is in equilibrium with the mineral phase. Therefore, when SI values are higher than zero, the groundwater is considered supersaturated with respect to the mineral and precipitation may occur. When SI values are lower than zero, the groundwater is considered undersaturated with respect to the mineral; the mineral cannot precipitate and dissolution may occur.

Statistical methods

As proposed by Güler and Thyne (2004), graphical methods were combined with multivariate statistical analyses to characterize the regional groundwater hydrogeochemistry. The combination of these methods permits the rapid determination of the geochemistry of various facies of water and allows one to explore and explain the various processes that control groundwater mineralization (Drever 1997; Cloutier et al. 2008). Graphical methods are also used to characterize the geochemical data.

In order to use the most representative parameters for a larger number of samples, only 11 of the 39 parameters analyzed were selected for multivariate analysis: HCO₃, Ba, Ca, Cl, Mg, Mn, K, Si, Na, Sr, and SO₄. NO₃ was not included in the analysis because its concentrations were below the detection limit for 50% of the samples. When a sample value of one of the 11 parameters was below the detection limit, it was replaced by a value corresponding to 55% of the detection limit, as proposed by Sandford et al. (1993). This replacement was performed for less than 1.3% of samples, with a maximum proportion of 9% for the manganese element. Finally, because concentrations among parameters were not homogeneous, they were standardized to make them comparable in the multivariate analysis.

Based on this matrix of 145 samples and 11 geochemical parameters, a hierarchical cluster analysis (HCA) was performed to group the samples. HCA is a semi-objective multivariate statistical method leading to the construction of a dendrogram that reveals groups of samples having similar parameter values. The Euclidean distance is used to measure the distance between samples, thereby identifying groups having the highest similarities. The distance between clusters was evaluated according to Ward’s method (Cloutier et al. 2008). Statistical analyses were produced using SPSS software as well as the ggbiplot and devtools packages from Rstudio.

Results and discussion

Geochemical characterization of shallow groundwater

Data from the 145 analyzed samples of groundwater were plotted on a Piper diagram, which revealed a wide range of geochemical compositions (Figure 2). A seawater value was also plotted as a potential source of recharge of the BSL aquifers. Seven geochemical facies can be identified based on major ions, from Ca-HCO₃ to mixed cations-HCO₃. The contribution of each facies is reported in Table 1. The Ca-HCO₃ facies dominated, representing 66% of the samples (Table 1). This facies dominated in almost all geological and hydrological settings and was characterized by low Cl concentrations (< 100 mg/L) and low TDS values (71 < TDS < 588 mg/L). The Na-HCO₃ and Na-Cl facies represent 20% and 7%, respectively, of all groundwater samples. These more mineralized groundwater samples were predominantly from the coastal plain (altitude between 105 to 165 m above mean sea level), which included 70% of the Na-HCO₃ and 100% of the Na-Cl samples. The Na-Cl facies was characterized by high concentrations of TDS (mg/L). One sample had concentrations of Na (1400 mg/L), Cl (2100 mg/L), TDS (8863 mg/L) and Br (8.8 mg/L) that were significantly higher than those of other samples, making it the most mineralized sample. It was collected in a deep confined bedrock aquifer (~68 m depth) located on the St. Lawrence shore near Trois-Pistoles. The Ca-Cl and mixed water (mixed cation-HCO₃, -SO₄) facies made up a small portion of the samples (4% and 1% for Ca-Cl and mixed, respectively). No statistically significant relationship among the geology, the aquifer confinement or the geochemical facies was revealed by the Piper diagram (i.e. there was no significant difference between the geochemical facies in granular and rock aquifers). Thus, it seems that factors other than the hydrological settings may control the chemical composition of groundwater in the study area.

Figure 2. Piper diagram showing water chemistry based on major ions in the granular and bedrock aquifers. The four geochemical groups identified by hierarchical cluster analysis (HCA) are reported. Characteristics of seawater (X) are presented in the diagram.

Table 1. Contribution of the different geochemical facies as revealed by the Piper diagram.

Geochemical facies	No. of samples	Contribution (%)
Ca-HCO3	96	66
Na-HCO3	29	20
Na-Cl	10	7
Ca-Cl	6	4
Mx cation-HCO3	2	1
Mx cation-SO4	1	< 1
Mx cation-anion	1	< 1
Total	145	100

The HCA performed on the 11 selected parameters revealed four geochemically distinct groups (C1, C2 C3, and C4; p < 0.05). This clustering explains 75% of the variance. The distances between C4, C2 and C1 are the smallest, indicating high similarity, while the distance between C3 and C4 is the largest, indicating high dissimilarity between the two groups. Table 2 presents the mean values of the geochemical parameters for each group. The mean values of pH, TDS, and well depths are also reported. The C3 group has the highest values for five ions (HCO₃, Cl, K, Na, SO₄) and the highest mean pH (~8.8), TDS (~792 mg/L) and total alkalinity (~259 mg CaCO₃/L) values. This group includes samples from the deepest wells (25–145 m). For five samples in that group, pH reached 9.5 with total alkalinity (TAlk) higher than 300 mg CaCO₃/L; these samples were collected in wells deeper than 40 m. The C3 group is classified as Na-HCO₃ water with high TDS and TAlk, but five Na-Cl water facies samples were also observed in this group. The ratio of Br to Cl ranged from 0.003 to 0.006 for five of the 16 samples. For the other samples, Br concentrations were below the limit of detection. Most samples belonging to the C3 group were collected in semi-confined to confined aquifers of the coastal plain, below the upper limit of the Goldthwait Sea transgression. Because of the elevated Cl concentrations and trace presence of Br, samples of C3 group could potentially result from salinization processes with the adjacent St. Lawrence Estuary or deep remnant seawater.

Table 2. Mean geochemical and physical parameter values characterizing the four geochemical groups identified by hierarchical cluster analysis (HCA).

Group	C1^*	C2	C3	C4
No. of samples	44	40	16	42
HCO3 (mg/L)	284.0 ± 76.8	159 ± 49	319.6 ± 156.9	133.0 ± 74.1
TAlk(mgCaCO3/ L)	232.5 ± 62	130.7 ± 41	259.2 ± 136	108.0 ± 61
Ba^2+ (mg/L)	0.02 ± 0.2	0.10 ± 0.1	0.5 ± 1.0	0.0
Ca^2+ (mg/L)	78.6 ± 33.6	42.3 ± 17.1	6.3 ± 6.4	37.8 ± 24.0
Cl^- (mg/L)	117 ± 313	22.3 ± 19.9	172.9 ± 267.2	5.6 ± 5.0
Na+ (mg/L)	91.1 ± 207.1	23.9 ± 22.1	236.4 ± 191.6	7.0 ± 5.4
Mg^2+ (mg/L)	16 ± 15	7.1 ± 2.5	4.2 ± 7.4	6.0 ± 3.9
Mn^2+ (mg/L)	0.1 ± 0.4	0.3 ± 0.6	0.0	0.0
K^+ (mg/L)	6.5 ± 11.6	4.0 ± 5.8	5.0 ± 3.7	1.3 ± 0.8
Si^2+ (mg/L)	5.1 ± 9.2	4.2 ± 1.2	2.9 ± 1.4	3.4 ± 1.1
Sr^2+ (mg/L)	0.7 ± 0.7	0.4 ± 0.3	0.1 ± 0.1	0.2 ± 0.1
SO4^2- (mg/L)	35.3 ± 47.6	21.9 ± 10.7	42.5 ± 40.8	10.6 ± 8.8
pH	7.4 ± 0.5	7.3 ± 0.7	8.8 ± 0.6	7.2 ± 0.6
TDS (mg/L)	641 ± 585	315 ± 93.1	792 ± 545	210.9 ± 98.5
Depth (m)	28 ± 22	23 ± 18	50 ± 20	21.2 ± 17.1
Geochemical facies	Ca-HCO3	Ca-HCO3	Na-HCO3	Ca-HCO3
δ^2H ‰^**	−84.6 ± 8.2	−83.8 ± 7.2	na	−85.5 ± 7.8
δ^18O ‰^**	−12.1 ± 1.1	−12.3 ± 0.8	na	−12.3 ± 1
δ^13C ‰^***	−19.7 ± 4.3	−21.4 ± 4.7	na	−17.5 ± 11
Confinement	unconfined to confined	unconfined	semi-confined toconfined	unconfined
Hydrogeological context	highlands	highlands	coastal plain	highlands

Note:

na: no data available.

* Mean values calculated by including the most salinized sample located on the St. Lawrence shore.

** Only 123 samples were analyzed.

*** Only 34 groundwater samples were analyzed.

In contrast to C3, the C4 group presented the lowest values for five ions (HCO₃, Ba, Cl, K, SO₄). C4 samples were collected from the shallowest wells (mean depth ~21 m), and they had the lowest pH (~7.2), TDS (~211 mg/L) and TAlk (~108 mg CaCO₃/L) values. The C4 group is classified as a Ca-HCO₃ water facies with low TDS and low pH. Samples from the C4 group were mainly collected in the unconfined bedrock aquifers of the highlands.

The C2 and C1 groups have characteristics that are intermediate to the C4 and C3 groups. Both C2 and C1 are characterized as Ca-HCO₃ waters, but C1 had higher cation contents (Na, Ca, Mg, K and Sr) and moderate TDS (~ 641 mg/L) and Cl contents (~117 mg/L). The C1 group also included the most mineralized sample, collected near Trois-Pistoles; it is characterized by a high TDS value (4318 mg/L; higher than values observed in C3 group samples), but with a medium pH (7.8) and high SO₄ concentration (310 mg/L), unlike other C3 group samples. This sample is the only Na-Cl water facies sample included in the C1 group. The locations of these water groups are reported in Figure 3.

Figure 3. Locations of the four geochemical water groups based on 11 parameters as revealed by a hierarchical cluster analysis (HCA). The two main hydrological contexts (the coastal plain and the highlands) and the degree of confinement of Bas-Saint-Laurent (BSL) aquifers are also reported.

Figure 4 shows the stable isotopic signatures of groundwater samples. The average composition of the groundwater samples was −12.4 ± 1.0‰ and −86.6 ± 8.3‰ for δ¹⁸O and δ²H, respectively. Groundwater samples plot along a line defined by δ²H = 7.55 × δ¹⁸O + 7.70, which is located slightly above the regression line between the rain and snow end members, presented as a preliminary meteoric water line by Buffin-Bélanger et al. (2015), and very close to the global meteoric water line (e.g. δ²H = 8.13 × δ¹⁸O + 10.38). The absence of a distinct isotopic composition in the groundwater samples suggests that infiltration and recharge of aquifers occurred under modern climatic conditions. Furthermore, the similarity between the groundwater isotopic distribution and the regression line linking modern precipitation (e.g. rain and snow) suggests a common meteoric origin. The average isotopic composition of BSL groundwater was slightly more negative than values reported by Cloutier et al. (2006) in a Paleozoic Basses-Laurentides sedimentary bedrock aquifer (located ~600 km southwest of BSL). The average signature of groundwater samples was also slightly more negative compared to the average signature of all the precipitation samples (−9.4 ± 4.6‰ for δ¹⁸O and −65.2 ± 31.0‰ for δ²H) and appeared to be mainly influenced by the isotopic composition of snow precipitation. Based on the mean signature of the two precipitation end members (e.g. assuming an average composition of −7.5 ± 2.6‰ δ¹⁸O and −53.2 ± 15.7‰ δ²H for rain and −17.2 ± 3.1‰ δ¹⁸O and −122.2 ± 21.9 ‰ δ²H for snow), the percentage of snow is calculated for the different samples using the following equation:

Figure 4. Distribution of δ¹⁸O and δ²H (V-SMOW in ‰) in groundwater facies. The dotted line represents the line linking the stable isotopes of solid and liquid precipitation, and the black line is the linear regression of groundwater samples. The mean values of δ¹⁸O and δ²H are also reported for the regional rain and snow precipitation and seawater (data from Buffin-Bélanger et al. 2015 and Chaillou et al. in press).

where X is the percentage of snow in the sample, C_sample is the composition of the sample, C_rain is the average composition of the same parameter in rain, and C_snow is the average composition of the same parameter in snow. The equation allows the estimation that snow contributed between 33 and 81% of the water sample signatures, with the highest contributions being for samples collected in unconfined granular aquifers. This highlights that snow cover plays a key role in the recharge of the regional aquifer, whatever the degree of confinement of the aquifers. The stable isotopes of all water facies show no difference in relation to the geochemical composition (p > 0.05). There is also no spatial trend toward the coast, suggesting little or limited direct intrusion of modern seawater along the shore.

Mineralization mechanisms: seawater intrusion or water–rock interactions?

Mixing between groundwater and modern or remnant seawater

The proximity of the BSL aquifers to seawater of the St. Lawrence Estuary creates unique issues with respect to groundwater sustainability. These issues are primarily those of seawater intrusion into freshwater aquifers, and changes in the quality of fresh groundwater resources. Seawater intrusion is the movement of saline water into freshwater aquifers; it is most often caused by groundwater pumping from coastal wells, as frequently observed along the east coast of the US (e.g. Richter and Kreitler 1993). The second issue is related to geological history and past seawater transgression, with the mixing of deep (and often stagnant) groundwater with remnant seawater or evaporated seawater still present in marine clay and basement aquifer (Aquilina et al. 2015). This has been observed in other regions of Eastern Canada, particularly in the St. Lawrence Lowlands region (Cloutier et al. 2008; Blanchette et al. 2010). Groundwater salinity was evaluated according to definitions proposed by Rhoades et al. (1992) using TDS and Cl values. In BSL aquifers, only one sample can be classified as saline groundwater (TDS > 1500 mg/L, Cl > 250 mg/L). Most samples from the C4 and C2 groups can be classified as fresh groundwater (e.g. TDS < 500 mg/L, Cl < 100 mg/L). All samples from the C1 and C3 groups, excluding the high-mineralized sample (Figure 3), are slightly to moderately saline groundwater (TDS between 150 and 1900 mg/L). Only six samples exhibited concentrations higher than the Canadian drinking water quality guideline for Cl of 250 mg/L (aesthetic objective of Health Canada 2014), reaching values as high as 800 mg Cl/L.

Since chloride is considered a conservative component (Hill 1984), the major ions are plotted against Cl in Figure 5. The theoretical mixing lines with current seawater of St. Lawrence Estuary (salinity = 29) are reported for binary diagrams of the major ions. In all the binary diagrams, as in the Piper diagram (Figure 2) and the stable isotope composition of water (Figure 4), the C4 group appears as a dilute pristine end member. It is considered as the less-mineralized group originated from the recharge. The distribution of samples from the C2, and C1 to a lesser extent, appears to be unrelated to modern seawater mixing along the St. Lawrence shore, except for the single highly mineralized sample of Na-Cl water facies.

Figure 5. Relationships between different ions and Cl content (in mg/L) in Bas-Saint-Laurent (BSL) groundwater. The grey lines represent the theoretical mixing line between fresh groundwater and seawater for Na (A), Mg (B), Ca (C), and NO₃ (D). The high-mineralized sample is reported in parentheses.

Cl and Na contents are, however, positively correlated (Figure 5A). Most of the groundwater samples have Na/Cl ratios higher than the theoretical molar Na/Cl seawater ratio of 0.85 (Quinby and Turekian 1983) and tend to reach a ratio of 1/1 or slightly higher (Figure 6A), probably induced by a combination of processes that includes cation exchange, silicate mineral weathering, and halite dissolution present in marine clays and bed rocks or evaporated seawater originating from past transgression. Whereas the Br concentrations of C4 group samples were below the limit of detection, a few samples of the C1, C2 and C3 groups had Br detected and exhibited Br/Cl ratios ranging from 0.003 to 0.0069 (Figure 6B). This range of Br/Cl ratios is smaller than that reported in the St. Lawrence Lowlands region by Blanchette et al. (2010) and Cloutier et al. (2010) and does not support halite dissolution. The seawater Br/Cl signature is clearly expressed in only a few samples of the C3 and C1 groups, with ratios close to 0.0035 (Whitemore 1988), but there is no clear relationship with proximity to the coast (except for the highly mineralized sample from C1 in which Br/Cl was ~0.0055; data not shown). Some samples of C1 and C3 groups were also characterized by SO₄/Cl mass ratios between 0.1 and 0.5, close to the theoretical marine mass ratio of 0.14 (Figure 6C). The other samples exhibited high SO₄/Cl concentrations, without a clear increase of Cl content, probably due to gypsum dissolution which is present in marine clays and Appalachian rocks (Buffin-Bélanger et al. 2015; B. Hétu, pers. comm.).

Figure 6. A, Sodium concentrations; B, bromide to chloride; and C, sulphate to chloride mass ratios reported as a function of chloride concentrations. Marine mass ratios are also reported. The high-mineralized sample is not reported.

These few samples and most samples belonging to the C3 group were collected below the upper limit of the Goldthwait Sea transgression, supporting the idea that high Cl and Na concentrations originated from evaporated or remnant seawater still recorded in marine clay and bed rock. There is limited knowledge on the impact of past climate change on modern groundwater resources. However, past marine transgressions, as far as 5 million years ago (Aquilina et al. 2015), could have introduced marine solutes throughout the aquifer domain by diffusion processes. These solutes could have induced high chloride concentrations in modern groundwater, as observed in the Armorican rock basement aquifer by Aquilina et al. (2015). Fluid signatures of basement aquifers could then have recorded hydrological changes at geological scales (Corcho Alvarado et al. 2011; Aquilina et al. 2015).

Water–rock interactions at regional scale

The geochemical evolution from Ca-HCO₃, with low-TDS and low-pH water (C4), to intermediate Ca-HCO₃ water (C2, C1), and finally to Na-HCO₃ water facies with high TDS and high pH (C3), suggests different modes of water–rock interactions. The low TDS, low pH and low TAlk of the C4 group (Table 1) reflect an early stage of groundwater recharge without much mineralization induced from water–rock interactions, while the intermediate Ca-HCO₃ water facies, with increasing Ca, Na, Mg and SO₄ concentrations and higher TDS and TAlk values, suggest dissolution of minerals as well as silicate weathering in the aquifer. In aquifers from the highlands to the coastal plain, groundwater samples evolved from Ca-HCO₃ to Na-HCO₃, with some Na-Cl facies in deep bedrock wells (see the above section), which is typical of water facies resulting from dissolution, leaching, cation-exchange reactions and mixing with mineralized end members (Bishop and Lloyd 1990; Cheung et al. 2010).

Mineral dissolution and leaching of carbonate phases (calcite, dolomite), gypsum and silicate would mobilize Ca, Mg and Na, and release HCO₃ and SO₄ in C4 and C2 groups, and to a lesser extent in C1 (Table 2). In addition, ion exchange processes where Ca and Mg in the aquifer matrix have been replaced by Na at favourable exchange sites may account for the dominance of Ca and Mg over Na for most of the groundwater samples with a Ca-HCO₃ water facies (Figure 7A). The high HCO₃ content nearly balances the sum of the major cations (Figure 8A, B) for the C4, C2 and, to a lesser extent, C1 groups, indicating a concomitant dissolution of Ca, Mg and HCO₃. The high HCO₃ could be generated by the dissolution of marine carbonate minerals such as calcite, by the dissolution of dolomite, or by the weathering of silicate minerals (e.g. dissociation of carbonic acid). δ¹³C_DIC values have been used to confirm the origin of HCO₃ concentrations in groundwater, since it reflects the isotopic fractionation between DIC species and provides relevant information about the origin and the evolution of HCO₃ in groundwater (Clark and Fritz 1997). If silicate weathering were occurring, the isotopic signature of the DIC species would be expected to be similar to the composition of the H₂CO₃ that controlled the weathering (~−22‰). The mean δ¹³C_DIC value of C2 group – and of Ca-HCO₃ water facies more generally – was close to this theoretical value (Table 1). However, the small concentrations of silica measured in BSL groundwater suggest that the dissolution of silicate minerals is not a major process.

Figure 7. Processes related to water–rock interaction reported for the four water groups: Ca (in meq/L) (A) and pH (B) in relation to Na (in meq/L). The high-mineralized sample is not reported.

Figure 8. A, Calcium content (in meq/L) in relation to bicarbonate (HCO₃ in meq/L). B, Sum of major cations (Ca + Mg in meq/L) in relation to bicarbonate (HCO₃ in meq/L). C, Sum of major ions (Ca + Mg − (HCO₃+SO₄) in meq/L) in relation to Na-Cl (in meq/L). The high-mineralized sample is not reported.

In the case of marine carbonate dissolution, HCO₃ with a δ¹³C _DIC of ~−9‰ is expected, assuming a closed system where calcite dissolution (~4.4 ± 0.5‰) and H₂CO₃ (~−22‰) occurred in equal proportions (Clark and Fritz 1997). Calcite dissolution would explain the enriched δ¹³C_DIC signature measured in at least five samples of the C1 group (data not shown). The SI calculation confirmed that the calcite is slightly undersaturated (SI < 0; Figure 9A) in groundwater samples of the C4 and C2 groups. Dissolution may occur and contribute to increases in Ca and HCO₃. The relatively depleted levels of δ¹³C_DIC (~−20‰) in most groundwater samples do not correspond to an evolution in a closed system in contact with the soil CO₂ (Clark and Fritz 1997). Despite the lack of traces suggesting significant anaerobic respiration in the BSL aquifer (i.e. high DOC content, low saturation of dissolved oxygen, and occurrence of reduced species), the influence of an organic source of DIC by microbial respiration cannot be excluded. The microbiological degradation of organic matter in the aquifer matrix could contribute to depleting the original signal.

The other source of HCO₃, Mg and Ca is the dissolution of (Mg/Ca) carbonate as dolomite. The behaviour of dolomite is quite similar to that of calcite: dolomite is undersaturated with respect to mineral phases for most of the Ca-HCO₃ water facies (Figure 9B). Favourable conditions for the formation of dolomite occur in marine-anoxic environments, as coastal sediment where organic-rich sediments are dominated by bacterial sulphate reduction (Backer and Burns 1985). In addition, the formation of dolomite may be favoured in seawater–groundwater mixtures, particularly in coastal permeable sediments of subterranean estuaries where the redox conditions support bacterial sulphate reduction (see Moore 1999). Geskeet et al. (2015) have recently shown that the isotopic dolomite signatures change as a function of their diagenetic settings and their environments of formation. Non-marine evaporative dolomite, such as lacustrine/palustrine dolomites, have distinct carbon and δ¹⁸O characteristics, with δ¹³C_DIC = −11 ± 12‰ and δ¹⁸O = −5.3 ± 7‰, while other dolomite minerals (e.g. marine unaltered and altered dolomite and hydrothermal dolomite) have enriched carbon signatures (Geskeet et al. 2015). The occurrence of the earliest diagenetic dolomites appears to be in agreement with the Quaternary setting of the studied area (i.e. occurrence of beaches, rivers and delta environments; see the section on the hydrological settings and geology). Examining all of these isotopic and geochemical constraints, it seems that neither of these mechanisms (e.g. carbonate and gypsum dissolution and silicate weathering) alone could dominate the groundwater chemistry of Ca-HCO₃ water facies; they probably act concomitantly to control the hydrogeochemistry of the dilute pristine end member (e.g. C4 and C2 groups).

Groundwater samples from the C1 group appear to be an evolution of the C4 and C2 groups and an intermediate between Ca-HCO₃ and Na-HCO₃ water facies. C1 group samples are also collected at intermediate depths, between the shallow C4/C2 and the deepest C3 groups (Table 1). The C1 group is characterized by higher Na compared to Ca and Mg and high TAlk values, as observed in samples of the C3 group, but pH is still nearly neutral, whereas it reaches values higher than 8.5 in C3 (Table 1; Figure 7B). It may result from the dilution of the mineralized C3 water group by the recharge or dilute pristine end members C2 and C4. This mix along the flow line leads to the progressive domination of Na on cation exchange sites and relatively elevated Cl and SO₄, also supported by the continuously dissolution of gypsum (SI gypsum < 0). Ca, and probably Mg, exchange for Na (Figure 7A), resulting in NaHCO₃-rich water. The Ca depletion drives incongruent carbonate dissolution where dolomite dissolution controls the precipitation/equilibrium of calcite, as shown in Figure 9B. The effect of this cation exchange induced by the mixing between evaporated and/or remnant seawater and diluted end member, besides changing the Ca, Mg and Na concentrations, is to increase the alkalinity and pH of groundwater with dolomite dissolution. The relationship between Ca + Mg-(HCO₃ + SO₄) and Na-Cl showed a negative linear trend with a slope near 1, supporting the idea that the cation exchange is coincidental with the mixing with a saline end member (Fisher and Mulligan 1997). The cation exchange seems to occur in all groundwater samples, but it predominates for samples of C3 and in some samples of C1 (Figure 8C). The ability to delineate the exact composition of the saline end member that generated the groundwater chemistry of the C3 group is limited, but the cation exchange reactions in addition to some mixing processes with evaporated or remnant seawater appear to control the mineralization of these deep groundwater samples. Further study is needed to characterize the residence time and the exact groundwater flow path in order to delineate the mechanisms responsible for the geochemical evolution of groundwater in BSL aquifers.

A conceptual model of the hydrogeochemical evolution in the BSL region along the regional flow line

Based on the distribution of major ions and on the above discussion, a conceptual model of the hydrogeochemical evolution of groundwater in the bedrock is proposed from the recharge to the discharge zone (Figure 10). This model is shown over a cross-section view, from the highlands at St. Gabriel to the St. Lawrence shore at Sainte-Luce (Figure 1B), and is based on real samples located along the target flow line.

-	Step 1: There is an infiltration of meteoric water to the bedrock through shallow till sediments and unconfined shale rock outcrops on Appalachian ridges. Congruent carbonate (dolomite and calcite) and gypsum dissolution produce Ca-HCO3 water facies (the C4 and C2 groups of groundwater).
-	Step 2: Under a semi-confined to confined setting, the groundwater moves downward to the deepest depths, where mechanisms of dissolution and, to a lesser extent, cation exchanges generate groundwater enriched in ions. Groundwater conserves the signature of the recharge area (Ca-HCO3 water facies), even though the dissolution of carbonate minerals is very active, probably due to the dilution with fresh recharge groundwater. In the valleys between Appalachian ridges, below the transgression limit, the succession of unconfined, semi-confined and confined conditions favours the generation of mixing processes between mineralized (C3) and fresh recharge groundwater (C4, C2) to form the intermediate C1 group of samples.
-	Step 3: During the seaward transport of regional groundwater, in the semi-confined to confined aquifer located below the upper limit of Goldthwaite Sea transgression, the mixing with a saline end member and the onset of ion exchange processes induce the production of Na-HCO3 water facies (C3 group) with some Na-Cl samples. The Ca depletion promoted carbonate mineral dissolution and the production of high levels of both bicarbonate and SO4, and increased the groundwater pH. The time needed for groundwater to flow to the coastal discharge zone, and the absence of a recharge zone all along the flow line, favour the evolution of highly mineralized and alkaline groundwater. Along the shore, direct salinization by modern seawater also occurs but appears to be a very local (and sporadic) mechanism.

Figure 9. Saturation indexes (SI) of carbonate minerals. A, SI of calcite in relation to pH; and B, relationship between calcite and dolomite SI. The high-mineralized sample is not reported.

Figure 10. Conceptual model of the hydrogeochemical evolution of groundwater in the Bas-Saint-Laurent (BSL) bedrock aquifers along cross section 1–2 (see location in Figure 1A and geological details in Figure 1B). Concentrations of major anions (A) and cations (B) along the section (concentrations are reported in mg/L). C, Geology, hydrological settings, and expected groundwater flow paths from the recharge area to the coastal discharge zone. The water groups are reported as a function of the ion concentrations.

Conclusions

The origin and the geochemical background of groundwater are key components to an adequate assessment of water quality at a regional scale. This study presents the first regional hydrogeochemical portrait of groundwater in the Bas-Saint-Laurent region, a region shaped by the Appalachians, a strong Quaternary glacial heritage, and coastal dynamics from the St. Lawrence Estuary. The proximity of BSL aquifers to the St. Lawrence Estuary, the geological history, and past seawater transgressions create unique issues with respect to groundwater mineralization and sustainability. The description presented here is based on an extensive groundwater sampling programme that was carried out to document the hydrochemistry and the stable isotope signatures of water and inorganic carbon.

Seven water facies, including the dominant Ca-HCO₃ and Na-HCO₃ facies (66 and 20%, respectively), were identified from the analysis of 145 water samples. No relationship among geology, aquifer confinement and geochemical facies emerged from the Piper diagram, suggesting that factors other than the hydrological setting may control the chemical composition of groundwater in the study area. Isotopic signatures from all samples suggest that infiltration to and recharge of the aquifers is contemporary to modern climatic conditions and that snow cover plays a key role in the region’s aquifer recharge, whatever the degree of confinement of the aquifers. A hierarchical cluster analysis was performed, which classified the water samples into four distinct geochemical groups showing a gradient from less mineralized in the highlands (the C4 and C2 groups, with a dominant Ca-HCO₃ facies) to more mineralized toward the coastal plain of the St. Lawrence (C1, with a Ca-HCO₃ facies, and C3, with a Na-HCO₃ facies).

Groundwater geochemical processes were investigated to determine the mechanism driving the evolution of groundwater mineralization. Analyzing the δ¹³C_DIC for the various groups led to the conclusion that calcite and more probably dolomite dissolution control the first steps of mineralization following the recharge in the highlands. The proximity to seawater from the St. Lawrence Estuary does not appear to be the main control on groundwater chemistry in the BSL aquifers. However, the Quaternary transgression events and the record of remnant or evaporated seawater in clay and basement aquifer appear to exert a major control on groundwater chemistry. The mixing between this undefined saline end member and fresh dilute pristine groundwater, in association with the onset of cation exchanges, drives the progressive replacement of Ca (and Mg) by Na in groundwater, inducing the progressive mineralization along the groundwater path.

The regional geochemical portrait allows the proposal of a conceptual model to explain the hydrogeochemical evolution of groundwater along the regional flow line. The model provides a new understanding of the origin of geochemical conditions and processes responsible for the mineralization of groundwater resources in the BSL region. This knowledge could help resource managers prevent the deterioration of the groundwater resource and anticipate future changes.

Funding

This project was funded by the Ministère du Développement durable, de l’Environnement et de la Lutte contre les changements climatiques of the Québec Government, as part of the Programme d’acquisition des connaissances sur les eaux souterraines (PACES) 2012–2015. This research was also supported by the Canada Research Chair Program (to GC), and the Université du Québec à Rimouski (GC, TBB).

Acknowledgements

The authors thank numerous field assistants for their support, the collaboration of the local population and municipal authorities, and Laure Devine for the English revision. The authors would like to thank Marie Larocque and Vincent Cloutier for the invitation to submit this contribution to the PACES special issue, and the two anonymous reviewers for their helpful comments that contributed to the improvement of the manuscript.