Skip to main content

The influence of surficial features in lava flow modelling


Lava flows can cause substantial and immediate damage to the built environment and affect the economy and society over days through to decades. Lava flow modelling can be undertaken to help stakeholders prepare for and respond to lava flow crises. Traditionally, lava flow modelling is conducted on a digital elevation model, but this type of representation of the surface may not be appropriate for all settings. Indeed, we suggest that in urban areas a digital surface model may more accurately capture all of the obstacles a lava flow would encounter. We use three effusive eruption scenarios in the well-studied Auckland Volcanic Field (New Zealand) to demonstrate the difference between modelling on an elevation model versus on a surface model. The influence of surficial features on lava flow modelling results is quantified using a modified Jaccard coefficient. For the scenario in the most urbanised environment, the Jaccard coefficient is 40%, indicating less than half of the footprints overlap, while for the scenario in the least urbanised environment, the Jaccard coefficient is 90%, indicating substantial overlap. We find that manmade surficial features can influence the hazard posed by lava flows and that a digital surface model may be more applicable in highly modified environments.


Recent effusive volcanic eruptions such as the 2018 Kīlauea eruption in Hawaii, USA, have reminded the global community how disruptive lava flows can be. Although few deaths are attributed to lava flows (Harris 2015; Brown et al. 2017), they can have severe, immediate impacts to the built environment and prolonged economic and societal consequences. For example, the 1973 Vestmannaeyjar Volcanic Field eruption on the Icelandic island of Heimaey threatened the local fishing harbour (e.g. Williams and Moore 1983; McPhee 1989) and had lasting financial impacts on the Icelandic economy (Morgan 2000). More recently, lava flows during the 2002–2003 eruption of Mt. Etna (Sicily, Italy) partially inundated tourist facilities and road networks on the volcano’s southern flank (Villa 2002; Bonaccorsco et al. 2016; Rongo et al. 2016). In yet another example, over the course of 3 months in 2018, the development of the lava flow field on Kīlauea’s Lower East Rift Zone (Hawaii, USA) destroyed over 700 houses (Neal et al. 2019) and the local government has reported small business closures over a year later (County of Hawai‘i 2019). These consequences have been recognised for decades. In fact, Thomas Jaggar devised ways to protect Hilo Harbour (Hawaii, USA) in future eruptions in the 1940s (Jaggar 1945). Such planning continues to this day (CDEM 2015). As with all hazards, pre-crisis planning can minimise the psychological and physical impacts and recovery costs for communities affected by lava flows (UNDDR 2019). Preparation and mitigation actions can take many forms, such as the development of plans and policies (e.g. the creation of the Auckland Volcanic Field contingency and evacuation plans (CDEM 2015; Wild et al. 2019)), lava flow modelling (e.g. the DOWNFLOW modelling undertaken by Favalli et al. (2006, 2009a) and Chirico et al. (2009) in Goma, Democratic Republic of Congo), volunteer and professional trainings (e.g. the electricity utility on Hawai’i Island creating trainings about how to strengthen their networks to withstand the high temperatures of lava flows (Tsang et al. 2019)), and exercises (e.g. Exercise Ruamoko, a New Zealand all-of-nation desktop exercise focused on unrest in the Auckland Volcanic Field (Brunsdon and Park 2009; Lindsay et al. 2010)). All these activities rely on understanding the potential hazard(s), in this case when and where lava flows can occur.

Assessing the potential lava flow hazard to a community can be difficult for a variety of reasons. Many volcanoes have not been studied in detail, and even well-studied volcanoes commonly lack extensive research on specific characteristics of previous individual eruptions; metrics such as recurrence interval or average lava flow episode volume are often unknown (e.g. Mt. Cameroon, Cameroon; Bonne et al. 2008; Favalli et al. 2012). Additionally, predicting the location of the next volcanic vent is critical to predicting lava flow inundation areas (Connor et al. 2012; Bilotta et al. 2019); this can be especially challenging in volcanic fields where the area under consideration for vent opening may be particularly large with few spatio-temporal trends (e.g. Allen and Smith 1994; Gallant et al. 2018).

In some volcanic regions the natural topography has been anthropologically altered during urbanisation (e.g. Al-Madinah, Saudi Arabia (Runge 2015), Auckland, New Zealand (von Hochstetter 1859; Searle 1964)). In these cases, the increasing number of obstacles (i.e. elements in the built environment) a lava flow may encounter may be increased; in other words, lava flows will not simply follow the natural topography (Tsang et al. 2019). Despite substantial anthropogenic alteration to some natural topographies, lava flow modelling is still normally conducted on digital elevation models (DEMs) that do not consider such alterations (e.g. Kereszturi et al. 2014).

