Lifetimes of Used Nuclear Fuel Containers Affected by Sulphate-Reducing Bacteria Reactions inside the Canadian Deep Geological Repository

: The present work aims at approximating the reduction of sulphate to sulphide caused by sulphate-reducing bacteria (SRB) inside the Canadian deep geological repository in order to calculate the expected lifetime of used nuclear fuel containers (UFCs). Previous studies have assumed a conservative constant concentration of sulphide at the host rock interface. The novelty of this study resides in the use of ﬁrst-order kinetics to explicitly account for the SRB-induced sulphide production. This reaction term is developed following an empirical approach using published results on actual sulphate reduction by SRB and included in a coupled reaction-diffusion system. Lifetimes of UFCs are subsequently calculated following the conditions of two scenarios: having SRB active only at the region closest to the host rock and having SRB active at the host rock and throughout the bentonite clay. This study shows that the mean lifetimes of UFCs in both cases are above one million years. However, more accurate results would require the characterization of the host rock and groundwater of the prospective emplacement, as well as additional experiments on growth and sulphide by the microbial communities from the site.


Introduction
Deep geological repositories (DGR) are considered the safest long-term nuclear waste disposal solution [1]. They involve the excavation of an underground site where used fuel containers (UFC) are placed for permanent storage. The basic idea of the design is to use several barriers of different materials to protect and contain the UFCs. A number of countries are currently developing DGRs, among them are Finland, Sweden, and Canada [2][3][4].
The Canadian DGR design is planned to host approximately 4.6 million CANDU (Canada Deuterium Uranium) used fuel bundles contained within 96,000 UFCs in a suitable rock formation, either crystalline or sedimentary. The vaults containing the UFCs will be located at a depth of about 500 m underground [5]. The design is based on a multiplebarrier system which is designed to contain and isolate the nuclear waste. Two important layers of protection in this system are provided by the UFCs and the highly compacted bentonite (HCB). The UFCs are 2.5 × 0.556 × 0.556 [m] canisters made of steel with a 47-mm-thick cylindrical shell with 30-mm-thick hemispherical caps and an external 3mm-thick copper coating [5]. The HCB is an absorbent aluminum phyllosilicate clay that has been compacted to achieve a designed dry density and is the material responsible for isolating the UFCs from the host rock of the vault and preventing any inflow or outflow of materials. Bentonite swells upon contact with water which makes it a natural water barrier and a sealer to the placement rooms. Buffer boxes made of HCB are created to enclose the UFCs, these boxes are then placed inside the placement rooms. The space between the the UFCs, these boxes are then placed inside the placement rooms. The space between the buffer boxes and the host rock is backfilled with bentonite pellets dubbed 'gapfill' in this work ( Figure 1).

Modelling Scope and Container Failure
Mass transport inside the DGR is an important process closely connected to the safety of its design, whether the inward transport of corrosive species [6,7] or the outward transport of radionuclides [8,9], the capacity of the engineered barriers to block the mobility of chemical species is tightly connected to the long term safety of the DGR. In the Canadian DGR design, the inward transport of sulphide has assumed a constant concentration of sulphide at the rock interface [6]. A value of 3 ppm (9.07 × 10 −5 M) has been deemed conservative as the stringent conditions of high pressure and low available water would limit bacterial growth [6]. However, so far, this constraining process has not been explicitly included in the mass transport models. Having a source term would provide information about the transient stage inside the placement room and be a step towards a more realistic and comprehensive model.
Three-dimensional models of the placement rooms show that the flux of sulphide from the host rock to the UFCs is higher at the hemispherical ends [6]. Therefore, in the present work, a one-dimensional model of the section that goes from the host rock to the endcap of the UFC is developed.
Container failure is defined as the through-wall penetration of the outer protective layer [10], i.e., the 3-mm-thick copper layer. Other definitions of failure exist, such as plastic deformation, brittle, or ductile fracture; however, since the focus of this study is the mass transport of corrosive species, the former definition is adopted.

Reaction-Diffusion System
The system proposed to describe the production of sulphide by the SRB and its dispersion over a one-dimensional space is a set of coupled reaction-diffusion equations: Initial and boundary conditions are given and presented in the next section; S stands for the concentration of sulphate, , and H for that of sulphide ; the reaction

