Models of geothermal reservoirs used for power generation need to simulate the production of multiphase fluid via a complex network of sources, as well as the reinjection of some of this fluid back into the system. The use of standard numerical methods to model such systems with interacting sources typically gives poor non-linear solver convergence and impractical limitations on time-step sizes. However, these issues can be overcome by modifying the Jacobian matrix used in the non-linear equation solution, adding extra terms to represent the interactions between sources. This approach has been incorporated into the Waiwera parallel, open-source geothermal flow simulator. The improved performance offered by this method is demonstrated through test problems and a real geothermal reservoir model.

A new geothermal reservoir modelling framework is discussed. The framework has two main objectives: first, all the geoscience and reservoir engineering data should be stored in a simple manner, not in any way dependent on the model grid to be used, and secondly, the data storage protocols should be easily transferable from one modelling project to the next. In our framework, some of the data are stored as part of a digital conceptual model created in Leapfrog Geothermal®, while most of the rest, including well-by-well reservoir engineering data, are stored in human- and machine-readable JSON files. Finally, some of the data related to the specification of production history and future scenario parameters are stored in control spreadsheets.

Waiwera is a new reservoir simulator designed to address the current and future challenges posed by modelling geothermal systems. It incorporates several innovations aimed at improving the robustness of phase transitions and steady-state model convergence, including a modified treatment of the multi-phase gravity term. The Waiwera code is parallelised and makes extensive use of the PETSc scientific computation library. It is written in object-oriented Fortran 2003 and released under a free, open-source software license. Results from demonstration problems show a high degree of consistency with comparable simulators such as TOUGH2, improved phase transition and steady-state convergence behaviour, good parallel scalability on large problems and demonstrate the suitability of Waiwera for complex models of real geothermal reservoirs.

Waiwera is a new open-source, parallelised geothermal flow simulator, developed initially as part of the Geothermal Supermodels research programme. Recent progress on Waiwera development has included implementation of a parallelised MINC method for simulating fractured media, a zone system for easier specification of model parameters over designated parts of the model mesh and usercustomizable simulation output including output from sources. Waiwera has adopted a new, modern build system and also uses a completely new unit-testing system, enabling easy installation and testing on a range of platforms. As well as describing these new features we also demonstrate Waiwera's performance on test problems and geothermal reservoir models.

Waiwera is an open-source, parallelized geothermal reservoir flow simulator which gives the ability to reliably run and calibrate complex geothermal reservoir models made up of millions of blocks. However, running such large models brings a number of new challenges that must be addressed in order to realise the benefits of the increased computational power. This paper describes these new challenges and our team's approaches for handling them. The challenges can be broadly categorised as either model setup or model post-processing. There are a small number of issues relating to simulation control that are also discussed. For model setup, our techniques for the following are discussed:  Grid definition and generation  Distributing rock properties to the model  Representing faults and fractures  Defining wells and their feedzones  Specifying boundary conditions including deep inputs  Upscaling/downscaling between different grid resolutions For model post-processing we make extensive use of python scripting tools and graphical user interfaces such as Leapfrog. We describe our techniques for extracting information from very large transient simulation output files and then presenting results in a range of meaningful formats to provide stakeholders with tools to support their decisions. Both model setup and model post-processing are demonstrated using examples from current real-world geothermal reservoir modelling projects.

