Remediation scenario of the abandoned Kettara mine site (Morocco): acid mine drainage (AMD) transport modeling

In this study, we present the efficiency of remediation scenario to attenuate the impact of acid mine drainage (AMD) contamination in the Kettara abandoned mine site. The study focuses on the AMD groundwater contamination of the Sarhlef shists aquifer. To predict the evolution of AMD groundwater contamination in the Kettara mine site under remediation scenario, a model of groundwater flow and AMD transport was performed. Piezometric heads were measured at the dry and wet periods from eleven wells located downstream of mine wastes. To elaborate a conceptual groundwater flow model, we faced with to the heterogeneity and anisotropy of fractured Sarhlef shists aquifer. Consequently, the study focused on the use of various approaches: 1. The inverse modeling by the CMA-ES algorithm is adopted as an alternative approach to determine hydraulic parameters indirectly, and 2. the model is treated as an equivalent porous media (EPM). The groundwater flow model was carried out in steady-state and transient conditions in the dry and wet periods using the PMWIN interface. The obtained results are satisfactory and show an excellent correlation between measured and computed heads. Contaminant transport model is used to solve the advection–dispersion equation and to generate the AMD concentration by MT3D via the PMWIN interface. A sensitivity analysis of the dispersivity coefficient is carried out. The AMD transport simulation was computed during periods of 1, 5 and 10 years, and the performed model indicates that the simulated concentrations under remediation scenario are reduced 1000 times comparing to the current concentrations. The study revealed a necessary approach in addressing an environmental issue for the AMD contamination. The results of the study will be a start-up for further research work in the study area and implementing it for the prevention of AMD propagation plume.