Modelling Scope and Container Failure
Mass transport inside the DGR is an important process closely connected to the safety of its design, whether the inward transport of corrosive species [6,7] or the outward transport of radionuclides [8,9], the capacity of the engineered barriers to block the mobility of chemical species is tightly connected to the long term safety of the DGR. In the Canadian DGR design, the inward transport of sulphide has assumed a constant concentration of sulphide at the rock interface [6]. A value of 3 ppm (9.07 × 10 −5 M) has been deemed conservative as the stringent conditions of high pressure and low available water would limit bacterial growth [6]. However, so far, this constraining process has not been explicitly included in the mass transport models. Having a source term would provide information about the transient stage inside the placement room and be a step towards a more realistic and comprehensive model.
Three-dimensional models of the placement rooms show that the flux of sulphide from the host rock to the UFCs is higher at the hemispherical ends [6]. Therefore, in the present work, a one-dimensional model of the section that goes from the host rock to the endcap of the UFC is developed.
Container failure is defined as the through-wall penetration of the outer protective layer [10], i.e., the 3-mm-thick copper layer. Other definitions of failure exist, such as plastic deformation, brittle, or ductile fracture; however, since the focus of this study is the mass transport of corrosive species, the former definition is adopted.

Reaction-Diffusion System
The system proposed to describe the production of sulphide by the SRB and its dispersion over a one-dimensional space is a set of coupled reaction-diffusion equations: Initial and boundary conditions are given and presented in the next section; S stands for the concentration of sulphate, SO 2− 4 , and H for that of sulphide [HS − ]; the reaction or source term R k ρ , S is a function of the bentonite dry density and the concentration of sulphate; and the diffusion coefficients are also a function of the dry density of the bentonite. The reaction term is assumed to be a first-order chemical reaction that approximates the reduction of sulphate to sulphide by the SRB and whose reaction rate constant is obtained from an optimization procedure applied to published experimental results as shown in Section 3.
The reaction-diffusion system of Equation (1) is implemented and solved using a one-dimensional finite-difference model in MATLAB ® . Since the model is intended to produce several thousand runs, an experiment is developed to find a balance between accuracy and computation time. The accuracy of the solution is mainly determined by the spatial increment ∆x; therefore, the experiment solves the reaction-diffusion model for different spatial increments and the criterion for selecting the optimal size of ∆x is the minimal computation time such that its accuracy does not deviate more than 0.5% with respect to the most accurate lifetime solution (see Appendix A).

Initial and Boundary Conditions
A sharp hydraulic gradient would exist between the pressure of the geosphere and the atmospheric pressure of the placement rooms of the DGR [11]. Groundwater would seep into a repository and saturate the gapfill and the HCB. Therefore, the initial and boundary conditions of the model are given by the composition and concentration of the chemical species in the groundwater of the host rock. Currently, no final location for the DGR has been selected and, therefore, no sample analysis has been performed. There exists, however, reference groundwater concentrations for crystalline and sedimentary rock at approximately 500 m underground [12]. The concentrations of sulphate and sulphide are taken as input values to calculate the boundary conditions and are shown in Table 1. Initial conditions for both chemicals are assumed to be zero, other saturation conditions were tested but they do not significantly affect the sulphide flux to the UFC or their lifetimes. The left-hand side (LHS) boundary condition for sulphate is kept constant at the given reference concentration of Table 1 and its right-hand side (RHS) boundary is given by the Neumann boundary condition: i.e., the RHS acts as a physical barrier. For sulphide, it is assumed that the host rock would act as a physical barrier preventing any sulphide from leaving the system and, therefore, the LHS boundary is defined using the Neumann boundary condition The RHS boundary for sulphide is kept constant at a concentration of zero, following the conservative assumptions in [6,13] that copper corrosion is mass-transport limited and that any sulphide molecules reaching the surface of the UFC will instantly corrode the copper layer.

Effective Diffusivity
The effective diffusivity of sulphide anions is taken from [14], which made the effective diffusivity depend on the bentonite dry density as follows D e = 5.30087 × 10 −10 exp −2.561 × 10 −3 ρ Dry (2) were the effective diffusivity (D e ) is in units [m 2 /s] and the dry density of the bentonite (ρ Dry ) has units [kg/m 3 ]. For the effective diffusivity of sulphate, this work follows [15] using half the value of the sulphide diffusivity. Table 2 shows the properties and diffusivities for the designed densities of the bentonite clay in the Canadian DGR.