The utility of marine protected areas (MPAs) as a means of protecting exploited species and conserving biodiversity within MPA boundaries is supported by strong empirical evidence. However, the potential contribution of MPAs to fished populations beyond their boundaries is still highly controversial; empirical measures are scarce and modelling studies have produced a range of predictions, including both positive and negative effects. Using a combination of genetic parentage and relatedness analysis, we measured larval subsidies to local fisheries replenishment for Australasian snapper ( Chrysophrys auratus : Sparidae) from a small (5.2 km 2 ), well-established, temperate, coastal MPA in northern New Zealand. Adult snapper within the MPA contributed an estimated 10.6% (95% CI: 5.5–18.1%) of newly settled juveniles to surrounding areas (approx. 400 km 2 ), with no decreasing trend in contributions up to 40 km away. Biophysical modelling of larval dispersal matched experimental data, showing larvae produced inside the MPA dispersed over a comparable distance. These results demonstrate that temperate MPAs have the potential to provide recruitment subsidies at magnitudes and spatial scales relevant to fisheries management. The validated biophysical model provides a cost-efficient opportunity to generalize these findings to other locations and climate conditions, and potentially informs the design of MPA networks for enhancing fisheries management.

In a companion paper we discuss the development of our new geothermal simulator, called Waiwera, and present results of tests on benchmark problems. In this paper we discuss the application of Waiwera to natural state simulations with models of Ngawha (New Zealand) and System B. These two models were chosen because they both involve mixtures of water and a non-condensible gas (carbon dioxide for Ngawha and air for System B) and they both have large contrasts in permeability. As a result the models run badly with a standard version of TOUGH2. Recently we have made changes in Waiwera and AUTOUGH2 (our version of TOUGH2) to the way some aspects of two-phase flow are handled. We found that these changes greatly improved the performance of the Ngawha and System B models. TOUGH2 and Waiwera solve the nonlinear mass and energy balance equations by using the Newton-Raphson (N-R) iterative procedure. At each N-R iteration a large system of linear equations has to be solved. Typically in a difficult natural state simulation as the time step increases the condition number of the linear equations degrades and the NR process may struggle to converge. With Waiwera we have found, by experiment, that careful selection of linear solvers and preconditioners can greatly improve numerical performance. Results for the Ngawha and System B models using different linear solvers and preconditioners are discussed below.

The Menengai geothermal system is hosted in a ring-like caldera a short distance north of the city of Nakuru in Kenya. In 2010 it became the second field in Kenya to be developed for energy generation, following Olkaria, which has proven to be successful and is currently undergoing expansion. In 2013 a preliminary numerical model of the Menengai field was developed (Kipyego et al., 2013) which highlighted the structural control of the system and provided initial insights into its behaviour. The model also demonstrated that the standard TOUGH2 simulator was not sufficient for modelling the supercritical conditions known to exist at depths below 3200m. In this work a new model is developed in which supercritical conditions are imposed in the upflow region of the bottom boundary. The model is run using the University of Auckland supercritical version of the TOUGH2 simulator (Croucher and O’Sullivan, 2008) with the air-water equation of state. Comparisons of modelled downhole temperatures with field data show that the supercritical model achieves a significantly better match than the preliminary model and that artificially low permeabilities are not required at the bottom of the model to replicate the high temperatures encountered. The possible production and reinjection scenarios investigated previously are implemented using the supercritical model and the differences discussed.

The “Geothermal Supermodels” project is a four-year New Zealand-based research programme, a major part of which is the development of a new open-source, parallelized geothermal reservoir flow simulator, called “Waiwera”. This paper describes recent progress on Waiwera. Features recently added include simulation of noncondensible gases, a flexible framework for modelling source and sink terms, support for radial and twodimensional models and capillary pressure functions. Numerical performance has been improved via nondimensionalising the primary thermodynamic variables, and better algorithms for polynomial evaluation. Thermodynamic behaviour has also been improved by implementing more sophisticated algorithms for phase transitions and for non-condensible gas energy of solution. Work has also been done to increase the efficiency of linear equation assembly and also solution, via improved preconditioning techniques. We demonstrate the performance of Waiwera via results from a range of benchmark tests.