Introduction
Groundwater constitutes a vital water resource in Morocco which is characterized by an arid climate with limited rainfall and an important rate of evaporation. In abandoned Kettara mine site, the groundwater is becoming undrinkable due to high sulfate concentration that exceeds 1800 mg/l [31]. The Kettara mine site concept remediation will consist of collecting and placing coarse tailings over tailing pond. The process requires placing a fine alkaline phosphate waste (APW) layer on the Kettara coarse tailing (Fig. 1e) [20]. The APW acts as a capillary barrier which must block the access of water to the mine wastes, and consequently, the generation of AMD will be limited. Remediation of abandoned mine areas is becoming a more significant environmental issue in the world [4,23,37,47,55].
The research work aims to evaluate the efficiency of remediation scenario at the Kettara mine site by predicting the evolution of AMD contamination in groundwater. By comparing the current AMD groundwater contamination to the predicted one under remediation scenario, we can quantify the rate of removal of contamination and the efficiency of remediation concept. Several studies on the prevention of environmental contamination related to mine waste have been carried out [2,9,17,26,33,37,38,45,[57][58][59].
The present work contains two steps: The first one is a conceptual groundwater flow model and the second is to simulate the AMD contamination. For elaborating the groundwater flow model, the primary purpose is how to describe the heterogeneity associated with fractures of the Sarhlef schists aquifer. Modeling processes in hard rocks and their associated fractured aquifers has been of high importance last decades [8,16,19,28,29,35,39,52,61]. One of the simplest approaches is basing on the concept of a representative elementary volume. They assume that the aquifers including fractures and conduit networks can be represented by an equivalent porous medium (EPM) with equivalent hydraulic conductivity. The EPM approach is commonly used by several studies for modeling groundwater flow and transport in heterogeneous aquifers [1,14,27,44]. In our present study, we have opted for the equivalent porous medium (EPM) approach to represent the Sarhlef schists aquifer.
Due to the unknown information of hydraulic parameters and limited financial resources for their measurement in situ, the determination of aquifer transmissivity from measured heads constitutes the second purpose. In heterogeneous aquifers, the groundwater flow model has been calibrated by the inverse approach. Inverse approach modeling is becoming a powerful method to approximate and attribute relevant values to unknown features in hydrological models [5,24,40,49,53,56]. Inverse approach combines two compounds: hydrodynamic model code and optimization algorithm. The CMA-ES algorithm, one of available the optimization algorithms, has been used in this research. The EPM approach and the CMA-ES algorithm are used to produce a model that accurately represents the real system of groundwater flow in Sarhlef schists aquifer. Processing Modflow (PMWIN) [7] is used to simulate the groundwater flow, and the results obtained are combined with MT3DMS mass transport model to predict the AMD transport in groundwater of the Kettara mine site under remediation scenario. The obtained results are promising and they show clearly the importance of remediation scenario to attenuate the AMD contamination. The remediation project of the abandoned Kettara mine site, under study, is the first in Morocco. Once tested and approved, the findings will be applicable to other mines located in semi-arid climate.

Study area
The extraction of pyrrhotite (FeS) in the Kettara mine is carried out between 1965 and 1982 and was mainly aimed at production of sulfuric acid. It has been estimated that during this mining activity, more than three million tons of sulfureous wastes were stored over an area of approximately 16 ha [20]. After the closure of the mine, these wastes are still without remediation and constitute the biggest pollution issue in the Kettara region. The abandoned Kettara mine is located at 30 km north-northwest of Marrakech, at the edge of the road connecting Marrakech to Safi (Fig. 1a). The climate of Kettara mine site is semiarid with an annual potential evaporation of 2500 and 250 mm of annual rainfall.

Geological and geophysical setting
The Kettara sulfide deposit is located in the Sarhlef series; this series belongs to the Hyrcynian Jebilet massif affected by the post-vesean metamorphism [25]. Metamorphic deformation is distinctly observed at the Kettara mine site (Fig. 2a). The microtectonic structure is the schistosity with an average orientation of N45° (Fig. 2c). Faults with quartzo-carbonated seams are the tectonic structure observed with principal direction of N75°, N95° and N110° (Fig. 2b).
An electrical resistivity tomography (ERT) profile situated 2 km downstream to the Kettara mine site was carried out [30]. The ERT shows the internal geological structure (Fig. 1d). The first layer with a low resistivity corresponds to Quaternary alluvium. The second layer with an irregular morphology corresponds to altered schists. In the last layer, we can observe the presence of very heterogeneous formations with lateral and vertical variations in resistivity. These variations make it possible to identify moderately the schist's formation as to highly resistant areas interspersed with more conductive areas (faults).

Hydrogeological setting
At the Jebilet massif, two major aquifer systems are identified [11]: a superficial aquifer located in the altered schists and granites, and a deep aquifer of discontinuous water flows at the level of faults and fractures of the crystalline basement (schists and granites). The majority of wells surveyed the groundwater located in the altered and fractured schists formation, in either quaternary alluvium along the Oued Kettara or other thalwegs. Other wells and boreholes are exploiting the deep aquifer sheltered in the schist's substratum.
The piezometric level of the altered schists aquifer is about 15 m, and the hydrodynamic parameters are estimated [11] with a transmissivity of 9.10 −4 m 2 / s and a coefficient storage of 5.10 −2 . This aquifer is subject to human exploitation by traditional wells. This recharge is provided by direct infiltration of meteoric water through fractures and permeable alluvium and by the Kettara Wadi.

Data collection
Eleven wells, located in the study area, were sampled and analyzed for physicochemical characterization [31]. These wells are located downstream of the Kettara mine and follow the general direction of the groundwater flow in this area. Sampling of well water and measuring of groundwater head are carried out during two campaigns covering dry and wet periods (March 2011 and June 2012).

Equivalent porous medium (EPM)
According to geological and hydrogeological studies, the Sarhlef schists aquifer is complex and highly heterogeneous. The presence of faults and fractures are probably responsible for the heterogeneity of this aquifer. Indeed, the notion of fractured media is based on the existence of cracks and/or faults influencing the fluids flow through these media. A fractured porous medium is imagined as an interconnected system of cracks dividing the medium into a series of porous blocks, called "porous matrices." The flow characteristics of a fractured medium depend on the degree of fracturing, the connectivity of the fracture, and the variation of porosity and permeability parameters [41]. The groundwater modeling in this medium requires to adapt a conceptual model. According to [47], there are several conceptual models: stochastic continuum (SC) model, triple porosity | https://doi.org/10.1007/s42452-021-04690-6 medium, equivalent porous medium (EPM), discrete fracture network (DFN). The EPM model consists of replacing the discontinuous values of porosity and permeability by equivalent mean values [14]. These mean values can be obtained by a homogenization procedure. This EPM model treats the fractured medium as a homogeneous medium with average hydraulic properties [48]. Various studies in fractured media have adopted the equivalent porous medium approach [1,8,14,16,27,29,47]. Owing to the heterogeneity of Sarhlef schists aquifer and the unknown mode of groundwater flow, we decided to adopt an equivalent porous medium (EPM) as a modeling approach. The objective is to reproduce the observed piezometric by calibrating equivalent hydraulic conductivity. Our model will be assimilated as a confined aquifer with a single layer (single-layer model).

Inverse approach modeling
Owing to the lack of hydrodynamic parameters data necessary for the characterization of the Sarhlef schists aquifer, the calibration of the model cannot be performed. Indeed, the capacity of a model to simulate the observed measurements depends on the satisfactory calibration (of one or more selected parameters). To obtain a deeper insight of aquifer features and improve groundwater modeling accuracy, the estimation of hydrodynamic parameters is necessary. In our study, we first tried to calibrate the model basing on the transmissivity values extracted from the study carried out by [11]. Due to the insufficient data, the calibration result after several tests was very poor. The foundations to carry out a conceptual model is in collecting the information and data [5]. In heterogeneous aquifers, the groundwater flow model has been calibrated by the inverse approach to obtain the optimal fits between the simulated and observed groundwater heads. Heterogeneity and limited features data about schists aquifer led us to adopt an inverse approach modeling. Indeed, the advantage of inverse approaches over direct approaches is that the formulation of the inverse problem applies to situations where the environment is heterogeneous and where the observations are few and poorly distributed [27]. Many issues of groundwater model flow have used an inverse approach in the last decades [5,49], and these methods generally lead to better solutions. The fundamental concept of inverse strategy is the optimization in order to minimize an objective function. When the objective function or the best fit in the least-squares sense minimizes the difference between the simulated value and the observed one, the parameter optimization is achieved.
An inverse method that combines the hydrodynamic model code and the optimization algorithm is followed for the identification of hydrodynamic parameters (transmissivity). A large number of optimization codes are used: Newton methods, artificial neural network (ANN), evolution strategies (ESs) and genetic algorithms (GAs) [49]. The covariance matrix adaptation evolution strategy (CMA-ES) algorithm, one of the most efficient optimization algorithms, is used as an identification approach in the present study [3,21,49]. We used the piezometric levels of the high-water period as input data. The boundary conditions correspond to the measured hydraulic heads imposed at each node. The average mesh size used is 4997.38 m 2 . The finite element method was selected to calculate the piezometric level. As shown in Fig. 3, the comparison between measured and calculated heads of the two maps is similar, so it can say that the CMA-ES algorithm reproduces the measuring head satisfactorily.

Groundwater flow modeling
The aquifer has been carried out as confined with a single hydrogeological layer, and the modeled domain covers an area of 16 km 2 (Fig. 4a). The model grid was discretized into 45 columns and 47 rows. A sensitivity study of the different cell sizes: 20, 60, 80, 100, 120, 140 and 160 m, was carried out. The top elevation of the aquifer was obtained from the DEM (digital elevation The volcanic outcrops and the basic tuffites constitute the lateral geological limits of the aquifer and are defined as no-flow boundaries. The limits with constant head boundaries correspond to the schistous outcrops and are localized in the northeast and the south of the study zone. The limit with fixed head boundaries corresponds to the Kettara wadi. This limit was taken in only for simulation in the wet period because the flow is temporary and taking place during this period. Groundwater flow simulation in Sharlef schists aquifer was carried out by the Processing Modflow for Windows code (PMWIN, version 5.3) [7]. This code is widely used in hydrogeology [2,15,26,33,43,47]. The finite difference is the method used by PMWIN code for solving the three-  The model of flow simulations was performed under steady-state and transient conditions. Water entries and exits were added to the model (Fig. 4c). The inflow comprises the recharge from precipitation (4,16 × 10 −4 m/j). We suppose that the recharge is homogeneous in all the studied zone. The water exits correspond to the pumping wells discharge (Q = 0,001 m 3 /S).
The main goal of calibration is to reproduce the field measured heads by an iterative process. Calibration results were assessed using the root mean square error (RMSE) of the hydraulic head, and a good fit indicates that the model is calibrated and therefore potentially representative of the aquifer flow [53]: where h sim : hydraulic simulated head, h obs : hydraulic measured head and n: number of calibration hydraulic heads used in error computations.
The calibration is an iterative process conducted in steady-state and transient-state simulation.
To reproduce the observed head in steady-state and transient-state simulations, the groundwater flow model was calibrated using measured head data from 11 wells during March 2011 (wet period) and during June 2012 (dry period). For this, the coordinates (X, Y), as well as the corresponding measured heads, of these wells have been integrated into the Boreholes and Observation module of PMWIN. The Kettara wadi works as drains during highwater periods and does not affect during dry periods.

Contaminant transport simulation
As an input to contaminant transport models, groundwater velocities are computed. Contaminant transport models are used to solve the advection-dispersion equation and to generate the contaminant concentration by MT3D via the PMWIN interface [7]. The equation of advection-dispersion contaminant transport in three dimensions is: where U: Darcy velocity and C: concentration.
The retardation factor R of the compound dissolved in water can be calculated from the estimation of the soil-water partition coefficient Kd of this compound [L 3 M −1 ], the effective porosity and the bulk density of the porous medium [ML −3 ] (Eq. 4). In our case, conservative solute, the retardation factor is equal to 1.
The solute is transported in the direction of groundwater velocities by the process of advection and convection.
Molecular diffusion results from molecular agitation of solute molecules. The result of this molecular agitation is a transfer flux of solute particles from high-concentration zones to low-concentration zones. For a fluid that circulates in a porous medium, the phenomena of advection and diffusion are easily combined by establishing again the conservation of the mass of the element transported in an elementary volume, by summing the two flows of matter.
The dispersion is linked to the heterogeneity of the porous medium on a small and large scale, and it is at the origin of the "spreading" of plume pollution and contributes to dilute the concentrations. There The dispersion coefficient is related to the effective velocity U. The spread of solutes is called longitudinal dispersion D L when it is in the direction of the main flow and transverse dispersion D T in perpendicular directions. In general, the longitudinal dispersion is much greater than the transverse dispersion. The ratio (α L / α T ) provides the form of the plume: The smaller this ratio, the larger the rising plume.

Inverse modeling approach
The transmissivity values repartition is heterogeneous (Fig. 5), the highest values (> 0.02 m 2 /s) are found in particular in the southeast and southwest parts and the values of 0.002 to 0.01 m 2 /s are distributed over the rest of the study site. A study carried out by [50] and using the CMA-ES code shows that the calculated transmissivity has a heterogenous distribution over the studied area. The transmissivity distribution was being explained by the aquifer lithology and thickness [50]. In our study, this distribution is related to the geometrical characteristic of the Sarhlef shists aquifer, especially owing to the important thickness and to the faults network in these areas. The highest values of transmissivity can be assigned to the important thickness of the Sarhlef shists aquifer, and the presence of faults affects transmissivity distribution by dividing the aquifer to hydrogeological blocks. In the absence of pumping well tests, the obtained results will be constituting a reference for future hydrogeological studies in the region.

Sensitivity analysis
To test the influence of the mesh size on the simulation modeling, a comparison between simulated and measured hydraulic head over eight different grid size: 20 × 20 m, 40 × 40 m, 60 × 60 m, 80 × 80 m, 100 × 100 m, 120 × 120 m, 140 × 140 m and 160 m × 160 m was made to establish the proposed grid size according to scatter plots (Fig. 6). The results show that the correlation coefficients between the measured values and the predicted values are between 47.5% and 84.7%. Thus, these findings show that the ideal grid size was chosen based on the highest correlation coefficient. Results propose a 100 m × 100 m grid size for the modeling process in the studied aquifer. Consequently, the model will be composed of 211,500 meshes of 100 × 100 m size. Several studies were adapting a grid refinement of contaminated area in numerical groundwater model [12,34]. The groundwater solute transport models may be constructed at a smaller scale with finer discretization than the flow models in order to accurately delineate the solute source and the modeled target [12]. A refinement of the grid with a mesh size of 50 × 50 m is located near the source of contamination (mine tailings).

Calibration and verification of groundwater flow model
By comparing the hydraulic heads obtained from observation wells (n = 11) and the computed hydraulic heads, steady-state calibration of model flow was achieved. Another way to display the calibration fit is the scatter plot. Furthermore, the calibrated model outputs were evaluated by the root mean square error (RMSE). The simulated piezometric map for the period of high water (March 2012) [31] is presented in Fig. 7a. On the whole of the modeled domain, we reproduce very well the general appearance of the piezometric map. It can be seen that most of the observation points are located around 45° straight line, illustrating an RMS value R 2 = 0.847. The difference between the calculated values and the measured values does not exceed 10% of the piezometric amplitude of the modeled area, which corresponds to a satisfactory calibration. The wells for which the deviation is large are P1, P2, P7, P10, P11 and P12. These deviations can be attributed, on the one hand, to the system for measuring the absolute altitudes of the wells and, on the other hand, to their situation in the zone of low transmissivity value.
The simulation was made based on the piezometric of June 2012 (Fig. 7c). In detail, when we compare the calculated and observed values, we see that the differences are generally acceptable (R 2 = 0.592), and most of them are less than 10 m. The piezometers for which the deviations are the greatest are all located at the extreme downstream of tailing mine (P9, P10, P11 and P12). Indeed, at the level of these wells, agricultural farms solicit and use water in a significant way.
During the dry period, most of the wells present a drop in the piezometric level and become exhausted under the effect of discharges. Therefore, the simulation is generated only during the high-water period. The transient  (Fig. 7b) and taking one year us a period of simulation. The setting is satisfactory and shows an excellent correlation between measured and computed heads (R 2 = 0.995). The well with the most significant deviation is P11, related to his localization in a low transmissivity value zone.
In order to confirm the model validation, calculation of water balance is necessary, and the finding is presented in Table 1. The analysis of these balances makes it possible to verify that the water inputs and outputs correspond well to the supposed behaving of the aquifer system. The balance of flows in wet period (steady-state and transient regime) and dry period (steady-state regime) is well balanced. The inputs represent the supply of the aquifer, and the outputs represent the drainage. We deduce that the water inputs in this aquifer derive from recharge by precipitations. Despite the limited rainfall in Kettara mine site, and in arid and semiarid areas, it can be concerned as the primary hydrogeological input. [60] assess very low groundwater recharge that may be less than 20 mm/a in semiarid and arid regions. The water outputs from the aquifer are due to i) the Kettara wadi and ii) the pumping well.