Digital elevation models are commonly used when modelling volcanic hazards. One modern method to create a DEM is using a LAS dataset, which is created using a Light (intensity) Detection and Ranging (LiDAR) system (e.g. White and Wang 2003; Favalli et al. 2009b). To create the dataset, a LiDAR system is flown over the area of interest while a laser in ultraviolet, visible, or infrared wavelengths is directed at the ground (e.g. Wagner et al. 2006; Heidemann 2018). As the LiDAR system flies over the area of interest, it measures the reflectance of the original wave emitted (e.g. Axelsson 1999; Wagner et al. 2006; Favalli et al. 2009b). Several returns are measured (e.g. Wagner et al. 2006), with each return providing different information about the ground below. The returns are numbered based on the order in which they are sensed (e.g. Heidemann 2018). The number of returns measured depends on the system being used, but the early returns detail surficial features while later returns describe lower features (e.g. the bare ground; Axelsson 1999). The first return is commonly used to create a digital surface model (DSM), which includes the bare ground and surficial features including trees, buildings, etc. (e.g. Wagner et al. 2006; Heidemann 2018). One of the later returns is used to create a DEM. In-between returns can be used to create hybrid surface models (e.g. a DSM that does not include trees) through the process of classification or filtration (e.g. Axelsson 1999; Wagner et al. 2006). A second modern method to create a DEM is using Structure-from-Motion photogrammetry (e.g. Westoby et al. 2012; Ouedraogo et al. 2014). In this method, a set of overlapping but offset photographs is taken (e.g. Westoby et al. 2012; Dietrich 2015). Algorithms can then identify features from multiple images at different angles to create a three-dimensional surface. Sufficient images are taken to blanket an area with images, enabling a DEM (or DSM) to be created (e.g. Westoby et al. 2012; Ouedraogo et al. 2014; Dietrich 2015).

Here, we present lava flow hazard modelling for three effusive eruption scenarios in the Auckland Volcanic Field (AVF), which underlies the city of Auckland, New Zealand, to highlight the effects of using different elevation models, namely DEMs and DSMs. In Section 2, we provide an overview of the AVF and the scenarios used, which were created in the context of the Determining Volcanic Risk in Auckland (DEVORA) research programme. In Section 3, we describe our methods, including how we selected the MOLASSES lava flow model and how we compared our modelling results. In Section 4, we present the MOLASSES results from modelling the DEVORA Scenarios on both a DEM and a DSM. Finally, in Section 5, we detail some of the implications of modelling lava flows on a DEM versus a DSM and compare our modelling results to the lava flow hazard footprints developed in an earlier version of the DEVORA scenarios.

Auckland volcanic field, North Island, New Zealand