Reservoir models have become important and increasingly widely-used tools for the development and management of geothermal systems. In many countries both regulators and financial institutions require that some level of numerical modelling is carried out as part of a feasibility study. The geothermal modelling team at the Geothermal Institute and Department of Engineering Science at the University of Auckland (UoA) have been developing reservoir models since the 1980s and are among the world leaders in the technology. Through the taught programme for the Postgraduate Certificate in Geothermal Technology and research degrees the UoA is also responsible for training many of the world’s geothermal reservoir modellers. The worked present here is a summary of the state-of-the-art in geothermal reservoir modelling and how it is used to aid strategic decision-making. Case studies from different countries are presented covering topics such as resource planning and optimisation, re-injection strategies, predicting subsidence, environmental sustainability and enhanced permeability.

The validity of subsurface flow simulations strongly depends on the accuracy of relevant rock property values and their distribution in space. In realistic simulations, this spatial distribution is based on two geological considerations: (1) the subsurface structural setting, and (2) smaller-scale heterogeneity within a hydrostratigraphic unit. Both aspects are subject to uncertainty, whereas techniques to address heterogeneity are well established, no general method exists to evaluate the influence of structural uncertainties. We present a method to include structural geological data (e.g. observations of geological contacts and faults) directly into an inversion framework, with the aim of enabling the inversion routine to adapt a full 3-D geological model with a set of geological parameters. In order to achieve this aim, we use a set of Python modules to combine several pre-existing codes into one workflow, to facilitate the consideration of a structural model in the typical model evaluation steps of sensitivity analysis, parameter estimation, and uncertainty propagation analysis. In a synthetic study, we then test the application of these three steps to analyse CO2 injection into an anticlinal structure with the potential of leakage through a fault zone. We consider several parts of the structural setting as uncertain, most importantly the position of the fault zone. We then evaluate (1) how sensitive CO2 arriving in several observation wells would be with respect to the geological parameters, (2) if it would be possible to determine the leak location from observations in shallow wells, and (3) how parametric uncertainty affects the expected CO2 leakage amount. In all these cases, our main focus is to consider the influence of the primary geological data on model outputs. We demonstrate that the integration of structural data into the iTOUGH2 framework enables the inversion routines to adapt the geological model, i.e. to re-generate the entire structural model based on changes in several sensitive geological parameters. Our workflow is a step towards a combined analysis of uncertainties not only in local heterogeneities but in the structural setting as well, for a more complete integration of geological knowledge into conceptual and numerical models.

A high-resolution 3D biophysical model was used to investigate the patterns of larval transport for an important commercial and recreational temperate fish, snapper Pagrus auratus, from a well-established marine reserve (Cape Rodney to Okakari Point marine reserve, CROP), and spawning ground. Our focus was to study the effects of local hydrodynamics, contrasting larval vertical behaviours and changing El Nino-Southern Oscillation (ENSO) cycles (via their effect on wind forcing) on the potential larval supply to adjacent non-protected areas. The model suggests the CROP marine reserve provides significant larval subsidies within a relatively small scale (<= 40 km), the details of which depend on larval behaviour and ENSO patterns. Changing ENSO patterns mostly affected the dispersal direction of larvae, while larval behaviour was a strong driver only under La Nina conditions. Modelling suggests that under El Nino conditions, snapper larvae are more likely to settle successfully and hence make a strong contribution to surrounding fished stocks. Understanding the contribution of MPAs to fisheries management will require multiple approaches. Modelling will help identify the strength of contributing physical and biological factors, and in due course enable site-and time-specific examples of larval subsidy to be generalised.