Acid mine drainage transport simulation
The groundwater contamination in the Kettara mine site is mainly related to AMD. Sulfates (SO 4 2− ), characteristic of AMD, constitute the element with the highest concentrations in water and therefore the main source of pollution of these waters [31]. For local scale transport simulations of AMD contamination in schists aquifer, it is necessary to describe the sulfates distributions within the source. The current situation of groundwater Kettara pollution by sulfates (SO 4 2− ) is shown in Fig. 8. Furthermore, the degree of contamination of groundwater depends on the interactions with the soil, sediment and rocks that constitute the path of the contaminant.
To determine the evolution and the future of AMD in the Sarhlef schist aquifer, we called up by the pollutant transport model. The modeling of pollutant transport is a tool to help predict risk.
As part of this modeling, the boundary conditions corresponding to a constant sulfate concentration of 200 g/ m 3 were imposed at the mine tailings site. This value corresponds to the standard sulfate concentration in groundwater [55]. It has been used as an initial concentration in the remediation scenario. For all the other cells of the model, the initial sulfate concentration is zero.
Unfortunately, no tracking tests have been carried out on the site. For this reason, dispersion coefficient values were estimated according to a commonly used method based on the scale of the study area (Eq. 7). This method consists in considering, as a first approach, a longitudinal dispersion coefficient of the order of 0.1 times the length of the pollutant plume or the length traveled by the pollutant [49], that is: α L : longitudinal dispersion coefficient (m) and Lp: length of the plume (m).
The dispersivity value obtained is only an estimate, which can give, at best, a large range of variation in dispersivity for a given plume length. However, this relationship represents a good starting point for modeling work. Besides, there are many empirical relationships between the longitudinal, transverse and vertical dispersion coefficients that have been described. Generally, α T is estimated at 0.1α L [12] and α V is estimated at 0.05 α L .
There is no information available on the value of the molecular diffusion coefficient, but in most cases, it is small and can be overlooked compared to other phenomena such as kinematic dispersion.
To estimate the impact of the variation in dispersivity values on the obtained results, a sensitivity analysis was carried out by varying the values of the dispersion coefficients according to four scenarios ( Table 2). In addition, the estimates of this parameter corroborate by [12] study which relates dispersion coefficient to the length of plume, obtained at porous and fractured aquifers, thus enabling it to apply to the Sarhlef schist aquifer. The predictive simulation times are of the order of one year, five years and ten years.
The comparison of the obtained maps after simulation makes it possible to estimate the sensitivity of the model to the variation of each of the above three parameters (α L , α T and α V ). Values assigned to the parameters for each simulation are shown in Figs. 9 and 10. Overall, the simulated periods (one and ten years), in scenarios 1 and 4, show a difference in terms of the shape and extent of the plume. The obtained results for the first scenario show that the pollutant transported mainly along the direction of the groundwater flow, the contours are moving away from the source (mine tailings) and the plume is getting more wider and longer at 10 years of simulation time. For scenario 4, we can show similar trend, with less wide and (7) L = 0.11Lp longer plume. This can be well explained by the values of the longitudinal and transversal dispersion coefficients which present the two borderline cases (maximum for scenario 1 and minimum for scenario 4). The concentration gradient is decreasing for both scenarios below the remediation assigned value. Scenarios 2 and 3 have the same overall plume shape over the three simulation periods. This plume grows over time and progresses in the direction of flow while moving downstream of the tailings. The results show that the longitudinal and transversal extension of the plume into sulfates downstream is greater in the case of a high dispersion coefficient and over a large simulation period. Investigations carried out by [6] confirm that dispersivity increases with distance or observation scale. Dispersion parameters, therefore, play an important role in the evolution of the contaminant plume.
The concentration versus time graphs is shown in Fig. 11. The concentrations at time level t = 0 (s) correspond to the current state. From t = 3.15 1007 s (time of the start of simulation), we note that the concentrations drop below the assigned remediation concentration value for all wells and increase slightly until they stabilize beyond five years. However, the wells of lower concentrations are P5, P6 and P7 located at the SW of mine tailings. We also note that for the low dispersion coefficients these concentrations are slightly higher (scenarios 3 and 4). The wells with the highest concentrations are those located downstream of the mine tailings.
The evolution of concentrations relative to the center of the plume and along the direction of flow for the four scenarios shows a gradual decrease while going downstream of the mine wastes (Fig. 12). The gradient is slightly higher for the minimum dispersion values. This evolution is mainly due to the effect of dispersion and advection.
In general, the transport of AMD depends on the flow conditions imposed on the hydrodynamic flow model. It is simulated during high-water period, and the aquifer is assuming to be equivalent and continuous-our explications was based on this simple approach to predict the evolution of AMD pollution under remediation scenario. However, the aquifer is very heterogeneous and fractured, so the results of the simulation of pollutant transport must be explained, taking into account the reel flow conditions. It would, therefore, be very interesting to carry out measurements in the field, of the tracing test type, to better constrain the model, and thus improve its predictive qualities.

Conclusion
In the current study, modeling approach was used to reproduce the flow and pollutant transport of an heterogenous aquifer by accordingly adopting an inverse approach to estimate the hydrodynamic parameters (transmissivity).
Results of numerical flow modeling are promising; however, it is recommended that the current groundwater flow model can in no way be considered as being definitively developed. In perspective, to improve and adjust the hydrodynamic model, pumping tests are necessary to get the hydrodynamic parameters of the aquifer. The simulation of contaminant transport shows that the pollutant concentrations are significantly reduced and the groundwater quality is improving. The obtained results indicate the future groundwater impact and constitute a valuable contribution of groundwater assessment under remediation scenario in the study area.
The remediation project will, therefore, create an appropriate solution that will make it possible to limit the generation of pollution by AMD and thus protect the groundwater against contamination. It is farther necessary to remove the uncertainties and to propose a definitive predictive model by setting up: 1-a system of piezometers for the Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.