2.6.2.1.2 Groundwater numerical modelling

The main geological domains of interest in the Galilee subregion are the Galilee Basin, the Eromanga Basin and the Cenozoic to Quaternary cover sediments and alluvial deposits (refer to companion product 1.1 for the Galilee subregion (Evans et al., 2014)). The Galilee Basin stratigraphic units crop out in the east of the subregion (Figure 3) and contain the main target formations for coal development: Betts Creek beds, Bandanna Formation and Colinlea Sandstone. For BAs these units are collectively known as the upper Permian coal measures. The Rewan Group, Clematis Group and Moolayember Formation separate the upper Permian coal measures from the Jurassic to Cretaceous strata of the Eromanga Basin. The upper Permian coal measures are underlain by the Joe Joe Group. The aquifers hosted in these stratigraphic units are mostly confined and feed several springs as well as providing water to public water supply and stock and domestic bores. The Cenozoic and alluvial sediments, not shown in Figure 3, contain watertable aquifers, some of which are perched (refer to companion product 2.1-2.2 for the Galilee subregion (Evans et al., 2018a).

Figure 3 Schematic representation of proposed coal resource development operations in the Galilee subregion

GAB = Great Artesian Basin; Fm = formation

Detailed surface geology is provided in Figure 4.

The coal developments included in the CRDP are all situated in the east of the Galilee Basin due to the proximity of the upper Permian coal measures to the surface (refer to companion product 2.3 for the Galilee subregion (Evans et al., 2018b)). The proposed mining operations are therefore limited to the Belyando river basin, a headwater catchment of the larger Burdekin river basin (companion product 1.1 for the Galilee subregion (Evans et al., 2014)).

The conceptual model for the Galilee subregion indicates that large coal mining operations included in the CRDP will result in dewatering of the upper Permian coal measures in the vicinity of the operations (Figure 3 and companion product 2.3 for the Galilee subregion (Evans et al., 2018b)). This dewatering has the potential to cumulatively impact regional groundwater systems due to lateral and vertical propagation of groundwater level and pressure reductions, both within the upper Permian coal measures and within adjacent hydrostratigraphic units.

Coal mining may potentially affect surface water runoff as rainfall is intercepted by mine workings and therefore no longer contributes to streamflow. Also, coal mining may cause a depressurisation of or drawdown in aquifers interacting with the surface water features, which then indirectly may alter streamflow.

As discussed in detail in Section 2.3.4.1 in companion product 2.3 (Evans et al., 2018b), the CRDP in the Galilee subregion consists of 17 proposed new coal and coal seam gas (CSG) resource development projects. For only seven of these projects is sufficient information available to be included in the modelling (Figure 4). These are all coal projects: Alpha, Carmichael, China First, China Stone, Kevin’s Corner, South Galilee and Hyde Park. No CSG developments are included in the modelling.

As outlined in submethodology M07 (as listed in Table 1) for groundwater modelling (Crosbie et al., 2016), different model types and model codes are chosen in BA, depending on the specific requirements of each subregion. The main goal of each groundwater model in BA remains, however, to deliver spatially explicit model outputs that are used as inputs to other BA models, including surface water modelling, uncertainty analysis and receptor impact modelling, and to directly evaluate change in water resources. Table 3 lists the criteria a groundwater model in BA needs to satisfy to be considered fit-for-purpose for BA. Beneath the table, these fit-for-purpose criteria are discussed briefly for the numerical modelling approach taken in the Galilee subregion. The remainder of this product describes in greater detail the numerical modelling, and the underlying assumptions and their effect on predictions.

Table 3 Capability requirements of the groundwater numerical modelling approach in the Galilee subregion

2.6.2.1.2.1 Prediction of hydrological response variables

The objective of the numerical modelling undertaken as part of a BA is to probabilistically assess hydrological changes arising from coal resource development at water-dependent assets and model nodes (see companion submethodology M07 (as listed in Table 1) for groundwater modelling (Crosbie et al., 2016)). The groundwater and surface water modelling predicts changes in hydrological response variables, the hydrological characteristics of the system or landscape class that potentially change due to coal resource development. These hydrological response variables are the input for receptor impact models that will evaluate how the change in hydrology and hydrogeology results in a change in the economic, social or ecological value of assets.

The model is required to estimate drawdowns caused by the coal resource developments outlined in the modelled CRDP at model nodes in the upper Permian coal measures, the Clematis Group and the Cenozoic and Quaternary cover sediments. Figure 4 shows the locations of these points and their hydrostratigraphic unit. Table 4 conveys the same information, together with a short description of each point, its water source and the model layer the point is assigned to.

Figure 4 Points for which drawdown is required to be estimated in the Galilee subregion (see Table 4 for point types)

Data: Bioregional Assessment Programme (Dataset 1)

Table 4 Summary of points for which drawdown due to coal resource development in the Galilee subregion are required

GDE = groundwater-dependent ecosystem; GW licence = groundwater extraction licence; TWS = town water supply

Model nodes were positioned to represent one or more identified water-dependent assets, such as clusters of groundwater bores, groundwater-dependent ecosystems, stream confluences and river reaches.

Most model nodes are located in the alluvium, the majority of which are associated with groundwater-dependent ecosystems with one model node assigned to a town water supply. The Ronlow beds, which are part of the Eromanga Basin, have two model nodes associated with springs. There are 13 model nodes associated with the Clematis Group and they include groundwater-dependent ecosystems, licensed groundwater extractions, monitoring bores, springs and the Jericho town water supply. The upper Permian coal measures have five model nodes, which are two groundwater extraction licences and three springs. There are no model nodes associated with the main aquitard units in the modelling domain, such as the Rewan Group, as the groundwater model used in this assessment (an analytic element model) does not provide estimates of drawdown for aquitards.

The hydrological response variables for groundwater are drawdown (dmax) and year of maximum change (tmax) at model nodes, where drawdown is defined as the difference in groundwater level between the baseline and CRDP, due to additional coal resource development.

For surface water, nine hydrological response variables are defined in companion submethodology M06 (as listed in Table 1) for surface water modelling (Viney, 2016) at 61 nodes along the stream network (companion product 2.6.1 for the Galilee subregion (Karim et al., 2018)).

Changes in the groundwater system can propagate to changes in the surface water system, which means that simulating the change in hydrological response variables at the various model nodes necessitates the development of an integrated surface water – groundwater modelling approach. The groundwater and surface water, however, operate at very different spatial and temporal scales. The surface water obviously is bound to river channels and floodplains. Streamflow is very responsive to individual rainfall events, requiring at least a daily temporal resolution to capture its ephemeral nature. Groundwater dynamics in the alluvial and Cenozoic deposits are mostly local and controlled through interactions with surface events, such as high rainfall or flooding (Section 2.1.3.2.2 in companion product 2.1-2.2 for the Galilee subregion (Evans et al., 2018a)). Capturing this dynamic in a numerical model necessitates at a minimum a monthly resolution.

The deeper hydrogeological units hosted in the Eromanga and Galilee basins are much more extensive, both horizontally and vertically. The groundwater dynamics are very slow. In the outcrop zones of the units, there are indications that groundwater levels are influenced by recharge events. In the deeper, confined parts of the hydrogeological units there is no indication that groundwater dynamics are affected by recharge and discharge processes. Simulating groundwater flow in the deeper hydrogeological units requires a spatially extensive model, but a high temporal resolution is not essential.

While fully coupled surface water – groundwater model codes are available (e.g. Hydrogeosphere, Brunner and Simmons, 2012), their use was not deemed to be justified within BA due to the high data requirements for parameterisation and due to operational constraints. The latter relates mainly to the general numerical instability of such models and long run times which would severely limit a probabilistic uncertainty analysis that requires the models to be evaluated 100s of times with vastly different parameter sets.

For the Galilee subregion, a pragmatic coupling of two models was developed, consisting of a regional groundwater model to simulate the change to the groundwater systems of the subregion and a rainfall-runoff model to simulate the change to the surface water systems of the subregion (Figure 5). The individual models have different spatial and temporal resolution which requires a set of customised processing steps to up or downscale model data to allow the models to be linked.

Figure 5 Model sequence for the Galilee subregion

AWRA-L = rainfall-runoff model; dmax = maximum difference in drawdown for one realisation within an ensemble of groundwater modelling runs results, obtained by choosing the maximum of the time series of differences between two futures; GW = groundwater; GW AEM = regional groundwater analytic element model; ΔHRV = change in hydrological response variable; tmax = year of maximum change; ΔQb = change in surface water – groundwater interaction flux; Qt_b = total streamflow baseline; Qt_c = total streamflow coal resource development pathway (CRDP); SW = surface water

The regional groundwater model is an analytic element model (referred to as GW AEM), designed to simulate the change in drawdown at the points shown in Figure 4 and the change in surface water – groundwater flux. As there is no coal resource development under baseline conditions, the drawdown and change in surface water – groundwater flux due to coal development is zero. There is therefore no need for a separate baseline conditions run for the groundwater model. The change in surface water – groundwater flux simulated with the CRDP run of the analytic element model, , is taken into account in the Australian Water Resources Assessment landscape model (AWRA-L) surface water model generated streamflow. The change in a number of hydrological response variables is modelled at the model nodes. The simulation of river management or routing of streamflow through the river network with a river model is not necessary as the salient features of streamflow can be simulated solely with a rainfall-runoff model (see companion submethodology M06 (as listed in Table 1) for surface water modelling (Viney, 2016)).

The AWRA-L baseline run simulates streamflow at surface water model nodes without any active mines. The AWRA-L CRDP run simulates streamflow at the surface water model nodes incorporating the effect of approved open-cut and underground coal mines in the CRDP. The total surface water streamflow under CRDP, is obtained as:

where is the total surface water flow under baseline and is the change in surface water – groundwater flux between baseline and CRDP computed by the groundwater model.

Due to the ephemeral nature of the river system and the difference in temporal resolution between the groundwater model and the surface water model it is possible that the change in surface water – groundwater flux is larger than the total available streamflow. Streamflow under those conditions is reduced to zero. This would typically occur during periods of low flows driven by baseflow. Any further loss of baseflow due to mine dewatering may result in the river running dry. Note that many of the smaller creeks and rivers are maximally losing. The surface water – groundwater flux for such systems is independent of the groundwater level and can therefore not be affected by drawdown.

The time series of and are summarised in the nine hydrological response variables to highlight different aspects of the hydrograph. These hydrological response variables will inform the receptor impact models for the model nodes associated with surface water.

2.6.2.1.2.2 Design and construction

According to the Australian groundwater modelling guidelines (Barnett et al., 2012), it is essential to design and construct the groundwater model in function of clearly stated objectives and to provide a model confidence level. The objective of the modelling is explicitly stated in the previous section. The model confidence level is an a priori categorisation of a groundwater model to reflect its predictive capability in function of the model complexity, prediction timeframe and data availability. As clarified in submethodology M07 (Crosbie et al., 2016), the groundwater models in the Bioregional Assessment Programme are all classified as level 1, the lowest level, as they are required to make predictions of unprecedented stresses at time frames longer than periods with data available to constrain the model.

The objectives of the modelling are not to simulate the state of groundwater in the future under baseline and coal resource development conditions, but to quantify the difference between those two futures. This is a very important nuance to the modelling objectives as it enables a number of simplifying assumptions based on the principle of superposition (Reilly et al., 1987). The principle of superposition means that for linear systems, the solution to a problem involving multiple inputs (or stresses) is equal to the sum of the solutions to a set of simpler individual problems that form the composite problem. To simulate the effect of change in stress, such as depressurisation and dewatering for coal resource development, it is therefore sufficient to only know the change in stress. It is not necessary to know the initial conditions in the aquifer or the other fluxes and stresses, provided these do not change due to the change in stress (Barlow and Leake, 2012).

The principle of superposition underpins most of the pumping test interpretations (Kruseman and de Ridder, 1994); aquifer parameters are inferred from the change in stress (pumping rate) and change in groundwater level (drawdown).

The principle of superposition is only valid for linear systems (i.e. systems where the response to a change in stress is proportional to the change in the stress). In other words, where a doubling of stress will result in a doubling of the response. In groundwater flow dynamics this condition is satisfied for confined aquifers. Unconfined aquifers are not strictly linear, as the transmissivity depends on the saturated thickness. Reilly et al. (1987) and Rassam et al. (2004) do show, however, that the concepts are still valid for mild violations of the linearity conditions. The assumption most likely to be violated when simulating dewatering of an unconfined aquifer is that the transmissivity is no longer constant as the saturated thickness decreases during dewatering. Singh and Atkins (1985) provide an overview of different analytic solutions that take into account this change in transmissivity in the context of mine dewatering. From these equations it is apparent that not accounting for changing transmissivity will lead to under estimating drawdowns in the vicinity of the mine. At greater distances from the mine, however, drawdowns will be overestimated. This is consistent with the interpretation of pumping tests, where high transmissivity results in a large, but shallow cone of depression, while a low transmissivity results in a small, but deep cone of depression (Kruseman and de Ridder, 1994). Barlow and Leake (2012) discuss the conditions for which analytical solutions are valid in the context of surface water – groundwater interaction. The principle of superposition is not valid if the connection status changes due to the stress (e.g. if the river changes from losing or gaining connected to losing disconnected). When the surface water and groundwater system is connected, the change in surface water – groundwater flux is proportional to the change in groundwater level. When the system becomes disconnected, the system is maximally losing and the flux is no longer proportional to the groundwater level (Lamontagne et al., 2014).

As such, the concept of superposition can be implemented in any groundwater modelling code, see, for example, Leake et al. (2008). The analytic element modelling code, implemented through TTim version 0.3 (Bakker, 2015), for this regional model is chosen because the scripted nature of the modelling code lends itself very well to automated changing of all aspects of the model, required for the sensitivity and uncertainty analysis (Bakker and Kelson, 2009). Another advantage is that due to the grid-independent nature of the analytic element method, the resolution of the results is not dependent on the grid discretisation as is the case for finite difference or element models. Outputs can thus be generated at any arbitrary locations within the modelling domain.

Further technical details of the conceptualisation, parameterisation and implementation are documented in the following sections of this product for the groundwater modelling and in companion product 2.6.1 for the Galilee subregion (Karim et al., 2018) for surface water modelling.

2.6.2.1.2.3 Integration with sensitivity and uncertainty analysis workflow

Companion submethodology M09 (as listed in Table 1) (Peeters et al., 2016) discusses in detail the propagation of uncertainty through numerical models in BAs. The workflow outlined in this product is tailored to numerical models with long run times and where observations of the groundwater system can be used to constrain the model parameters.

The analytic element model for the Galilee subregion has a very short run time and, because only the change in the system is simulated and not the state variables, it is not possible to use state observations, such as fluxes and groundwater levels, to constrain the model parameters through Approximate Bayesian Computation Markov chain Monte Carlo sampling of prior parameter ensembles.

These ensembles capture the range of each parameter that is deemed likely, based on the available local information and international literature. This is discussed in greater detail in the parameterisation section (Section 2.6.2.6).

The uncertainty analysis then consists of a direct sampling of these ensembles of parameter values with the analytic element model to arrive at an ensemble of predictions of additional drawdown, year of maximum change and change in surface water exchange flux.

2.6.2.1.2.4 Water balance components

A secondary objective of the numerical models is to inform the water balance assessment (companion product 2.5). The AWRA-L model produces surface water estimates of the water balance under baseline and coal resource development futures and can therefore be used in that assessment. The analytic element model, however, only simulates the change in stress due to coal resource development. Its model output therefore has no information on other components of the regional water balance such as recharge or lateral exchange fluxes.

2.6.2.1.2.5 Transparent and reproducible model outputs

An over-arching requirement of the BAs is for all model outputs to be transparent and reproducible.

The model files, including the pre- and post-processing scripts and executables, will be made available through the BA data repository.

As the evaluation of the model chain is a highly automated and scripted process, it is possible to reproduce the results reported in this document using the scripts and executables in the repository, provided the computational resources are available. Scripting the modelling process has the added advantage that it minimises the risk of human error in the data analysis.