Numerical modelling has become an important tool in managing geothermal systems and planning their exploitation for renewable energy production. The TOUGH2 simulator has been the industry standard tool for developing numerical models for many years. It includes several different equation-of-state modules and thus can be used for modelling many different kinds of geothermal fields as well as other complex sub-surface flow problems such as carbon sequestration and nuclear waste storage. Under certain conditions TOUGH2 simulations stall at a relatively small time step size and are unable run up to the very large times required for a natural state simulation. This behaviour leads to slow model development and also poses a significant obstacle to inverse modelling using iTOUGH2 or PEST as forward simulations are more computationally expensive and may not finish. In previous work the authors identified and analysed conditions leading to stalled simulations and proposed corrections to the air-water and CO2-water equation-of-state modules to eliminate the behaviour.While these improvements prevent simulations from stalling at small time steps they do not prevent a limiting time step size from occurring as simulations approach a large simulation time. This behaviour is observed often, affects all equation-of-state modules and can add significant computational time to simulations. Improvements to the convergence criteria used by the TOUGH2 algorithm during its internal Newton-Raphson iteration are proposed following an analysis of the algorithm's behaviour as the time step size becomes very large. The results presented show that this improvement can dramatically reduce the total simulation time for all equation-of-state modules. In particular the improvement can reduce the computational time required for calculating numerical derivatives during inverse modelling by an order of magnitude. This has a profound effect upon our ability to carry out inverse modelling and uncertainty analysis of numerical models of large-scale, real geothermal systems.

Numerical modelling has become an important tool in managing geothermal systems and planning their exploitation. The TOUGH2 simulator is the industry standard tool for developing numerical models. It includes several different equation-of-state modules and thus can be used for modelling many different kinds of geothermal fields as well as other complex sub-surface flow problems. Under certain conditions TOUGH2 has difficulty in running a model up to the very large times required for a natural state simulation as it stalls at a relatively small time step size. This problem leads to slow model development and also poses a significant obstacle to inverse modelling using iTOUGH2 or PEST as forward simulations are more computationally expensive and may not finish. In this paper two conditions leading to stalled simulations are identified and analysed. A correction to the air-water and CO2-water equation-of-state modules eliminates the first problem. The second is eliminated by an adjustment to the saturation temperature calculation. Results are presented showing a dramatic improvement in the simulation time required to achieve a steady state for large models of real geothermal systems.

The TOUGH2 simulator and its parallel version TOUGH2-MP continue to be the industry standards for the development of numerical models of geothermal systems. Increasing processing power enables us to simulate larger, more complex systems, while improved data collection and remote sensing techniques provide an ever greater suite of observations against which to calibrate the model. As a consequence, the raw text input files that control the simulations grow increasingly cumbersome to construct, while post-processing of model output becomes more challenging. We use the Python scripting language and the PyTOUGH library to interact with TOUGH2 in a number of novel ways to control simulation tasks accurately and efficiently. In this paper several examples of the use of PyTOUGH are described. A new method for automatically generating geothermal production wells is described that allows well information to be stored in concise, readable and easily updatable files. Examples are given of sequential modification of model geometries to represent effects such as eruptions and excavation. The generation of atmosphere blocks that vary in space and time due to changes in lake levels and extreme altitude are presented. Methods for controlling sequential simulations with both TOUGH2 and TOUGH2-MP are described including techniques for dealing with numerical nonconvergence. Finally, PyTOUGH's ability to combine and post-process results from multiple simulations and different types of output files is discussed.

The subsurface structural setting, as a part of the conceptual model, can have a strong influence on the outcome of a flow simulation. The structural geological model should therefore be considered in uncertainty studies or parameter estimations. However, in practice this inclusion is not easily possible because the generation of an input mesh for a flow simulation, populated with properties according to a geological model, requires at least some manual work. We present here a method that removes this gap between geological data and flow simulation and propose that it will enable insights into the significance of structural uncertainties in the context of inversion modeling that is not yet possible. Practical applications of the methods are demonstrated with two simple examples: the inversion of the distance to an impermeable boundary in an aquifer, and the significance of uncertainties in the exact position of structural geological observations in a CO 2 injection scenario with possible leakage over a fault with an uncertain position. The presented examples show that the integration of geological structures into an inversion framework enables a consideration of structural uncertainties that is not easily possible with the methods available to date. Geological data can be highly uncertain and limited, and our examples emphasize the fact that these uncertainties can strongly influence the outcomes of a flow simulation.