The AVF is a monogenetic, basaltic volcanic field situated on the North Island of New Zealand (Fig. 1). It has been active for approximately 190kyr with the most recent eruption 550 yr BP (Leonard et al. 2017). Since the last eruption, New Zealand’s largest city, Auckland, has been built on top of the volcanic field, and consequently, Auckland could be severely impacted by a future AVF eruption (Deligne et al. 2015; Hayes et al. 2017). There are neither strong spatial-temporal nor spatial-volumetric trends in the AVF (Bebbington and Cronin 2011; Bebbington 2015), although there are trends correlating bulk rock chemistry with distance to neighbouring vents, volume, and time (Le Corvec et al. 2013; McGee et al. 2015). Without clear spatial-temporal or spatial-volumetric trends (Bebbington and Cronin 2011; Bebbington 2015), it is difficult to assess probabilistic hazard for the city, including location and probability of the next vent opening (Lindsay et al. 2010; Sandri et al. 2012; Hayes et al. 2018). The Determining Volcanic Risk in Auckland (DEVORA) research programme ( was established in 2008 to study the volcanic risk in the city from local and distal volcanoes. To aid stakeholders in preparing for an eruption, the DEVORA research team developed eight hypothetical AVF eruption scenarios (Hayes et al. 2018, 2020), four of which include lava flows (Fig. 1).

Fig. 1
figure 1

Map showing Auckland with the volcanic centres of the AVF. The black oval shows the possible boundary of the AVF based on the known eruptive centres, shown as black triangles (Runge 2015; Kereszturi et al. 2014). The DEVORA scenario vents that we use as case studies are denoted with solid, red triangles; the letters in each triangle indicate the scenario. B is for Birkenhead; the Birkenhead scenario has two vents, both of which fall within the northernmost red triangle. M is for Mount Eden; O is for Ōtāhuhu. Inset shows the North Island of New Zealand and major cities; Auckland is boxed. Figure modified from Bebbington (2015) and Hayes et al. (2018). Data from Kermode (1992) and Land Information New Zealand

DEVORA scenarios

The DEVORA scenarios are hypothetical eruption sequences, incorporating multi-hazard modelling (Hayes et al. 2018, 2020). The hazards were modelled to represent the full range of possible eruptive phenomena and hazard intensities (Hayes et al. 2018, 2020). Not all hazards are included in every scenario; rather, hazards are included based on how frequently they have occurred in past AVF eruptions, and the local environmental conditions at the scenario vent area (e.g. onshore/offshore; depth to water table). The resulting DEVORA scenarios have been or will be used in a variety of projects conducted by DEVORA researchers and their stakeholders, including calculating the impact of an AVF eruption on electricity transmission networks (Tsang SWR, Lindsay JM, Kennedy B, Deligne: Thermal impacts of basaltic lava flows to buried infrastructure: workflow to determine the hazard, submitted), reviewing the local government volcanic contingency plan (A. Doherty, pers. comm.), assessing waste disposal during the recovery phase of an AVF eruption (Hayes 2019), considering how long it will take roads to clear during an evacuation (e.g. Wild et al. 2019; Wild A: Quantitative hazard and risk modelling approaches for volcanic crisis management, in preparation), and modelling multi-year economic impacts to the area (Cardwell RJ, McDonald GW, Wotherspoon LM: Simulation of post eruption time variant land use and economic impacts of a volcanic eruption scenario in the Auckland region of New Zealand, submitted).

Four of the DEVORA scenarios in Hayes et al. (2018) include lava flows: Mt. Eden, Birkenhead, Ōtāhuhu, and Rangitoto Island. In these scenarios, lava flow hazard was represented by a temporal series of hand-drawn areal lava flow footprints, following the precedent of hand-drawing an areal footprint set by Deligne et al. (2015). In Hayes et al. (2018), topography maps were used to estimate the steepest gradient at the flow front’s location. The flow front was then drawn to advance a short direction in the direction of the steepest gradient. Once the footprint had been drawn, the volume represented on the map was calculated using the area of the flow’s footprint and an average (i.e. mean) thickness based on analogue lava flows. This process was used until the total volume prescribed by the scenario had been distributed. One of the motivations of our study was to update the hand-drawn lava flow hazard footprints in these scenarios using numerical modelling.


To analyse the influence of the built environment on lava flow modelling, we quantitatively model the lava flow hazard for three of the DEVORA Scenarios (Mt Eden, Birkenhead, and Ōtāhuhu; Table 1) on both a DEM and a DSM. Lava flow modelling of the Rangitoto Island scenario was not undertaken because the scenario vent is within 25 m of the Hauraki Gulf, minimising the possibility of built environment impacts. Although the vent in the Ōtāhuhu scenario also lies at the edge of a tidal area, the closest body of water is very shallow.

Table 1 Table summarising the rationale and scenario storyline for the three DEVORA scenarios being used here as case studies, after Hayes et al. (2018)

Lava flow hazard model selection

As our modelling results are intended to replace the hand-drawn hazard footprints presented in Hayes et al. (2018), our model selection was guided by the lava flow characteristics included in the original scenarios together with DEVORA stakeholder and research programme guidelines in place for selecting volcanic hazard models for the Auckland Volcanic Field.

The selection of the most appropriate lava flow model(s) is vital to ensure the necessary outputs for the intended application. Hayes et al. (2018) includes quantitative descriptions of the lava flow’s presence (areal footprint) and thickness. Seven published and peer-reviewed models produce both of these outputs: LAVASIM, MAGFLOW, MOLASSES, MULTIFLOW, SCIARA, and VOLCFLOW. Additionally, the DEVORA research team use the following criteria to select hazard models used for the AVF:

  • Ideally hazard models should be explored that can be run in-house in a reasonably short timeframe should an eruption begin.

  • Codes should ideally be open access, i.e. not dependent on a paid software license. If necessary, models that are dependent on proprietary software (e.g. ArcGIS) can be considered if the software is widely available to researchers and stakeholders.

  • Since we are modelling processes in a highly topographically modified and densely populated area, codes that can handle data with high spatial resolutions are preferred.

  • Preference is given to codes that do not need to be significantly modified to be applied to Auckland context.

  • The code(s) (with minimal or no modification) must be executable on a desktop computer.

  • Finally, there is a preference for codes that have already been validated on natural case studies.

Based on these guidelines, we selected the MOLASSES lava flow simulation code (Gallant et al. 2018). MOLASSES determines an areal footprint and the thickness of the lava flow at given locations, yielding the required outputs for our case studies, and can be run locally on both a DSM and DEM without requiring a proprietary license. For a more detailed explanation of how and why MOLASSES was selected, see Tsang (2020).

Lava flow hazard modelling

We have run MOLASSES six times here to illustrate how results can differ if a surface model is used instead of an elevation model. All three scenarios were run twice in MOLASSES, once on a DEM and once on a DSM. The DEM and DSM used in this project are publicly available from Land Information New Zealand (LINZ) and were created using the same LAS dataset captured by Auckland Council during the second half of 2013. They both have a resolution of 1 m and a vertical accuracy of ±0.2 m and a horizontal accuracy of ±0.6 m (both at a 95% confidence level). The same source DEM and DSM were used for all model runs although the DEM and DSM were cropped to the area of interest in each scenario. No special processing was applied to either product; thus, the DSM includes all surficial features, including trees. In our DSM, many of the buildings in the localised areas of interest are used as housing and generally vary from approximately the lava flow’s mean thickness to many metres above the lava flow’s upper crust. The input parameters for the MOLASSES modelling are included in Table 2. Since the lava flow modelling results presented in this paper effectively update the qualitative footprints in Hayes et al. (2018), the lava flow hazard outputs presented should be referred to as the DEVORA scenarios, version 1.1.

Table 2 Input parameters for the MOLASSES lava flow modelling for three DEVORA scenarios. Coordinates are provided in New Zealand Transverse Mercator (2000)

Quantifying the influence of surficial features on lava flow modelling results

We are unaware of published work presenting lava flow modelling using a DSM. The closest equivalent is Charbonnier et al. (2018)‘s study that models lahars. Charbonnier et al. (2018) found that when the lahar model FLO-2D is calibrated with geological data, a modelled lahar can be highly sensitive to sub-meter obstacles. We present lava flow modelling in an urban area and provide a quantitative measure of how the built environment alters the lava flow model results. There are locations where a DSM is unnecessary, and indeed its inclusion could needlessly increase processing times. For example, a DSM is unlikely to be required if the area modelled coincides with a park; indeed, a DSM of a forested park would include combustible trees which would be treated as permanent barriers, which likely does not reflect reality. In order to determine whether the results from a DSM or DEM are more appropriate to use given the wider hazard assessment context of the modelling, we propose using two measures to compare the DEM and DSM results. The first measure is a modified version of the Jaccard coefficient:

$$ Overlap\ \left(\%\right)=100\ast \frac{A_{DEM}\cap {A}_{DSM}}{A_{DEM}\cup {A}_{DSM}} $$

where ADEM is the area covered by the lava flow when modelled on the DEM and ADSM is the area covered by the lava flow when modelled on the DSM. When the resulting overlap is high, then the DSM is not adding value to the simulations run on the DEM. Conversely, when the overlap is low, the implication is that human modifications to the land are resulting in a difference in the modelled lava flow areal footprint. Thus, the lower the percent overlap, the more important it is to consider what is causing the difference between the two footprints. If structures included in the DSM are likely to behave as barriers (e.g. concrete buildings), then the modeller should critically consider if modelling on a DSM is possible. The second measure compares the maximum run out lengths of the footprints:

$$ Length\ \left(\%\right)=100\ast \frac{L_{DSM}}{L_{DEM}} $$

where LDSM is the runout length of the DSM results in the direction of the maximum runout length of the DEM results and LDEM is the maximum runout length of the DEM results. The lower the percentage is, the less overlap of the results there is in that direction.


The lava flow model inputs are arranged by scenario in Table 2 while the model outputs are visually shown in Fig. 2 and compared in Table 3.

Fig. 2
figure 2

Results of lava flow modelling using MOLASSES on a DEM (results shown in pink) and a DSM (results outlined in white) with scenario vent indicated with a red triangle for a) the Mt Eden scenario. b the Ōtāhuhu scenario. c the first lava flow in the Birkenhead scenario and d) the second lava flow in the Birkenhead scenario. Note different scales in each figure. Percentages refer to the modified Jaccard coefficient calculated for each scenario

