CC BY
15
RUSSIAN JOURNAL OF EARTH SCIENCES, VOL. 17, ES5003, doi:10.2205/2017ES000612, 2017
Improving the system of numerical simulation of volcanic ash propagation using the PUFF model
Sergey I. Malkovsky1, Aleksei A. Sorokin1, and Sergey P. Korolev1
Received 11 December 2017; accepted 11 December 2017; published 30 December 2017.
The ash clouds and plumes produced during explosive volcano eruptions are extremely hazardous for the population and economic activities and, for this reason, their detection and control over their movements are considered an urgent and vital problem. The solution of the problem is primarily linked to the satellite data analysis and modeling volcanic ash cloud movement trajectories. One of the most accepted models widely used in these investigations is the PUFF model and its Puff-UAF software implementation. The experience of applying the model in different regions of the world for more than a decade to investigate volcano explosive eruptions has shown a sufficient correlation established between the forecasted and observed volcanic ash cloud movement trajectories obtained using remote sensing methods. The paper describes the existing approaches and program solutions applied to meteorological dataset preprocessing for numerical simulation of the ash cloud movement and dispersion from explosive volcano eruptions using the Puff-UAF software suite. The detected problems of wgrib2-generated NetCDF file usage for the puff program operation, as well as their solutions, are discussed. The test numerical simulation results derived using different Puff-UAF versions, the original Puff-UAF 2.2.2 and the Puff-UAF with the authors' puff source code modifications are presented. Their comparison with remote sensing data obtained for the Shiveluch volcano eruption is made and the correctness of the puff source code modifications is concluded. KEYWORDS: Volcano monitoring; forecasting; computational models; algorithms; software tools and services; remote sensing.
Citation: Malkovsky, Sergey I., Aleksei A. Sorokin, and Sergey P. Korolev (2017), Improving the system of numerical simulation of volcanic ash propagation using the PUFF model, Russ. J. Earth. Sci., 17, ES5003, doi:10.2205/2017ES000612.
1. Introduction
PUFF is a numerical volcanic ash-tracking model capable of forecasting the ash plume location after an explosive volcano eruption. Tanaka [1994] has originally proposed the model based on the three-dimensional La-grangian formulation of the pollutant dispersion. Searcy et al. [1998] (the University of Alaska, Fairbanks, UAF) have developed the software implementation of the PUFF model in C++, a high-level computer programming language. This software implementation has become the basis for Puff-UAF [Peterson R. A., 2003, Puff UAF User's Manual, http://puff.images.alaska.edu/doc.shtml (accessed August 20, 2017)], which is a freely distributed software suite applied for forecasting the airborne ash cloud movement and dispersion and the ash fallout numerical simulation. The
1Computing Center of the Far Eastern Branch, Russian Academy of Sciences, Khabarovsk, Russia
Copyright 2017 by the Geophysical Center RAS. http://elpub.wdcb.ru/journals/rjes/doi/2017ES000612-res.html
aforementioned open source software (GNU General Public License version 2) is available for download from this website: https://sourceforge.net/projects/puff-uaf/ (actual version 2.2.2).
PUFF is widely used for operational forecasting of the ash plume trajectories from explosive volcano eruptions [Scollo et al., 2009; Webley et al., 2008, 2009], and for studying the previous events in order to determine or refine their characteristics [Aloisi et al., 2002; Daniele et al., 2009; Fero et al., 2008; Tanaka and Yamamoto, 2002].
Applying Puff-UAF for numerical simulation of the Kamchatka volcano ash cloud movement and dispersion, we detected the incorrectness of the results obtained by its puff component when loading meteorological data in the NetCDF file format [Rew and Davis, 1990], prepared using the wgrib2 utility (http://www.cpc.ncep.noaa.gov/products/wesley/wgrib2/ index.html).
The incorrectness implies inconsistency between the forecasted and observed volcanic ash cloud movement trajectories. It appeared to be caused by incomplete support of the NetCDF files, whose time variable is expressed in seconds.
ES5003
1 of 6
Figure 1. An overview of a conceptual scheme illustrating the Puff-UAF software operation.
When loading the data, the puff consistently uses only part of these, which, eventually, influences the performed calculation results.
The paper briefly describes the meteorological data preprocessing algorithm using the wgrib2 utility and our Puff-UAF source code modifications that allow for avoiding the detected shortcomings and ensure the implementation of the correct program operation on the similar NetCDF files.
2. Puff-UAF Architecture
The Puff-UAF software source code is distributed with a set of scripts, assigned for its compilation. After the Puff-UAF compilation and installation, a set of executable files including puff, ashdump, ashgmt and ashxp is generated. The puff program provides the PUFF model implementation and numerical simulation of the ash cloud movement and dispersion. The ashdump displays the simulation results (latitude, longitude, height, size and age of individual ash particles) as a plain text. The ashgmt and ashxp utilities allow us to visualize the calculation results in the PNG and JPEG graphic file formats, and the ashgmt utility enables to create the GIF animation of the simulated ash plume trajectories.
Figure 1 illustrates a conceptual schematic overview of the Puff-UAF software operation. The primary data for the puff calculation include the meteorological data files in the NetCDF format, the program options (output mode, verbosity level, etc.) and the modeled event parameters (the volcano name and location, the eruption origin time and duration, the initial ash cloud height, etc.). After the puff program execution, a set of the NetCDF files, containing the information on the ash particle position, size and age or ash concentration, is generated.
3. Meteorological Data Preprocessing for Puff-UAF Operation
The basic input data for computer calculations are the zonal and meridional components of wind velocity (u- and ^-components) at a set of pressure levels, covering the modeling area. The auxiliary data may include the temperature fields at the same pressure levels (T, used at variable fall dynamics) and geopotential height data (Z). These meteorological data are transferred to the program as a NetCDF file of a 4D structure (time, level, latitude and longitude).
To produce these NetCDF files, both the analysis and forecast data can be applied that are derived from various numerical weather forecast models, including the Global Forecasting System (GFS) [Environmental Modeling Center, 2003], the Integrated Forecasting System of the European Centre for Medium-range Weather Forecasts (IFS ECMWF) [Persson, 2015], etc. These forecast data are usually distributed in the GRIB format and therefore require preprocessing prior to their use.
The meteorological data files can be prepared in a variety of ways. The Puff-UAF manual [Peterson R. A., 2003, Puff UAF User's Manual, http://puff.images.alaska.edu/doc.shtml (accessed August 20, 2017)] suggests using the gribtonc utility (http://www.unidata.ucar.edu/software/decoders/), which supports files of GRIB Edition 1 format (grib), and generates a NetCDF file of the required structure on output. However, the numerical weather forecast model products are nowadays usually distributed in GRIB Edition 2 format (grib2). Therefore, if using gribtonc, then the downloaded meteorological data first need to be converted from grib2 to grib format, applying, for example, the cnvgrib utility (http://www.nco.ncep.noaa.gov/pmb/codes/GRIB2/) (Figure 2). This procedure involves the additional technical operations that lead to an increase in the total time needed for the input data preprocessing, and results in sufficient disk space usage to store temporary files. For instance,
Figure 2. Conversion of grib2 files into NetCDF format using the cnvgrib and gribtonc utilities.
Table 1. Resource Usage
Utility Conversion Input file
time (sec) volume (MB)
cnvgrib + gribtonc 22.6 193
wgrib2 4.4 152
erations or leading to the incorrect simulation results. However, in terms of using Puff-UAF for automated batch processing, it is more important to reduce the number of the auxiliary data preprocessing steps and resources allocated. Therefore, the puff runtime execution has been analyzed in detail to reveal the causes of the previously detected errors.
when preparing the input NetCDF file on a computer system with Intel Xeon Processor E5-2650 v2 2.6 GHz and Intel solid-state drive DC S3500 Series of 120 GB volume to simulate the event considered in Section 5 (eight grib2 files obtained using the GFS weather forecast model: 0.5° global grid, 26 pressure levels; u and v - the horizontal wind components, T - temperature and Z - geopotential height), the conversion time needed for the gribtonc and cnvgrib utilities exceeds that for the wgrib2 utility by a factor of 5. The storage space usage by employing these utilities attains 21% exceeding that using the wgrib2 utility (see Table 1). Besides, the development and support of the gribtonc utility, included in NetCDF Decoders package, have been terminated. Taking all the above into consideration, using the wgrib2 utility becomes more convenient.
Two problems have been detected during puff calculations performed with wgrib2-converted input files. The first problem is that in some cases, the simulation results appear to be incorrect. The second problem is concerned with the erroneous simulation results when using geopotential height values contained in the input files.
To summarize, several meteorological data preprocessing approaches exist for numerical simulation of the ash cloud movement and dispersion using Puff-UAF. Each of these approaches has its own shortcomings either associated with performing the additional resource-consuming technical op-
4. The Problem Analysis and the Suggested Solutions
During the investigation of the first problem, we have clarified that the cause lies in the difference associated with time representation in the NetCDF files prepared in the two aforementioned ways. In wgrib2-generated NetCDF files, the "time" coordinate variable, which stores the time dimension of each record variable, has the values expressed in the POSIX format. This implies that it contains the number of seconds that have elapsed since the start of the UNIX epoch (00:00:00 UTC, January 1, 1970). At the same time, in the NetCDF files, prepared using the method recommended in the official Puff-UAF manual, the corresponding coordinate variable has the values expressed as the number of hours elapsed since January 1, 1992. Further investigation has shown that the puff incorrectly processes the coordinate variable values expressed in seconds.
The algorithm of loading the "time" coordinate variable values from the NetCDF files to the puff program is as follows (Figure 3, the prefix "nc" denotes the NetCDF file, whereas "puff" stands for the program):
1. Determine the number of the "time" coordinate variable values in the NetCDF file (corresponds to the number of records).
Figure 3. An overview of the scheme of work with the "time" coordinate variable.
2. Load the "time" coordinate variable values to the float-type internal array using the "as_float" method of the NcValues class of the NetCDF library [http://www.unidata.ucar.edu/software/netcdf/docu mentation/historic/netcdf/ (accessed August 20, 2017)]. If the coordinate variable values are expressed in seconds, then convert them to hours (because the puff has all time values stored in hours).
3. Remove the out-of-simulation-time-range values from the array.
4. Use the resultant array values to identify the indices of the necessary NetCDF file records and then read them.
As seen from the algorithm description, the "time" variable values in the POSIX format are incorrectly processed because a reverse conversion from hours to seconds is not performed during the identification of the needed file record indices. This has resulted in a situation when the program is always loading the variable values from the NetCDF file that refer to the very first record. Likewise, because the time values are loaded using the "as_float" method of the NcValues class, the loss of accuracy and incorrect program operation may occur in case of loading the time value greater than 224. This is due to floating point number representation in the IEEE 754 standard. The loss of accuracy does not occur when the NetCDF files have the time expressed in hours because of small time values.
If the time in the file is stored in the POSIX format, then the correct program operation is ensured (after correction of the record index selecting process) for the dates ranging between 19:39:44 UTC June 20, 1969 and 04:20:16 UTC July 14, 1970 only.
In order to resolve the aforestated problems, we have modified the part of the puff source code related to the NetCDF file reading, whose timings of records are stored in the POSIX format. The gist of the introduced modifications is as follows:
1. At the second step, the temporary "double" variable has been introduced that stores the time values loaded using the "as_double" method of the NcValues class of the NetCDF library. If the time value is expressed in seconds, then it is converted to hours. Once this is done, the variable value is assigned to the appropriate element of the array.
2. At the third step, when indices of the loaded records are selected, a reverse conversion from hours to seconds is performed, if the time in source NetCDF files is expressed in seconds.
The second problem appears to be caused due to the NetCDF file generation by the wgrib2 utility in the format compatible with the GrADS program (http://cola.gmu.edu/grads/). Using "m" instead of "gp m" in the "units" attribute value of the variable storing the geopotential height data is one of the features of this format. These values cannot be correctly processed by the puff utility which leads to erroneous determination of the pressure level heights, and, consequently, to the invalid simulation
Table 2. Main Characteristics of the Shiveluch Volcano Eruption
Date and time Duration Minimum Maximum
height of the height of the
ash cloud ash cloud
2015.03.03, 80 min 8 km 12 km
22:50 UTC
results. To prevent incorrect processing of the values, we have modified the corresponding utility source code. The patch changes the logic of processing the "units" attribute values of variables in the "PtoH" and "pressureGridFromZ" methods of the puff's Grid class.
The implemented puff source code modifications have provided a complete support of the NetCDF files created with the wgrib2 utility.
5. A Comparative Analysis of the Obtained Results
To verify the correctness of the puff source code modifications, the experiments on numerical simulation of the ash cloud movement and dispersion have been performed. The parameters of the explosive event occurred at one of Kamchatka volcanoes, for which the remote sensing data (RSD) were available providing a comparative analysis of the obtained results, were used as the basic input data for performing numerical calculations.
The study investigates one of the events occurred at Shiveluch volcano during 1-6 March 2015 (VONA/KVERT WEEKLY INFORMATION RELEASE 10-2015, http://www.kscnet.ru/ivs/kvert/van/?n=10-2015) [Gordeev and Girina, 2014]. The main parameters of the event are listed in Table 2. We used the parameters obtained from Advisory No. 2015/73 of the Tokyo Volcanic Ash Advisory Center (VAAC) of the Japan Meteorological Agency (http://ds.data.jma.go.jp/svd/vaac/data/TextData/2015/ 20150303_SHEV_0073_Text.html), and those derived from the information release issued by the Kamchatka Branch of the Federal Research Center "Geophysical Survey of the Russian Academy of Sciences" on 3 March, 2015 (http://www.emsd.ru/ ssl/monitoring/arhiv/2015/Mar/03_ Mar.htm). Because the volcano was not visible due to cloudiness during the eruption, all the eruption parameters were recovered from seismological data. This is the reason why the listed data on the initial ash cloud height are broadly ranging.
The numerical forecast products generated from the 18:00 UTC, March 3, 2015 run of the GFS weather forecast model (0.5°global grid, 26 pressure levels; u and v - the horizontal wind components, T - temperature and Z - geopotential height) have been used for numerical simulation as the input meteorological data (http://www.nco.ncep.noaa.gov/pmb/ products/gfs/). The input NetCDF file prepared using the
Figure 4. A comparison of the numerical simulation results on dynamics of the ash cloud movement obtained using the original (a) and modified (b) puff utility versions (color scale showing a height range of ash particle occurrence) and the satellite data.
wgrib2 utility contained 3 to 24 hours with 3-hour time step forecast records. Using the cnvgrib and gribtonc utilities, the reference NetCDF file has also been prepared that contained the same data and was used for control.
To perform numerical simulation, the following parameters have been used: -seed 129 -volc shiveluch -sedimentation reynolds -dtMins 1 -saveHours 373min -eruptHours 1.33 -runHours 7 -nAsh 2000 -plumeMin 8000 -plumeMax 12000 -eruptDate "2015 03 03 22:50" -model gfs. In order to obtain the ash plume location for the time that corresponds to the time of the satellite image capture, the "-saveHours" parameter was assumed to be equal to 373 minutes.
To estimate the numerical simulation results, the remote sensing data recorded by the AVHRR sensor operating on board of NOAA-18 satellite on 4 March 2015 (05:03 UTC) and stored in the Information System "VolSatView" were used [Gordeev et al., 2016; Sorokin et al., 2016a]. Figure 4a and Figure 4b illustrates the results of the satellite data processing (a difference of the 11 and 12 micron channels) combined with the numerical simulation results obtained using the original (a) and modified (b) puff utility versions making
use of wgrib2-generated NetCDF input files. The ash cloud (white color) and its outline (solid red line) constructed by the experts of the Kamchatka Volcanic Eruption Response Team (KVERT) using "VolSatView" and the simulation results (dots) can be seen in the figures.
As can be seen from the Figure 4, the results obtained using the modified puff utility version are more consistent with the RSD data as compared to the results obtained by the original puff utility version employing the similar simulation parameters and meteorological datasets. These results are also consistent with those derived from the original puff utility version using the control NetCDF file. This gives grounds for validation of the correctness of the puff source code modifications.
The complementary experiments on validation of the modified Puff-UAF version for explosive events occurred at various volcanoes were performed in the normal mode of operation of the Information System "VolSatView". They have also shown good agreement between the modeling results obtained using the NetCDF input files generated by the wgrib2 utility and the factual data on volcanic ash cloud movements.
The modified Puff-UAF software source code is available at the website: http://lits.ccfebras.ru/assets/files/puff-uaf-2.2.2-1.tar.gz.
All calculations performed within the framework of the study were carried out in the Shared Facility Center "Data Center of FEB RAS" (Khabarovsk) [Sorokin et al., 2016b]
Conclusion
The modifications of the source code of puff component of the Puff-UAF software suite provide a complete support of the input NetCDF files, prepared using the wgrib2 utility. This allows for avoiding a range of technical operations and reducing the total time expenditure for meteorological dataset preparation and disk space usage necessary for their storage.
This also enables using the high resolution meteorological data for Puff-UAF simulations that are obtained from the GFS weather forecast model (0.25 degree by latitude and longitude, 1-hour time step) and are not supported by the cnvgrib and gribtonc utilities to convert from grib2 format. Improving the accuracy and precision of the performed calculations becomes possible.
The above characteristics are critically important for monitoring volcanic activity which requires high accuracy, precision and operability during the performance of numerical calculations. Using the modified Puff-UAF software suite, the Information Systems "Signal" and "VolSatView" provides new possibilities to create the tools capable of simulating the ash cloud movement trajectories and to perform a joint analysis of the obtained results and the remote sensing data [Gordeev et al., 2016].
Acknowledgments. This work was funded by the Russian Science Foundation (grant no. 16-17-00042). The authors used the computation systems and data storage technologies created within
the framework of the projects of the Russian Foundation for Basic Research (RFBR) (grants nos. 15-29-07953 and 16-37-00026) and the Far Eastern Branch of the Russian Academy of Sciences (grant no. 15-I-4-072) .
References
Aloisi, M., M. D'Agostino, K. Dean, A. Mostaccio, G. Neri (2002), Satellite analysis and PUFF simulation of
the eruptive cloud generated by the Mount Etna paroxysm of 22 July 1998, J. Geophys. Res., 107, No. B12, 2373, doi:10.1029/2001JB000630 Daniele, P., L. Lirer, P. Petrosino, N. Spinelli, R. Peterson (2009), Applications of the PUFF model to forecasts of
volcanic clouds dispersal from Etna and Vesuvio, Comp, and Geosci, 35, No. 5, 1035-1049, doi:10.1016/j.cageo.2008.06.002 Environmental Modeling Center (2003), GFS atmospheric
model. NCEP Office Note 442, 14 pp., Global Climate and Weather Modeling Branch, EMC, Camp Springs, Maryland. Fero, J., S. N. Carey, J. T. Merrill (2008), Simulation of the 1980 eruption of Mount St. Helens using the ash-tracking model PUFF, J. Volcan. Geotherm. Res., 175, No. 3, 355366, doi:10.1016/j.jvolgeores.2008.03.029 Gordeev, E. I., O. A. Girina (2014), Volcanoes and
their hazard to aviation, Her. Russ. Acad. Sci., 84, 1, doi:10.1134/S1019331614010079 Gordeev, E. I., et al. (2016), The VolSatView information system for Monitoring the Volcanic Activity in Kamchatka and on the Kuril Islands, J. Volcanol. Seismol., 10, No. 6, 382394, doi:10.1134/S074204631606004X Persson, A. (2015), User guide to ECMWF forecast products (ECMWF, version 1.2), 129 pp., ECMWF, United Kingdom. Rew, R., G. Davis (1990), NetCDF: An interface for
scientific data access, IEEE Comp. Graph. and Appl., 10, No. 4, 76-82, doi:10.1109/38.56302 Scollo, S., M. Prestifilippo, G. Spata, M. D'Agostino, M. Coltelli
(2009), Monitoring and forecasting Etna volcanic plumes, Nat. Hazards Earth Syst. Sci., 9, No. 5, 1573-1585, doi:10.5194/nhess-9-1573-2009 Searcy, C., K. Dean, W. Stringer (1998), PUFF: a highresolution volcanic ash tracking model, J. Volcanol. Geotherm. Res., 80, No. 1/2, 1-16, doi:10.1016/S0377-0273(97)00037-1 Sorokin, A. A., S. P. Korolev, O. A. Girina, I. V. Balashov, V. Yu. Efremov, I. M. Romanova, S. I. Malkovsky (2016a), The integrated software platform for a comprehensive analysis of ash plume propagation from explosive eruptions of Kamchatka volcanoes, Current problems in remote sensing of the Earth from space, 13, No. 4, 9-19, (in Russian)doi:10.21046/ 2070-7401-2016-13-12-9-19 Sorokin, A. A., S. P. Korolev, A. N. Polyakov (2016b), Development of Information Technologies for Storage of Data of Instrumental Observation Networks of the Far Eastern Branch of the Russian Academy of Sciences, Procedia Computer Science. Elsevier Science B.V., 66, 584-591. Tanaka, H. L., K. Yamamoto (2002), Numerical simulation of volcanic plume dispersal from Usu volcano in Japan on 31 March 2000 using PUFF model, Earth Planets Space, 54, No. 7, 743-752, doi: 10.1186/BF03351727 Tanaka, H. L. (1994), Development of a prediction scheme for volcanic ash fall from Redoubt volcano, Alaska, Proc. First Int. Symp. on Volcanic Ash and Aviation Safety p. 283-291, Geol. Survey Bull. 2047, Seattle, WA, U.S. Webley, P. W., et al. (2008), Predicting and validating the tracking of a volcanic ash cloud during the 2006 eruption of Mt. Augustine volcano, Bull. Amer. Meteorol. Soc., 89, No. 11, 1647-1658, doi:10.1175/2008BAMS2579.1 Webley, P. W., K. Dean, J. E. Bailey, J. Dehn, R. Peterson (2009), Automated forecasting of volcanic ash dispersion utilizing Virtual Globes, Nat. Haz., 51, No. 2, 345-361, doi:10.1007/s11069-008-9246-2
Sergey P. Korolev, Sergey I. Malkovsky, and Aleksei A. Sorokin, Computing Center of the Far Eastern Branch, Russian Academy of Sciences, Khabarovsk, 680000, Russia. ([email protected])