Sulphate-Reducing Bacteria Activity
The chemical reaction describing the sulphate-reducing mechanism [16] carried out by sulphate-reducing bacteria (SRB) is the following This chemical reaction is likely non-elementary; however, under the assumption of an abundance of hydrogen gas (H 2 ) and hydrogen ions (H + ), it can be approximated by an elementary first-order reaction. In [17] is shown that the reduction of sulphide into sulphide can indeed be approximated by first-order kinetics.
The procedure adopted in this work is to use a two-step process to model SRB activity. First, a first-order kinetics is found from experimental results about sulphate reduction under favorable laboratory conditions using SRB microorganisms extracted from a deep bedrock [18,19]. Then, the kinetic constant is adjusted to match experimental results made under conditions similar to those present inside the DGR [20].

SRB-Driven Sulphate Reduction
Data was used from [18,19], who performed an experiment to study copper corrosion on the presence and absence of SRB. Copper specimens were exposed to biotic and abiotic artificial anoxic groundwaters for periods of 4 and 10 months and then examined their chemical and electrochemical changes. The microorganisms used in the experiment came from the deep bedrock of Olkiluoto, Finland.
A particular featured of the results is that pH decreased after the experiment, meaning that the concentration of hydrogen ions increased, since H + = 10 −pH (from pH = −log 10 [H + ]). A decrease in hydrogen ions was expected in reaction (3) under controlled conditions, indicating that likely other reactions took place which increased the hydrogen ion concentration.
Therefore, since hydrogen ions were not limiting the reaction, it was assumed that hydrogen ions would be in abundance for groundwater with similar characteristics. This is a conservative, bounding scenario as the real rate of sulphide production will be significantly less than this. Additionally, the electrolysis of water into hydrogen gas (H 2 ) and oxygen (O 2 ) was assumed to be faster than the sulphate reduction to represent a conservative scenario; therefore, hydrogen gas was also considered in abundance inside the DGR. However, the use of HCB would prevent the formation of a biofilm on the surface of the UFCs as it did in the experiments; therefore, this sulphide production rate is an extreme bounding case representing the worst-case scenario in which every engineered barrier system fails [20] The conservative bounding upper limit of the reaction rate is described by where k ρ is the reaction rate coefficient as a function of the dry density ρ dry , χ is the efficiency of the reaction, and S and H are concentrations in [mol/L]. If concentrations were in [mg/L], then the right-hand side of Equation (4b) must be multiplied by the factor 33.07/96.06 to account for the mole conversion.

Reaction Rate Constant
The analytical solution of (4a) describes an exponential decay: where a 0 can be determined from the initial conditions and for k ρ an optimization procedure was used to find the optimal value of k ρ * that minimizes the error of the predicted reduction with respect to the observed data. Two objective functions are tested, the squared differences and the square of log differences.
The optimization model (6) is minimized using the bracketing method, which is an iterative approach and where d i is the concentration of SO 2− 4 after time t i obtained from the experiment. This model yields an optimal kinetic constant of k ρ * = 5.9631 × 10 −9 [1/s]. The second optimization model: has a closed optimal solution given by and the optimal kinetic constant in this model is k ρ * = 6.1614 × 10 −9 [1/s]. This rate constant is adopted since it refers to the most unsafe situation and because it yielded a smaller error in two out of three error measures (see Appendix A).
Therefore, the optimal reduction of sulphate describing the results from experiments in [18,19] is

Efficiency of the Reaction
Using the same experimental results from [18,19] and reaction (3), the theoretical sulphide production from the stoichiometry is calculated and compared with the actual sulphide production to determine the efficiency of the reaction. The theoretical production is 34.6 [mg/L] and the actual production is approximately 14.14 [mg/L]; therefore, the efficiency of the kinetics is