Table 3 Comparison of the DEM and DSM modelling results, including median and thicknesses in metres, the Jaccard coefficient, and the runout length. Thicknesses have been rounded to the nearest quarter metre


Lava flow hazard models and computing power have greatly improved over recent decades, allowing for more detailed lava flow hazard assessments to be conducted on a variety of surface models.

Importance and implications of surface model selection

Although the selection of a lava flow method is important, choosing the correct elevation model for the project also greatly influences the accuracy of the model’s outputs. When modelling hazards, it is standard to use a DEM (e.g. Barca et al. 1993; Scifoni et al. 2010; Favalli et al. 2012). Occasionally, the DEM may be modified to incorporate a manmade barrier to determine if the barrier could protect an area of significance (Fujita et al. 2008; Scifoni et al. 2010). Although using a modified DEM is not common practice, the topography in some volcanic areas has been extensively modified (e.g. Charbonnier et al. 2018; Hayes et al. 2018). Thus, DEMs may not be the most accurate representation of the obstacles a lava flow may encounter (Kereszturi et al. 2012), and a DSM may be more appropriate.

In our AVF case studies, we presented lava flow hazard footprints on both DEMs and DSMs. In each case, the hazard footprint generated on the DEM had a longer runout length (≥2% to 24% longer) than the corresponding DSM results (Fig. 2). The DSM results also frequently had a corresponding equal or thinner lava flow thicknesses the DEM results (Table 3). The shorter runout length and thinner footprints in the DSM results both have implications for evacuation zones and expected building damage. Shorter runout lengths may necessitate smaller evacuation or exclusions zones. However, thicker lava flows are more likely to overtop single or double storied buildings, increasing expected damage. Additionally, thicker lava flows would have a larger surface area in contact with large structures, which in turn would increase the pressure exerted on the structure by the lava flow (likely increasing expected damage) and potentially increase combustion rates.

However, lava flows frequently damage the built environment they encounter (e.g. Harris 2015; Neal et al. 2019), so obstacles such as timber-framed buildings, represented in the DSM (Fig. 3), may be destroyed during the eruption and no longer affect the lava flow path (e.g. Williams and Moore 1983; Jenkins et al. 2017). Conversely, natural, topographic obstacles are unlikely to disappear during an eruption. Therefore, for lava flow modelling, a surface model that combines a DEM with aspects of a DSM (i.e. a DSM that does not include trees, or buildings where the primary building material has a low combustion point) may be the most accurate representation of a modified terrain. The original LAS dataset would be needed to develop such a hybrid DSM, and the endeavour is likely to be time consuming and to require familiarity with local construction practices. Ideally, lava flow models would allow for changes of the surface model from one iteration to the next to represent progressive damage to the built environment. This will require further research about lava flow impacts. Thus, we suggest that DSMs should be strongly considered when simulating lava flows advancing through heavily modified areas. In lightly modified areas and in areas where the dominant construction material has a low combustion point (i.e. timber-framed buildings), a DEM could still be more appropriate, especially if the final areal footprint is important.

Fig. 3
figure 3

A closer look at MOLASSES lava flow modelling results from a portion of the Birkenhead flow 1 from Fig. 2c, illustrating how the DSM results differ from the DEM results (in blue). White filled in areas indicate obstacles that the lava would potentially flow around, according to the MOLASSES run on the DSM. It is important to note that not all of the buildings would necessarily survive the eruption as many are residential timber-framed structures that have relatively low combustion points

Using a DSM can have important planning implications, but using a DSM also comes with a computational cost. DSMs are frequently produced using much higher resolution data than DEMS as surficial features tend to be on metre scales. Most lava flow models were not created with metre scale surface models in mind (Mossoux et al. 2016). Thus, processing times are often substantially longer when modelling on a high spatial resolution surface model than the low-resolution surface models traditionally used. Depending on the model used and computing power available, rapidly generating lava flow model outputs, as is required during a crisis, may not be possible on a DSM due to the higher resolution required of the topography model (or a DEM with equivalent spatial resolution). Despite these drawbacks, DSMs may be more appropriate input data in certain circumstances, such as when planning exercises.

Quantitative evaluation of surface model selection

Previous applications of the Jaccard coefficient in the context of hazard modelling has been used to assess how accurately a hazard model’s output matches the scenario the model is attempting to replicate (Cordonnier et al. 2016; Mossoux et al. 2016; Favalli et al. 2012). We have adapted it to analyse the overlap between the DEM areal footprint and the DSM areal footprint. Following the lead of Cordonnier et al. (2016) and Dietterich et al. (2017), we do not provide guideline Jaccard coefficient percentages to determine when to use a DSM versus a DEM, as the appropriate cut-off percentage undoubtedly varies from one project to the next.

The Jaccard coefficient indiscriminately calculates the overlap between two areas. If one is modelling on a DSM, it is important to note that the size of the kipuka (any area that is not inundated but is surrounded by a lava flow) in the results. If a building survives an eruption (forming an anthropogenic kipuka (e.g. the kipuka formed by a church during the eruption that began at Paricutin in 1943 (Luhr and Simkin 1993))) but is surrounded by lava, it will be difficult to access. Thus, some may consider such buildings to be useless. A kipuka can also be large enough to enclose multiple buildings, though (Neal et al. 2019). In such cases, the buildings may not sustain any damage although would still be difficult to access during and after the crisis. Both natural and man-made diversion structures could likely create such kipuka. Additionally, buildings surrounded by lava that survive an eruption and can still be entered may preserve smaller assets that can be evacuated when the lava flow is safe to cross. No matter the kipuka size, the kipuka will influence the Jaccard coefficient results (e.g. Fig. 2). For the aforementioned reasons, we feel that the Jaccard coefficient calculation should include the DSM areal footprint as is, i.e. with the kipuka created by the buildings. While including the diversion structure and building kipuka will decrease the total overlap compared to using only the perimeter of the areal footprints to calculate the Jaccard coefficient, the kipuka have planning implications. In some instances such as evacuation planning, the outer perimeter of the areal footprint is of more interest than the entire footprint.

Comparison of results with qualitative lava flow hazard footprints in Hayes et al. (2018)

Most of the input parameters used by Hayes et al. (2018) are the same as those used here. In fact, Hayes et al. (2018) used a surface model, rather than an elevation model, to draw their qualitative lava flow footprints. The volume pulse is the only input parameter that is not shared between Hayes et al. (2018) and this paper. Despite very similar starting points, the resulting areal extents of the lava flows in Hayes et al. (2018) and this paper are quite different. First, the extents in Hayes et al. (2018) are notably longer (Fig. 4). This is likely because effusion rates, represented by the pulse volume parameter in MOLASSES, are not directly considered in Hayes et al. (2018). Discontinuous lava effusion was considered in Hayes et al. (2018) and represented by the stalling of the lava flow fronts; variations in effusion rates were, otherwise, not considered. In order to create MOLASSES footprints with similar extents to Hayes et al. (2018), a lower pulse volume would be necessary. Second, the surface models used in this paper had a higher spatial resolution than the surface model used by Hayes et al. (2018) enabling more precise topographic gradient calculations than were used in Hayes et al. (2018). Assuming accurate representation of the physical process, the paths lava followed in our case studies (Fig. 2; Fig. 4b) are more realistic than those drawn in Hayes et al. (2018; Fig. 4a). Although both the effusion rate and spatial resolution discrepancies suggest that the areal extents presented here consider more nuanced data than those in Hayes et al. (2018), inputs to the MOLASSES modelling undertaken here could be improved. For example, we used the AVF lava flows’ mean thickness, rather than the requested modal thickness, due to availability of data, but this assumes the AVF lava flow thickness distribution is Gaussian, which we deem unlikely. Thus, the lava flow modelling in this paper can be further updated as we learn more about AVF lava flows.

Fig. 4
figure 4

Comparison of the Mt Eden lava flows created a) for Hayes et al. (2018) where the lava flow effuses from the bottom of the cone (a large area represented by the brown circle) and b) using MOLASSES where the lava flow effuses from a point source vent (orange star)

The case studies presented here provide more information (i.e. thickness maps on a metre scale) on the lava flows in the DEVORA scenarios than previously existed, so our model results can be used in a variety of preparedness activities including stakeholder workshops and evacuation and contingency modelling. The results of modelling conducted in this paper are aligned with Hayes et al. (2018)‘s goal of using quantitative models for all of the considered hazards, although we acknowledge that qualitative methods may be more appropriate depending on the required outputs.