Suppressing Effect of Bentonite Dry Density
Since the extreme pressure conditions inside the DGR differ from those on the experiments [18,19], the reaction rate constant of Equation (8) is adjusted using the dry density of the bentonite as a suppressive factor because it determines the pore size of the clay and the swelling pressure which are two major factors limiting SRB activity.
An exponential functional form is selected to attenuate the reaction rate constant as experiments have shown that the swelling pressure of bentonite increases exponentially with respect to the dry density of the bentonite clay [11]. The exponential form developed must meet the following conditions , the reaction constant under the baseline density of the experiment of Section 3.2 must equal the optimal reaction constant of Equation (8).
, if the dry density equals the specific density of the bentonite clay ρ S = 2760 kg/m 3 , then the reaction constant becomes zero as there is no pore space available for bacteria to remain active.
The exponential form developed is where k A ρ Dry is the attenuated or adjusted rate, ρ Dry is the dry density of the clay in [kg/m 3 ] and C f is a value at which the rate constant is taken to be zero (details on the derivation of this formula are presented in [21]). Sulphate to sulphide reduction rates estimated by [22] are used to determine the value of C f that best approximate the sulphide production on those experiments. Therefore, an optimization model is created to find the optimal value of C f matches the sulphide reduction presented in [22]. The model is the following Gradient descent, which is an iterative approach to solve optimization problems, is used to find the optimal value for C f yielding an optimal value of C f * = 1.37 × 10 −17 [1/s]. The attenuated reaction rate constants are shown in Table 3.

Lifetimes Calculations
Sulphide flux to the UFC is calculated following Fick's first law of diffusion at the RHS boundary and time t of the simulation using: in units [mol/s]. This work follows [6,7] assuming that copper corrosion reaction occurs instantaneously and that the RHS sulphide boundary condition has a constant value of zero at any time. Therefore, the depth of the corrosion per unit of time is given by the formulation in [13]: in units [m/s], where f HS is the unitless coefficient reflecting the stoichiometry of the copper corrosion reaction, M Cu is the molar mass of copper, ρ Cu is the density of copper. Table 4 shows the values used for these parameters. The total depth of corrosion up to time T is given by the accumulation of corrosion during the simulation time: where τ Corr is the cumulative corrosion depth from time t = t 0 = 0 to t = T in units [m], and T is the number of time steps simulated in the model. The lifetime of the canister is then computed as follows where θ Cu is the thickness of the copper layer (3 mm) and d SS Corr is the speed of corrosion once the simulation has reached steady state [m/t] which is normally after 3000 years.

Scenarios
Scenario A explores the assumption that the bentonite clay will suppress any SRB activity within the emplacement rooms and that it would only be active at the host rock. In this scenario, the rough nonuniform surface of the rock produced by the excavation would allow SRB activity on a thin layer of the gapfill bentonite closest to the rock interface. Consequently, sulphide production is restricted to within one centimeter from the rock interface. The gapfill and the HCB are kept in their designed densities shown in Table 2; therefore, a sharp change in density is present between the gapfill and the HCB layer. The concentrations of sulphate and sulphide in groundwater are those of Table 1.
Scenario B test the assumption that the SRB would be active, yet highly suppressed, inside the bentonite clay. This scenario also incorporates the design idea of increasing the dry density of the gapfill to match that of the HCB. Therefore, in this setup, the bentonite has a homogenous density throughout, and no distinction is made between the gapfill and HCB layers.

Uncertainty in Groundwater Concentration
In both scenarios, the concentration of sulphate is given an interval of uncertainty around the reference concentrations to explore a wider range of input values. A lognormal distribution is selected with parameters: µ SO 4 = ln SO 2− 4 0 and σ SO 4 = ln(2)/3; where is the sulphate groundwater reference concentration of Table 1. This selection of parameters corresponds to having a lognormal standard deviation of approximately 25% of the mean of the reference concentration.

Scenario A
Lifetimes under this scenario are above one million years in the deterministic model for both groundwater sulphate concentrations ( Table 5    The steady-state concentration profile ( Figure A3 in the Appendix A) shows that sulphide decreases linearly and has a slope change at the layer between the filling gap bentonite (gapfill) and the HCB. This result is expected as steady-state diffusion is linear. The concentration of sulphide at the host rock reaches a steady-state value of approximately 0.38% of the groundwater sulphate concentration, i.e., 3.8 and 1.2 ppm of sulphide for 1000 and 310 ppm of sulphate, respectively.
Results from the stochastic model also yielded lifetimes above or around one million years. The coefficient of variation is about the same (0.23) in both sulphate concentration runs, which follows from the specification of the lognormal distribution. This result implies that larger mean lifetime values entail larger standard deviations as noticeable in Figure 3.

Scenario B
This scenario, similar to A, has lifetimes above one million years for the deterministic model using both groundwater sulphate concentrations ( Table 6). The flux of sulphide to the UFC is slightly higher than the previous scenario, yielding values of approximately 1.4 and 0.4 [mol year −1 m −2 ] for groundwater sulphate concentrations of 1000 and 310 ppm respectively (Figure 4).
The steady-state sulphide concentrations at the rock interface are similar to those from scenario A, yielding values of approximately 4 and 1.2 ppm of sulphide for the sulphate concentrations of 1000 and 310 ppm respectively ( Figure A4 in Appendix A). Unlike scenario A, here the profile of sulphide concentration at a steady-state decreases nonlinearly as it reaches the UFC.
The mean and medium lifetimes values are above one million years for the stochastic model. However, for the higher sulphate concentration simulation, lower than one-million-year lifetimes are possible. The minimum lifetime sampled is around 630,000 years and because the following relation is observed in the results: * 1000 = * 310, the sulphate concentration necessary to yield this minimum lifetime is 2300 ppm, which is an unlikely high value. A histogram of lifetimes of the UFCs for this scenario is shown in Figure 5.

Scenario B
This scenario, similar to A, has lifetimes above one million years for the deterministic model using both groundwater sulphate concentrations ( Table 6). The flux of sulphide to the UFC is slightly higher than the previous scenario, yielding values of approximately 1.4 and 0.4 [mol year −1 m −2 ] for groundwater sulphate concentrations of 1000 and 310 ppm respectively (Figure 4).    The steady-state sulphide concentrations at the rock interface are similar to those from scenario A, yielding values of approximately 4 and 1.2 ppm of sulphide for the sulphate concentrations of 1000 and 310 ppm respectively ( Figure A4 in Appendix A). Unlike scenario A, here the profile of sulphide concentration at a steady-state decreases nonlinearly as it reaches the UFC.
The mean and medium lifetimes values are above one million years for the stochastic model. However, for the higher sulphate concentration simulation, lower than onemillion-year lifetimes are possible. The minimum lifetime sampled is around 630,000 years and because the following relation is observed in the results: Li f etime 1000ppm * 1000 = Li f etime 310ppm * 310, the sulphate concentration necessary to yield this minimum lifetime is 2300 ppm, which is an unlikely high value. A histogram of lifetimes of the UFCs for this scenario is shown in Figure 5.

Sensitivity Analysis
A sensitivity analysis is performed on the dry density of the gapfill bentonite and the optimal parameter of Equation (10). These parameters are selected because they encompass a level of uncertainty and their true value has a significant effect on the lifetimes of the UFCs.
The gapfill bentonite is the material directly in contact to the host rock; therefore, its dry density determines the reduction rate at the rock interface and the mass transport to the HCB buffer box. For scenario A, the gapfill dry density is allowed to vary from 1300 to 1700 [kg/m 3 ] and the density of the HCB is kept at the design value of 1700 [kg/m 3 ].

Sensitivity Analysis
A sensitivity analysis is performed on the dry density of the gapfill bentonite and the optimal parameter C f of Equation (10). These parameters are selected because they encompass a level of uncertainty and their true value has a significant effect on the lifetimes of the UFCs.
The gapfill bentonite is the material directly in contact to the host rock; therefore, its dry density determines the reduction rate at the rock interface and the mass transport to the HCB buffer box. For scenario A, the gapfill dry density is allowed to vary from 1300 to 1700 [kg/m 3 ] and the density of the HCB is kept at the design value of 1700 [kg/m 3 ]. Since scenario B assumes a uniform density throughout the gapfill and HCB, for this scenario the density of both materials is changed jointly following the range from 1300 to 2000 [kg/m 3 ].
The C f parameter, which controls the sulphate reduction inside of the bentonite clays (gapfill and HCB), is allowed to vary from four orders of magnitude lower to four orders of magnitude higher than its optimal value.
The highest concentration of sulphate (1000 ppm) is used as a baseline value for these experiments.

Discussion and Conclusions
This analysis, the first of its kind on the Canadian deep geological repository (DGR) design, has concentrated on the production of sulphide due to sulphate-reducing bacteria (SRB) activity and the subsequent effect on the time to failure or lifetime of the used fuel containers (UFCs). A first-order kinetics was approximated to describe the reduction of sulphate into sulphide carried out by the SRB. The reaction rate constant is then exponentially adjusted to reflect the high-pressure conditions inside the DGR. The two scenarios created yielded similar sulphide fluxes to the UFC and their lifetimes are above one million years.
The concentration of sulphide at the LHS boundary eventually reaches a steady-state value of between 1 and 4 ppm in both scenarios. This result is in agreement and supports previous studies that have assumed a constant sulphide concentration of 3 ppm [6].
Lifetimes of UFC are inversely proportional to the groundwater sulphate concentration in both scenarios. This relation is expected for scenario A which may, roughly speaking, be characterized as a diffusion system with a source term at the LHS boundary because its steady-state determines the rate of corrosion. However, this relation is somewhat unexpected for scenario B which has a source term throughout the bentonite clay.
Despite that the scenarios entail different assumptions on the behavior of SRB, their results are of the same order of magnitude, which bounds the uncertainty of the effect of having active SRB inside the DGR under the conditions tested. The assumption for scenario A is that SRB would be active only at the rock interface and the region closest to it. While scenario B assumes a more conservative situation, namely, that bacterial activity would be present across the bentonite layers. In both cases, the suppressing effect of high-density bentonite significantly reduces the mass transport of sulphide and corrosion to the UFCs.
The sensitivity analysis shows that more important than the true groundwater chemical concentration of corrosive species, the dry density of the gapfill bentonite or the HCB has a much more significant impact on the mean lifetime of the nuclear waste containers (Figures 6 and 7). lion years.
The concentration of sulphide at the LHS boundary eventually reaches a steady-state value of between 1 and 4 ppm in both scenarios. This result is in agreement and supports previous studies that have assumed a constant sulphide concentration of 3 ppm [6].
Lifetimes of UFC are inversely proportional to the groundwater sulphate concentration in both scenarios. This relation is expected for scenario A which may, roughly speaking, be characterized as a diffusion system with a source term at the LHS boundary because its steady-state determines the rate of corrosion. However, this relation is somewhat unexpected for scenario B which has a source term throughout the bentonite clay.
Despite that the scenarios entail different assumptions on the behavior of SRB, their results are of the same order of magnitude, which bounds the uncertainty of the effect of having active SRB inside the DGR under the conditions tested. The assumption for scenario A is that SRB would be active only at the rock interface and the region closest to it. While scenario B assumes a more conservative situation, namely, that bacterial activity would be present across the bentonite layers. In both cases, the suppressing effect of highdensity bentonite significantly reduces the mass transport of sulphide and corrosion to the UFCs.
The sensitivity analysis shows that more important than the true groundwater chemical concentration of corrosive species, the dry density of the gapfill bentonite or the HCB has a much more significant impact on the mean lifetime of the nuclear waste containers (Figures 6 and 7).  Equally significant is an adequate description of the SRB activity at the walls of the emplacement rooms and inside the bentonite clay. The sensitivity analysis on the parameter , that controls the sulphide production inside the bentonite clays, shows that lifetimes are quite sensitive to changes in this value. However, in order to have a more accurate description of the sulphide production, additional experiments on microbial growth and sulphide production would be needed, but this is hindered by the fact that carrying Equally significant is an adequate description of the SRB activity at the walls of the emplacement rooms and inside the bentonite clay. The sensitivity analysis on the parameter C f , that controls the sulphide production inside the bentonite clays, shows that lifetimes are quite sensitive to changes in this value. However, in order to have a more accurate description of the sulphide production, additional experiments on microbial growth and sulphide production would be needed, but this is hindered by the fact that carrying out these physical experiments takes months or longer periods of time. Equally important is the characterization of the rock type, groundwater, and microbial communities of the prospective location of the DGR so that the modeling parameters reflect the particularities of the chosen site.
Future work will focus on including temperature-dependent coefficients to account for the effect of the thermal energy output coming from the UFCs.

Conflicts of Interest:
The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A
To find the spatial increment ∆x in order to solve the reaction-diffusion model (1), the setup model of scenario B is solved for different spatial segments ranging from 10 to 10,000. The resulting UFC lifetimes and the computation time of each run are shown in Figure A1. The increase in percentage accuracy gained by increasing the number of spatial segments is shown in Figure A2, a line shows the threshold where the accuracy gain is less than 0.1%.    Error measure comparison of the optimization models (6) and (7):  Error measure comparison of the optimization models (6) and (7):