As settlements continue to encroach upon volcanoes, it is important to start or continue to prepare for future eruptions. Lava flow modelling has become increasingly sophisticated over the past few decades, but with the proliferation of models, selection of the most suitable model for a specific purpose becomes more difficult. Although most lava flow models have been developed to run on DEMs, DSMs capture all obstacles (irrespective of their physical properties) the flow may encounter. Thus, when lava flow hazard models are run on a DSM, the results will have considered the obstacles encountered by the flow. We illustrate the difference a DSM can make in three case studies in the AVF. In our case study scenarios, we found that lava flows modelled on a DEM generally have a longer runout length compared to lava flows modelled on a DSM. We use a modified Jaccard coefficient to quantify the difference between using a DEM and using a DSM. The most overlap (90%) occurred in the Ōtāhuhu scenario which is located in an area with minimal anthropogenic modifications. Thus, a DEM adequately portrays the topography in this scenario. Since there is much less overlap in terrestrial areas with generally high levels of urban development (e.g. the Birkenhead and Mt Eden scenarios), we would recommend using the DSM results in these cases. A similar decision would likely be made in most of the terrestrial Auckland Volcanic Field.

Availability of data and materials

The dataset supporting the conclusions of this article is available in the University of Auckland’s Figshare repository, doi: in



Auckland Volcanic Field


Digital Elevation Model


Determining Volcanic Risk in Auckland


Digital Surface Model


Light intensity detection and ranging


  • Allen SR, Smith IEM (1994) Eruption styles and volcanic hazard in the Auckland volcanic field, New Zealand. Geosci Reports Shizuoka Uni 20:5–14

    Google Scholar 

  • Axelsson P (1999) Processing of laser scanner data- algorithms and applications. J Photogrammetry Remote Sens 54:138–147.

    Article  Google Scholar 

  • Barca D, Crisci GM, Di Gregorio S, Nicoletta F (1993) Cellular automata methods for modelling lava flows: simulation of the 1986-1987 eruption, Mount Etna, Sicily. In: Kilburn CRJ, Luongo G (eds) Active lavas: monitoring and Modelling. Routledge, New York

    Google Scholar 

  • Bebbington MS (2015) Spatio-volumetric hazard estimation in the Auckland volcanic field. Bull Volcanol 77.

  • Bebbington MS, Cronin SJ (2011) Spatio-temporal hazard estimation in the Auckland volcanic field, New Zealand, with a new event-order model. Bull Volcanol 73:55–72

    Article  Google Scholar 

  • Bilotta G, Cappello A, Herault A, del Negro C (2019) Influence of topographic data uncertainties and model resolution on the numerical simulation of lava flows. Environ Model Softw 112:1–15

    Article  Google Scholar 

  • Bonaccorsco A, Calvari S, Boschi E (2016) Hazard mitigation and crisis management during major flank eruptions at Etna volcano: reporting on real experience.

  • Bonne K, Kervyn M, Cascone L, Njome S, Van Ranst E, Suh E et al (2008) A new approach to assess long-term lava flow hazard and risk using GIS and low-cost remote sensing: the case of Mount Cameroon, West Africa. Int J Remote Sens 29:6539–6564.

    Article  Google Scholar 

  • Brown SK, Jenkins SF, Sparks RSJ, Odbert H, Auker MR (2017) Volcanic fatalities database: analysis of volcanic threat with distance and victim classification. J App Volcanol 6.

  • Brunsdon D, Park B (2009) Lifeline vulnerability to volcanic eruption: Learnings from a national simulation exercise. Technical council on lifeline earthquake engineering conference, Oakland, p 2009

    Google Scholar 

  • CDEM (2015) Auckland volcanic field contingency plan. Auckland Emergency Management.

    Google Scholar 

  • Charbonnier SJ, Connor CB, Connor LJ, Sheridan MF, Oliva Hernández JP, Richardson JA (2018) Modeling the October 2005 lahars at Panabaj (Guatemala). Bull Volcanol 80.

  • Chirico GD, Favalli M, Papale P, Boschi E, Pareschi MT (2009) Lava flow Hazard at Nyiragongo volcano, D.R.C. 2. Hazard reduction in urban areas. Bull Volcanol 71:375–387.

    Article  Google Scholar 

  • Connor JL, Connor CB, Meliksetian K, Savov I (2012) Probabilistic approach to modelling lava flow inundation: a lava flow hazard assessment for a nuclear facility in Armenia. J App Volcanol 1:3

    Article  Google Scholar 

  • Cordonnier B, Lev E, Garel F (2016) Benchmarking lava-flow models. In: Harris AJL, De Groeve T, Garel F, Carn SA (eds) Detecting, Modelling, and responding to effusive eruptions. Geological Society of London, London

    Google Scholar 

  • County of Hawai’i (2019) Kīlauea Eruption Recovery. Accessed 25 Sept 2019

    Google Scholar 

  • Deligne NI, Blake DM, Davies AJ, Grace ES, Hayes JL, Potter S et al (2015) Economics of resilient infrastructure Auckland volcanic scenario. ERI research report 2015/03, p 153

    Google Scholar 

  • Dietrich JT (2015) Riverscape mapping with helicopter-based structure-from-motion photogrammetry. Geomorph 252:144–157.

    Article  Google Scholar 

  • Dietterich HR, Lev E, Chen J, Richardson JA, Cashman KV (2017) Benchmarking computational fluid dynamics models of lava flow simulation for hazard assessment, forecasting, and risk management. J Appl Volcanol 6.

  • Favalli M, Chirico GD, Papale P, Pareschi MT, Boschi E (2009a) Lava flow hazard at Nyiragongo volcano, D.R.C. 1. Model calibration and hazard mapping. Bull Volcanol 71:363–374.

    Article  Google Scholar 

  • Favalli M, Chirico GD, Papale P, Pareschi MT, Coltelli M, Lucana N, Boschi E (2006) Computer simulations of lava flow paths in the town of Goma, Nyiragongo volcano, −Democratic Republic of Congo. J Geophys Res Solid Earth 111.

  • Favalli M, Fornaiciai A, Pareschi MT (2009b) LIDAR strip adjustment: application to volcanic areas. Geomorph 111:123–135.

    Article  Google Scholar 

  • Favalli M, Tarquini S, Papale P, Fornaciai A, Boschi E (2012) Lava flow hazard and risk at Mt. Cameroon volcano. Bull Volcanol 74:423–439.

    Article  Google Scholar 

  • Fujita E, Hidaka M, Goto A, Umino S (2008) Simulations of measures to control lava flows. Bull Volcanol 71:401–408.

    Article  Google Scholar 

  • Gallant E, Richardson J, Connor C, Wetmore P, Connor L (2018) A new approach to probabilistic lava flow hazard assessments, applied to the Idaho National Laboratory, eastern Snake River. Geol 46:895–898.

    Article  Google Scholar 

  • Harris AJL (2015) Volcanic hazards, risks and disasters. In: Shroder JF, Papale P (eds) Volcanic hazards, risks and disasters. Elsevier, Waltham

    Google Scholar 

  • Hayes J, Wilson TM, Deligne NI, Cole J, Hughes M (2017) A model to assess tephra clean-up requirements in urban environments. J App Volcanol 6.

  • Hayes JL (2019) Clean-up and restoration of urban environments after volcanic eruption. Dissertation. University of Canterbury, Christchurch

  • Hayes JL, Tsang SW, Blake DM, Lindsay JM et al (2018) The DEVORA scenarios: multi-hazard eruption scenarios for the Auckland volcanic field. GNS Science, Avalon, p 138.

    Book  Google Scholar 

  • Hayes JL, Wilson TM, Deligne NI, Lindsay JM, Leonard GS, Tsang SWR, Fitzgerald RH (2020) Developing a multi-hazard eruption scenario suite using an interdisciplinary approach. J Volcanol Geotherm Res 392.

  • Heidemann HK (2018) Lidar base specification (ver. 1.3). In: U.S. Geological Survey techniques and methods. USGS, Reston.

    Chapter  Google Scholar 

  • Jaggar TA (1945) Protection of harbors from lava flow. Am J Sci 243-A:333–351

    Google Scholar 

  • Jenkins SF, Day SJ, Faria BE, Fonseca JFBD (2017) Damage from lava flows: insights from the 2014-2015 eruption of Fogo, Cape Verde. J App Volcanol 6.

  • Kereszturi G, Cappello A, Ganci G, Procter J, Németh K, del Negro C, Cronin SJ (2014) Numerical simulation of basaltic lava flows in the Auckland volcanic field, New Zealand- Implication for volcanic hazard assessment. Bull Volcanol 76.

  • Kereszturi G, Procter J, Cronin SJ, Nemeth K, Bebbington M, Lindsay J (2012) LiDAR-based quantification of lava flow susceptibility in the City of Auckland (New Zealand). Remote Sens Environ 125:198–213.

    Article  Google Scholar 

  • Kermode LO (1992) Geology of the Auckland urban area 1: 50 000. GNS Science geological map 2(1). 1 sheet + 63p

    Google Scholar 

  • Le Corvec N, Bebbington MS, Lindsay JM, McGee LE (2013) Age, distance, and geochemical evolution within a monogenetic volcanic field: analyzing patterns in the Auckland volcanic field eruption sequence. Geochem Geophys Geosys 14:3648–3665.

    Article  Google Scholar 

  • Leonard GS, Calvert AT, Hopkins JL, Wilson CJN, Smid ER, Lindsay JM et al (2017) High-precision40Ar/39Ar dating of quaternary basalts from Auckland volcanic field, New Zealand, with implications for eruption rates and paleomagnetic correlations. J Volcanol Geotherm Res 343:60–74.

    Article  Google Scholar 

  • Lindsay J, Marzocchi W, Jolly G, Constantinescu R, Selva J, Sandri L (2010) Towards real-time eruption forecasting in the Auckland volcanic field: application of BET_EF during the new Zealand National Disaster Exercise ‘Ruaumoko’. Bull Volcanol 72:185–204.

    Article  Google Scholar 

  • Luhr JF, Simkin T (1993) Paricutin: the volcano born in a Mexican cornfield. Geoscience Press, Inc, Phoenix

  • McGee LE, Millet M-A, Beier C, Smith IE, Lindsay JM (2015) Mantle heterogeneity controls on small-volume basaltic volcanism. Geology 43(6):551–554.

    Article  Google Scholar 

  • McPhee J (1989) The control of nature. Farrar, Straus and Giroux, New York

    Google Scholar 

  • Morgan AV (2000) The Eldfell eruption, Heimaey, Iceland: a 25-year retrospective. Geosci Canada 27:11–18

  • Mossoux S, Bartolini S, Poppe S, Canters F, Kervyn M (2016) Q-LavHA : a flexible GIS plugin to simulate lava flows. Comput Geosci 97:98–109.

    Article  Google Scholar 

  • Neal CA, Brantley SR, Antolik L, Babb JL, Burgess M, Calles K et al (2019) The 2018 rift eruption and summit collapse of Kīlauea volcano. Science 363:367–374.

    Article  Google Scholar 

  • Ouedraogo MM, Degre A, Debouche C, Lisein J (2014) The evaluation of unmanned aerial system-based photogrammetry and terrestrial laser scanning to generate DEMs of agricultural watersheds. Geomorph 214:339–355.

    Article  Google Scholar 

  • Rongo R, Lupiano V, Spataro W, D’ambrosio D, Iovine G, Crisci GM (2016) SCIARA: cellular automata lava flow modelling and applications in hazard prediction and mitigation. In: Harris AJL, De Groeve T, Garel F, Carn SA (eds) Detecting, Modelling, and responding to effusive eruptions. Geological Society of London, London

    Google Scholar 

  • Runge MG (2015) Improving probabilistic hazard forecasts in volcanic fields Application to Harrat Rahat, Kingdom of Saudi Arabia. Dissertation. University of Auckland, Auckland

  • Sandri L, Jolly G, Lindsay J, Howe T, Marzocchi W (2012) Combining long- and short-term probabilistic volcanic hazard assessment with cost-benefit analysis to support decision making in a volcanic crisis from the Auckland volcanic field, New Zealand. Bull Volcanol 74:705–723.

    Article  Google Scholar 

  • Scifoni S, Coltelli M, Marsella M, Proietti C, Napoleoni Q, Vicari A, Del Negro C (2010) Mitigation of lava flow invasion hazard through optimized barrier configuration aided by numerical simulation: the case of the 2001 Etna eruption. J Volcanol Geotherm Res 192:16–26.

    Article  Google Scholar 

  • Searle EJ (1964) City of volcanoes: a geology of Auckland. Paul Longman, Auckland

    Google Scholar 

  • Tsang SW, Lindsay JM, Deligne NI (2019) Short-term preparation for and response to an impending lava flow: lessons from the June 27th lava flow (2014-2015), Hawaii, USA. GNS Science, Avalon, p 87.

    Book  Google Scholar 

  • Tsang SWR (2020) Modelling the hazard footprint and consequences of lava flows in an urban environment. Dissertation. University of Auckland, Auckland

  • UNDDR (2019) Global assessment report on disaster risk reduction. United Nations Office for Disaster Risk Reduction, Geneva

    Google Scholar 

  • Villa F (2002) Mount Etna huffs, puffs, but causes no injuries; lava flow is far from settlements. The Record, Bergen County New Jersey

    Google Scholar 

  • von Hochstetter F (1859) The isthmus of Auckland with its extinct volcanoes. T. Delattre, Auckland

    Google Scholar 

  • Wagner W, Ullrich A, Ducic V, Melzer T, Studnicka N (2006) Gaussian decomposition and calibration of a novel small-footprint full-waveform digitising airborne laser scanner. J Photogrammetry Remote Sens 60:100–112.

    Article  Google Scholar 

  • Westoby MJ, Brasington J, Glasser NF, Hambrey MJ, Reynolds JM (2012) ‘Structure-from-motion’ photogrammetry: a low-cost, effective tool for geoscience applications. Geomorph 179:300–314.

    Article  Google Scholar 

  • White SA, Wang Y (2003) Utilizing DEMs derived from LIDAR data to analyse morphologic change in the North Carolina coastline. Remote Sens Env 85:39–47 10.1016.S0034-4257(02)00185-2

    Article  Google Scholar 

  • Wild A, Lindsay JM, Costello SB, Deligne NI, Doherty A, Leonard GS, Maxwell K, Rollin J, Wilson TM (2019) Auckland volcanic field eruption crisis management decision-making pilot workshop, vol 44. GNS Science, Avalon.

    Book  Google Scholar 

  • Williams RS, Moore JG (1983) Man against volcano: the eruption on Heimaey, Vestmennaeyjar, Iceland. US Geological Survey, Reston

    Google Scholar 

Download references


We would like to thank Laura Connor for her help setting up MOLASSES and Dr. Graham Leonard for discussions of this work. We are grateful to Prof. Chuck Connor and an anonymous reviewer for helpful comments on this manuscript.


We are grateful for funding from the Determining Volcanic Risk in Auckland (DEVORA) research programme, the New Zealand Earthquake Commission (EQC), and the GNS Science Strategic Investment Fund.

Author information

Authors and Affiliations



Authors SWRT and JML conceived this study. SWRT undertook the lava flow modelling. GC advised SWRT on the lava flow modelling. All authors contributed to the preparation and editing of the manuscript. The author(s) read and approved the final manuscript.

Corresponding author

Correspondence to Sophia W. R. Tsang.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

We declare we do not have any competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

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 The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Tsang, S.W.R., Lindsay, J.M., Coco, G. et al. The influence of surficial features in lava flow modelling. J Appl. Volcanol. 9, 6 (2020).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: