The following table contains the list of all contributions accepted for CMG 2010 and the day of their presentation. Abstracts can be viewed as plain text or downloaded as Portable Document Format (.pdf) by clicking on the displayed links in the Title and PDF columns. Please note that some misspellings could be present in the database and the automated pdfs. Revised versions will be published as soon as possible.

All database fields can be browsed by using the "QuickSearch" function (please notice that pressing the Enter button will clear the form). To search also INSIDE the abstracts text, please check the appropriate box into the Search Settings menu.

QuickSearch:   Number of matching entries: 166/166.

Search Settings

AuthorsTitleSessionSchedulePresentation typePDF
Abe, S., Mair, K. & Urai, J.L. Roughness Evolution in DEM simulations of fault shear. Session 4 Monday 7  oral file 
Abstract: The surfaces of most known faults are not flat planes but show roughness at all scales. Due to many field measurements, mainly using LIDAR, the structure of this roughness is relatively well known. However, the mechanisms which generate the fault roughness are not yet fully understood. Fault slip processes are difficult to observe directly in nature and laboratory experiments designed to investigate the details of such processes have some limitations. We therefore use a numerical model to augment the field and laboratory observations. Due to the inherently discrete nature of the investigated processes, which are dominated by brittle fracturing and frictional interactions, we use a Discrete Element Method (DEM) model to investigate the process of fault surface abrasion. The simulation models consist of two blocks of solid material (made of bonded particles that can break apart) with a fault in between which is initially smooth at the macroscopic scale. A constant normal stress is applied to the outer surfaces of the model and the upper and lower blocks are moved at a constant speed so that the fault is sheared. At the edges of the model, in the shear direction, periodic boundary conditions are applied so that the amount of shear displacement is not limited by the model size. During fault shearing particles break away from the intact fault surfaces so that the roughness of the surfaces evolves and a fault gouge develops. Our initial investigation of the roughness spectra of the abraded surfaces produced shows that for sufficiently large displacement and normal stress the spectra evolve towards a power-law relationship. Furthermore, we observe a strong anisotropy of the roughness parallel and perpendicular to the slip direction. In some models slip-parallel striations resembling those commonly seen in natural faults are clearly recognisable. The development of these striations, their influence on mechanical properties of the fault and the role of fault gouge on roughness development are now being investigated.
Aharonov, E., Katsman, R. & Laronne Ben-Itzhak, L. How do stylolite networks and stylolite-fracture networks form: insights from modeling Session 2 Tuesday 8  poster file 
Abstract: Stylolites and solution seams are surfaces of localized dissolution, commonly observed in many sedimentary rocks, especially in carbonates. Many stylolites, from diverse origins, locations, and tectonic settings, share a common morphological feature: the amount of dissolved rock in the center of the seam is proportional to the stylolite length. However, kilometer-long stylolites with uniform dissolution along the seam are also observed. In addition we observe anastomosing networks of stylolites, and systems of fractures and veins interconnected with the stylolites. The question addressed here is how to explain these different morphologies and the conditions for their formation. The mechanism of stylolite formation has long been under intensive discussion, and it is currently thought to consist of two coupled processes: pressure solution and clay-enhanced dissolution. Our numerical elasto-plastic spring-network model, that couples these two processes, indeed reproduces growth of all the types of stylolites observed in the field: single self-similar stylolites, infinitely long parallel stylolites, anastomosing network of stylolites, and systems of fractures and veins interconnected with the stylolites. Our results suggest that the initial distribution of clay (bedding planes or dispersed defects) dictates which of the network geometries will form. To further understand how the interconnected systems of stylolites and fractures form, we apply the Eshelby Transformation method and obtain an analytical solution revealing the tensile/shear stress concentrations and predicting the angles of brittle fracturing in the vicinity of stylolites.
Ailleres, L., Lindsay, M., Jessell, M. & deKemp, E. Geological uncertainty: a mean to reduce potential field ambiguity? Session 7 Tuesday 8  oral file 
Abstract: We have identified the need for a better integration of geological into 3D geological modelling and especially during potential field inversions. To be reliable, 3D models need to reconcile all relevant & available geological and geophysical data. Current practices involve building a 3D geological model from geological and some geophysical data and assume that the geometries and more importantly the topology of the model produced are correct. This reference model is then inverted against potential field data and in the end no check is proposed with respect to the input data. The uncertainty attached to geological data is not taken into account. In this paper, we investigate the idea of combining geological uncertainty with potential field data ambiguity in order to produce 3D geological models with a reduced geophysical ambiguity and more importantly an estimated geological uncertainty. We propose to calculate multiple reference geological models using geological data variability; some of which derives from varying models explaining the structural evolution of the considered area (with implications on the model topology). We need to emphasize that as structural geologists, we can utilise structural information (such as fold axes, axial surfaces, cleavages, etc.) to predict where a given lithology should be in 3D space, therefore predicting the spatial variation of a given rock property (density or magnetic susceptibility). These structural elements derive directly from the considered structural evolution model. In other words, we are testing ``what if" scenarios as proposed in the literature against geological datasets and consequently geophysical datasets during potential field inversions. The robustness of the models is going to be estimated with respect to a (low) deviation of the modelled geology from the input data. This requires the development of geological penalty functions that will estimate this deviation and the potential variability of the models. Our proposed process allows for the production of a 3D variability map similar to a sensitivity analysis of the input data on the geometry/topology predicted by the models. This process helps investigate a larger part of the parameter space defined by geology (geometry of rock packages), petrophysical properties and geophysics during inversion. We will present some of the preliminary uncertainty cubes produced out of simple scenarios.
Albarello, D., D'Amico, V. & Mucciarelli, M. Empirical Testing of Probabilistic Seismic Hazard Estimates Session 8 Thursday 10  poster file 
Abstract: Several probabilistic procedures are presently available for seismic hazard assessment (PSHA). These result in a number of different outcomes (hazard maps), each generally compatible with available observations and supported by plausible physical models. To take into account this inherent uncertainty (``epistemic"), outcomes of these alternative procedures are combined in the frame of logic-tree approaches by scoring each procedure as a function of the respective reliability. This is deduced by evaluating ex-ante (by expert judgements) each element concurring in the relevant PSH computational procedure. This approach appears unsatisfactory also because the ``value" of each procedure depends both on the reliability of each concurring element and on that of their combination: thus, checking the correctness of single elements does not allow evaluating the correctness of the procedure as a whole. An alternative approach to scoring is here presented that is based on the ex-post empirical testing of the considered PSH computational models. This is performed by comparing the probabilistic ``forecasts" provided by each model with empirical evidence relative to seismic occurrences (e.g., strong-motion data or macroseismic intensity evaluations) during some selected control periods of dimension comparable with the relevant exposure time. In order to take into account the inherent probabilistic character of hazard estimates, formally coherent procedures have been developed that are based on Likelihood estimates and Counting protocols. Some results obtained by the application of these testing procedures in Italy will be shortly outlined.
Angheluta, L., Mathiesen, J. & Aharonov, E. Teething stylolites and the ache of an unknown instability Session 2   withdrawn file 
Abstract: Pressure solution is a major ductile deformation mechanism operating in sedimentary rocks, with practical implications to oil, gas, and water flow, as well as to controlling rates of basic geological processes such as rock compaction, folding, and healing between earthquakes. Pressure solution occurs as a mass transfer and chemical compaction process, where material dissolves in regions of high stress, diffuses through the fluid-saturated pore space and precipitates in regions of low stress. This process is often observed to occur in a localized manner in `stylolites' - irregular interfaces that fluctuate around a plane, with the stylolite `teeth' pointing in the direction of maximum compressive stress. The formation of stylolites has been a longstanding issue and up to now there has been no consensus on the underlying mechanisms. Here we consider two dynamical instabilities which both can generate irregular solution interfaces. The first is based on the accumulation of residual insoluble minerals (i.e. clays) at the interface of the stylolites. We assume that the clay plays an important role in pressure solution by enhancing the local solubility, as previously suggested by independent experiments and observations. That is, there are two different equilibrium concentrations, one for the solute in the bulk and one where the pore fluid is in contact with the clay layer. We show that this situation initiates a reactive transport away from the clay layer which then becomes morphologically unstable. The second instability we consider, appears when there are gradients in material properties. Such gradients give rise to corresponding gradients in the chemical potential and in turn drive the roughening process. Both instabilities cause roughening of an initially flat dissolution surface, possibly explaining growth of stylolite teeth.
Asef, M.R., Farrokhrouz, M. & Haghi, H. Evaluation of P-wave and S-wave Correlations Session 4 Tuesday 8  poster file 
Abstract: Although a number of equations have been suggested for prediction of S-wave without direct measurement, but these prediction equations can be divided into three main groups. The main predicting equations define Vs from Vp which are usually used in rock mechanics and petroleum industry with limited number of Vs measurements. The correlation between such data points are used as a based for the other set in nearby geological sequences. The second group of estimations defines Vs based on other logging data such as resistivity logs, gamma ray, NPHI and etc. Sometimes a mixture of these logging data is used. This method is basically used in petroleum industry where sonic logging costs more than other logging operations and the reason for such predictions are economical concerns. The final series of predicting equations correlate Vp or Vs with some physical properties like porosity, density, quarts content, cementation degree and etc. The purpose of such equations is to identify effective parameters in soil or rock velocity and measure their degree of influence rather than regressions between shear and compressional velocities. Based on this grouping, a number of different correlations were gathered and a new series of data from deep petroleum reservoir was evaluated with this correlations. It was decided to evaluate the degree of accuracy of these equations in new data set and define uncertainty limits. The main shortcoming of previous suggested equations was limitation to rock types. Main prediction equations were suggesting for sandstone which has a simple and predictable framework and grain structure and predicting equations are usually good estimators; but when it comes to other sedimentary rocks like limestone and shale or metamorphic and igneous rocks, major discrepancies occurred in predicted values. It was decided to evaluate predicted equations based on rock types and define the range of accuracy for new data set. The results showed that based on rock origins, predicted equations have better results for sedimentary rocks compared with igneous and metamorphic rocks. Between sedimentary rock types, sandstone showed better results of estimation in comparison with limestone and shale. It seems that fractures and secondary porosity in limestone and bedding and fissility in shale are the main causes of uncertainty in wave velocity predictions. In metamorphic or igneous rocks, porosity seems to be the most important parameter which prevents good estimation.
Ashkenazy, Y., Eisenman, I., Gildor, H. & Tziperman, E. The effect of Milankovitch variations in insolation on equatorial seasonality Session 3 Tuesday 8  poster file 
Abstract: Although the sun crosses the equator twice a year, at the equinoxes, at times in the past the equatorial insolation has had only one maximum and one minimum throughout the seasonal cycle due to Milankovitch orbital variations. Here we use a state-of-the-art coupled atmosphere-ocean general circulation model to study the effect of such insolation forcing on equatorial surface properties including air and sea temperature, salinity, winds, and currents. We show that the equatorial seasonality is altered according to the insolation, with, for example, either maximal sea surface temperature (SST) close to the vernal equinox and minimum SST close to the autumnal equinox, or vice versa. Our results may have important implications for understanding tropical climate, as well as for the interpretation of proxy data collected from equatorial regions.
Ashkenazy, Y., Paldor, N. & Zarmi, Y. On the meridional structure of extra-tropical Rossby waves Session 3 Tuesday 8  poster file 
Abstract: The most common derivation of Rossby waves is based on the quasi-geostrophic approximation, which is usually assumed to be valid for planetary meridional scales of the order of 1000 kilometers. Several decades-long numerical simulations using a general circulation model are used to demonstrate that initial conditions of the quasi-geostrophic solution do not propagate westward with a uniform phase speed, which implies that these harmonic solutions are not eigensolutions of the shallow water equations in channels of width exceeding several hundred kilometers. It is shown that the latitudinal variation of the phase speed of these numerical solutions results from the fact that the eigensolutions are the parabolic cylinder functions and not sinusoidal functions, where the higher modes are slower and extend further poleward. To explain these numerical results a simple non-harmonic approximation for extratropical Rossby waves on the sphere is proposed, where the meridional coordinate turns out to be a parameter instead of as a continuous variable. It is shown that, in contrast to the quasi-geostrophic solution, the meridional structure of Rossby waves is irrelevant to a first order approximation. The proposed approximation closely reproduces the numerical results and captures the latitudinal dependence of the Rossby-waves-phase speed, when starting from arbitrary initial meridional structure. We propose a similar approximation for topographic Rossby waves that accounts for the wave structure and the phase speed for arbitrary depth profile in one direction and cross bathymetry initial conditions. In all these cases the theory yields explicit expressions both for the dispersion relations and the non-harmonic wave structure.
Aspinall, W. Evaluating uncertainty and probabilistic forecast skill in volcanic hazard and risk assessments Session 7 Thursday 10  poster file 
Abstract: The Soufriere Hills volcano in Montserrat has been erupting, with timevarying levels of activity and intensity, for nearly 15 years. Since late 1997, a consistent risk estimation approach has been used in more than 21 successive biannual meetings, with risk calculated from probabilistic forecasts of event scenarios and population numbers. The outcomes of these forecasts are amenable to quantitative analysis, the first time such an appraisal has been undertaken for scientific advice during a volcanic crisis. The assessment cycle involves forecasts over the following 6 or 12 months. We cannot predict sensu stricto the exact timing, locations and magnitudes of dangerous volcanic events at the Soufriere Hills volcano, except in very exceptional circumstances, and certainly not for years ahead. Nevertheless the capabilities of the Montserrat Volcano Observatory have improved with effort and experience, and computer models of hazardous processes have advanced during the eruption. Confirming that these improvements have resulted in better forecasts is difficult because their effects cannot be separated easily from the temporal and intensity variations in activity and other changeable elements of the hazard/risk assessment process. Expert judgments concerning potential hazardous scenarios and their uncertainties are the main source of information that go into the risk assessment process, reflecting experience of similar ruptions elsewhere whilst conflating this in a Bayesian manner with the properties of the Montserrat volcano. Important scenarios for assigning risk relate to vulnerable populated areas, involving incursion of pyroclastic flows or directed lateral blasts; neither have happened to an unevacuated area, so these forecasts cannot be tested directly. However, 110 scenario probabilities can be examined against actual outcomes using the Brier Skill Score. Of these, 75 might be termed safetycritical had the areas been occupied, and for 83% of these a positive forecast skill score was achieved. Thus the experts have outperformed an ''uninformed'' baseline probabilistic forecast in the majority of cases, underlining the value for decision support for civil protection. This forecast skill can be expressed also in terms of an equivalent financial investment gain, which make the benefits of the scientists' forecast performance more readily understandable to decision makers, and the public generally.
Avdev, S.N. & Tzankov, C.V. Presenting of Earth's gravity field with optimized models of point masses Session 5 Thursday 10  poster file 
Abstract: The task of the report is to present the most important results obtained by calculation of an Earth's gravity field model by using optimized point masses. The used gravity acceleration data solves the inverse problem by seeking field sources in Earth's core, mantle and crust. The main idea involved in the report is through gravimetrically equivalent distribution of optimized point masses to build a stable Earth's density model and by calculation of its potential to determinate the surfaces of Earth's geoid, core and other geospheres. The proposed approach for solving the problem gives the possibility for complete exhaustion of Earth's internal structure information contained in the used gravimetrical data sample. The method is tested on $5^rctimes 5^rc$ grid data of absolute gravity acceleration values which includes 2592 data points. The geoid heights and absolute values for g used in the model are taken from the GRACE Gravity Model 02 (GGM02C) released October 29, 2004 and published to the public on As a base the mean-tide system model calculated to degree 360 is taken. A number of more and more involved geoid models are constructed by founding of unique and sustainable solutions for gradually increased number of point masses. Comparing the geoids with the base model indicates that with increasing the complexity of the point models, the differences gradually lowering. The underlying purpose is to show that the average of the absolute deviations between the base geoid and the point masses geoid become negligible when the secondary model consist of several hundred point masses.
Barsotti, S., Coltelli, M., Costa, A., Folch, A., Macedonio, G., Nannipieri, L., Neri, A., Prestifilippo, M., Scollo, S., Spata, G. & Boschi, E. Investigating the model-dependent uncertainty of volcanic ash dispersal forecasts Session 7 Thursday 10  poster file 
Abstract: Accurate modelling of volcanic ash dispersal and deposition is a challenging goal for volcanology. In particular, forecasts of tephra dispersal have a key role in the assessment of the associated hazard and impact mitigation. In order to produce robust forecasts of this phenomenon in case of explosive events occurring at Mt. Etna (Sicily, Italy), the Istituto Nazionale di Geofisica e Vulcanologia has set up an automatic web-based procedure based on a multi-model simulation approach. Each day five numerical models, named FALL3D, HAZMAP, PUFF, TEPHRA and VOL-CALPUFF, produce forecasting maps of aerial ash concentration and tephra deposit generated by hypothetical explosive events at Mt. Etna. By adopting the same initial meteorological dataset provided by an external mesoscale model, the numerical codes produce, for each eruptive scenario assumed, time-varying maps of ash concentration at different vertical atmospheric levels, maps of ash concentration vertically integrated over the domain, and maps of cumulative ground deposition. Model outcomes are then integrated to produce synthetic maps to be used as early-warning information for aviation and airport operations. Such a multi-model approach allows to quantify the model-dependent uncertainty characterizing the modelling capability. In fact, the five numerical codes, based on different modelling approach (e.g. Eulerian vs Lagrangian), different physical formulations of the process, as well as different numerical algorithms, allow us to obtain independent representations of the phenomenon. Quantitative comparisons between different model results were performed estimating key variables of ash dispersal and deposit. Some of the most important variables are the plume height, the main dispersal axis of the cloud, the ash cloud aerial extension and shape, as well as the tephra ground deposition at specific locations. Comparisons were extended over a time period of several months and allowed analytical comparisons between the models as a function of the eruptive and meteorological conditions assumed. Statistical analyses were also carried out to discriminate between casual random and systematic differences between these models. Such differences were then interpreted in terms of the specific features of each model so to gain crucial insights on the strengths and weaknesses of different models.
Benites, R.A. & Ben-Zion, Y. Seismic wavefield generated by shear dislocations in a structure with a bimaterial interface and near-by cracks Session 4 Tuesday 8  poster file 
Abstract: We develop a computational framework to calculate seismic wave propagation in a damaged fault zone structure, consisting of two welded quarter spaces of different elastic properties with various distributions of cracks in the slower medium near the interface. The incident wave corresponds to a strike-slip shear dislocation, and the total wavefield is calculated at observation points along lines perpendicular to the strike of the fault. The full waveform simulations incorporate the multiple scattering among the cracks, a vertical bimaterial interface and the free-surface. We use a Boundary Integral numerical scheme based on artificial wave-source distribution around the boundaries of all cracks $[$1$]$, for which both the displacement and stress Greenś functions due to each source are computed exactly using a kernel analytical solution derived by the Cagniard-de Hoop method $[$2$]$. Examined case studies include several realizations of crack distributions, crack sizes and aspect ratios. The results show head and body waves, expected from previous related model of a bimaterial interface separating two homogeneous quarter spaces, and scattered wavefield generated by the cracks. For wavelengths much larger than the cracks size, the presence of cracks reduces the average wave speeds and increase the seismic attenuation compared to the properties of the hosting rocks.\br> $[$1$]$ R. Benites, et al., J. Acoust. Soc. Am., 86 (1992).\ $[$2$]$ Y. Ben-Zion, Geophys. J. Int., 98 (1989). \
Bergantz, G.W. & Burgisser, A. Challenges in model validation of multiphase mixtures: scaling requirements for application to volcanic and magmatic systems Session 8 Friday 11  oral file 
Abstract: Validation of models for volcanic/magmatic systems requires assessment of the appropriate kinematic and dynamic scaling of multiphase mixtures. We present scaling relations that provides a more complete framework for scaling dilute magmatic systems by explicit consideration of the complexity caused by the feedback between particles (crystal, bubble, or pyroclast) and the continuous phase (liquid or gas). We consider three canonical igneous multiphase systems: magma chambers, volcanic plumes, and pyroclastic surges, and we provide estimates of the proposed scaling relations for published experiments on those systems.
Dilute magmatic mixtures can display a range of distinct dynamical regimes that we characterize with a combination of average (Eulerian) properties and instantaneous (Lagrangian) variables. The Eulerian properties of the mixture yields the Reynolds number (Re), which indicates the level of unsteadiness in the continuous phase. The Lagrangian acceleration of particles is a function of the viscous drag and gravity forces, and from these two forces are derived the Stokes number (St) and the stability number (Rt), two dimensionless numbers that describe the dynamic behavior of the particles within the mixture. The compilation of 17 experimental studies relevant for pyroclastic surges and volcanic plumes indicates that there is a need for experiments above the mixing transition (Re $>$ 104) and for scaling St and Rt. Among the particle dynamic regimes present in surges and plumes, some deserve special attention, such as the role of mesoscale structures on transport and sedimentary processes or the consequences of the transition to turbulence on particle gathering and dispersal.
The compilation of seven experimental studies relevant to magma bodies indicates that, in the laminar regime, crystals mostly follow the motion of the melt, and thus the physical state of the system can be approximated an ensemble single phase. In the transition to turbulence, magmas can feature spatially heterogeneous distributions of laminar regions and important velocity gradients. This heterogeneity has a strong potential for crystal sorting. In conclusion, the Re-St-Rt framework demonstrates that, despite numerous experimental studies on processes relevant to magmatic systems, some and perhaps most, geologically important parameter ranges still need to be addressed at the laboratory scale.
Bisson, M., Sulpizio, R., Zanchetta, G., Demi, F. & Santacroce, R. Assessment of collapsing colluvial cover volume as crucial input for volcaniclastic flow modelling in the Circumvesuvian area Session 1 Thursday 10  poster file 
Abstract: The volcaniclastic flows are among the most recurrent and dangerous natural hazards in volcanic terrains characterized by recurrent pyroclastic deposition. Their recurrence usually increases during and shortly after an eruptive event but they can also be generated by storms or earthquakes during volcanic quiescence. Damage and casualties relate to these events are common in the history of recent (AD 1631-1944) activity of Vesuvius, both near the volcanic complex and in the sub-Apennine basins affected by pyroclastic fallout. An example is the catastrophic event occurred on 5 and 6 May 1998 in the area of Sarno (Southern Italy), which highlights the destructive potential of debris flows, even when they are of relatively low magnitude. More than 130 people were killed and severe property damage took place when volcaniclastic debris flows triggered by heavy rainfall inundated various towns located in piedmont areas. With the aim to mitigate and prevent the volcaniclastic flow hazard, several models (physical or empirical) are proposed to delineate the areas affected by volcaniclastic flow-inundation and different methodologies are currently applied to individuate the source basins more prone to disruption for triggering volcaniclastic flows. Whatever is the model used, the more important input parameter is the volume of remobilised deposit. With this aim, we present here a method for the assessment of the volume of collapsing material under different rainfall conditions and thickness of the colluvial cover using the software SINMAP. The proposed procedure allows the calculation of volume of collapsing material for each drainage basin, and therefore to assess the initial volume feeding volcaniclastic flows. SINMAP allows to classify the terrain stability basing on the infinite slope stability model with wetness (pore pressures) obtained from a topographically based steady state model of hydrology.
Bobrov, A.M. & Baranov, A.A. The role of variable viscosity in the modeling of global stress fields in the Earth's mantle and in floating continents Session 9 Tuesday 8  poster file 
Abstract: In numerical two-dimensional experiments we investigate the evolution of spatial field of viscous overlithostatic horizontal stresses in the mantle and in moving continent. A continent moves self-consistently with the time-varying mantle forces which act from viscous mantle. A continent is modeled by the set of numerous markers, which transfer the high viscosity area. We compare two model laws for viscosity: the first, simplest variant - isoviscous mantle case; and the second, exponential p,T-dependent viscosity case. For these two models we analyze how a form of viscosity law in the mantle can change the horizontal stress fields in the mantle and continent and the stress evolution in the course of movement of a continent. Inclusion of variable viscosity gives the possibility to take into account, in given framework of purely viscous statement, the presence of oceanic lithosphere and the difference in viscosity values of upper and lower mantle. It is found at 10e7 Rayleigh number that the horizontal tensile stresses observed at the initial stages (above ascending mantle flow) are succeeded by the compressive stresses, especially at the leading continental edge (up to 40 Mpa when the continent thrusts onto subduction zone). It is obtained also that the distribution of horizontal stresses along the continent distinctly points to the locations of upward and downward flows in the subcontinental mantle. The stresses in the major portion of the mantle are in the range plus minus 5 MPa (50 bars). This work was supported by the Russian Foundation for Basic Research, Project No. 09-05-00319-a.
Brazier, R.A. & Boomer, K. Empirically based seismic location criteria for simple to complex geological structures: A jacknife approach for local/regional networks Session 7 Thursday 10  poster file 
Abstract: Determining accurate seismic locations with representative uncertainty estimates is of fundamental importance to ground-based nuclear explosion monitoring, including the assignment of accurate ground truth (GT) levels. The monitoring community relies on selection criteria for classifying seismic events at the GT5 level, which specifies the absolute location and depth errors as being less than 5 kilometers. A new approach to obtaining empirically based (EB) criteria for estimating GT levels of seismic events recorded on a regional network has been developed using a jackknife resampling method applied to carefully picked phase arrival times for GT reference events.
A criterion for the Archaean Kaapvaal craton using reference events from several South African gold mines and recorded locally on the Southern African Seismic Experiment (SASE) network has been developed. This criteria developed on a relatively simple crustal structure is compared to a preliminary criterion in the more complex region of Ethiopia and the global criterion. In addition we explore transporting criteria to similar geological structures such as the Archean Tanzanian craton where no reference event is available.
Brazier, R. & Boomer, K. Stochastic velocity Modeling of the Kaapvaal Craton beneath the BASO station. Session 7 Thursday 10  poster file 
Abstract: From a statistical view point, geophysical parameters are interpreted as random variables. For example, arrival times for P and S waves have a random component, attributed to sources such as distant-dependent measurement error, human error in picking arrival phases, and random noise attributed to travel path. Stochastic modeling of such phenomenon quantifies these random processes. Genetic algorithms are iterative stochastic models that evaluate progressively improved mathematical models. In this presentation, we will provide a background on how genetic algorithms in general, and the Non-dominated Sorting Genetic Algorithm (NSGA-II) in particular, model the velocity structure, with specific reference to the BOSA station on the Kaapvaal Craton in Southern Africa.
Brovkin, V. Limitations of global terrestrial biosphere models used for future climate projections Special Session Thursday 10  oral file
Abstract: Processes in terrestrial ecosystems, to large extent, are controlled by changes in climate and CO2 concentration. In turn, geographical distribution of ecosystems and plant physiology strongly affect heat, moisture, and momentum fluxes between land surface and atmosphere. These interactions form a feedback loop between terrestrial biosphere and climate which modulates substantially the Earth system dynamics on different time scales. Experiments with coupled atmosphere-vegetation models suggest that these interactions could be strong enough to support multiple steady states in the climate-vegetation system. However, how robust are the models of terrestrial biosphere, and how good are constraints on the terrestrial biosphere processes?
Intercomparison of land surface models included into the Earth system models used for future climate and CO2 projections reveals a remarkable dissimilarity among different models. One of key uncertainties is our limited understanding of the ecosystem response to the increase in the atmospheric CO2 concentration. The CO2 fertilization of terrestrial productivity is included into all models; however the availability of nutrients such as nitrogen and phosphate necessary for the biomass buildup is usually neglected. Another example of poorly constrained process is a response of soil organic carbon to elevated temperature, which is especially essential for permafrost regions. A proper representation of the land surface heterogeneity requires running the global models on very high spatial resolution, and this is still beyond the computational capacity of the Earth system models. These and other limitations of current generation of terrestrial biosphere models as well perspectives of the model development will be discussed.
Bykov, P.L. & Gordin, V.A. The Analysis of Three-Dimensional Geometry of Atmospheric Fronts Session 9 Tuesday 8  poster file 
Abstract: The atmospheric fronts' geometry is evaluated on standard baric levels according to data of the objective analysis (or of the forecast) of basic meteorological fields on global regular grids. The traditional predictors of the fronts are 1) vertical component of the wind's vorticity; 2) Laplacian of pressure upon a sea level; 3) temperature gradient. We use (instead of the first ones) maximal eigenvalues of the following matrices: a symmetric part of matrix Jacobi horizontal wind and the second horizontal derivatives (matrix Hesse) geopotential (or pressure upon a sea level), correspondingly. Then we search the optimal weight for these three predictors' combination.
We use high order approximations (Fourier transform and compact scheme approaches).
Heuristic methods of construction of the frontal thin zones along ``crests'' in the predictor's 2D field are developed. Also we use some algorithms for a vertical adjustment of the frontal lines.
A validation of the geometry's evaluation is considered, too. We calculate two kinds of correlation functions (CF) for the basic meteorological fields. CF are functions of the distance between to points. Two clusters of points' pair are used: the points are divided by a frontal zone and points not are divided by a frontal zone. Our algorithm (with optimal weights) provides significant difference (in the integral norm) between these CF.
Thus, we demonstrate an algorithm and codes that construct frontal surfaces with complex geometry in 3D space.
V.A.Gordin. Mathematical Problems and Methods in Hydrodynamical Weather Forecasting. - Amsterdam etc.: Gordon & Breach Publ. House, 2000.
V.A.Gordin. Mathematics, Computer, Weather Forecasting, and Other Mathmatical Physics' Scenarios. - Moscow, Fizmatlit, 2010 (in Russian, in press).
Ph.L.Bykov, V.A.Gordin. The Objective Analysis of Three-Dimensional Geometry of Atmospheric Fronts. Research Activities in Atmospheric and Oceanic Modeling, 2009\br> documents/ 02$$Gordin$$Vladimir$$fronts.pdf P.2.11-2.12.
Campillo, M., Froment, B., Hadziioannou, C., Larose, E., Rivet, D., Roux, P., Brenguier, F. & Shapiro, N. Noise based seismic speed monitoring: co-seismic and post seismic speed drops Session 4 Monday 7  oral file 
Abstract: Continuous recordings of ambient seismic noise can be used to perform a continuous monitoring of the elastic properties of a medium. We present the theoretical background and laboratory tests which demonstrate the feasability of the approach. We then consider different cases including a moderate earthquake, a large earthquake and a silent earthquake. Our results indicate that earthquakes are associated with relative velocity drops of the order of 10$^-3$. The velocity drop extends to a regional scale for large (M8) earthquake. We observe that the velocity change is associated both with strong motion effect on the surficial materials and to the deformation at depth. This last point is illustrated in the case of a silent earthquake. The velocity change seems to be instantaneous for the upper ''elastic'' crust but delayed for the deep crust. These results will be confronted with conceptual models of deformation.
Cardona, Y. & Bracco, A. Mesoscale variability, high frequency winds and their impact on the vertical velocity fields of the South China Sea Session 3 Tuesday 8  poster file 
Abstract: The South China Sea (SCS) is a marginal basin with a complex dynamics influenced by the East Asian Monsoon, river discharge and intricate bathymetry. As a result, both the mesoscale eddy field and the near-inertial energy distribution display large temporal and spatial variability, and they strongly influence the biogeochemistry response of this region. Here, with an ensemble of numerical integrations using ROMS, we investigate how the temporal resolution of the forcing wind field modifies the vertical velocity patterns and impacts vertical mixing. We evaluate the response of the ocean circulation in the SCS under three different forcing conditions: monthly mean, daily mean and six-hourly winds from NCEP/QUICKSCATT for the period 2000 and 2007. While the surface circulation does not display significant differences, the vertical velocity field shows high sensitivity to the frequency of the wind forcing in correspondence to mesoscale vortices and filaments. The high frequency wind energy injected at the surface is transferred at depth through the generation and subsequent straining effect of Vortex Rossby Waves (VRWs) immediately below the surface, and through near-inertial internal oscillations trapped inside anticyclonic vortices below ~ 50m. We analyze the physical mechanism responsible for those changes and the details of the interplay between the near-inertial and mesoscale eddy fields under the various forcing fields, and we quantify the impact of those changes on the vertical mixing. Our results are consistent with theoretical approximations which have been tested over analytical domains.
Carrassi, A. & Vannitsem, S. Accounting for model error in data assimilation. A deterministic Formulation Session 8 Friday 11  oral file 
Abstract: Data assimilation is confronted with the presence of model errors arising from the imperfect description of atmospheric dynamics. Due to the lack of a unique framework for the model error treatment and to the difficulty to accumulate reliable statistics, these errors are usually either ignored or modeled on the basis of simple assumptions such as white noise or first order Markov process.
By using a deterministic formulation, an approach to account for the model error in data assimilation is proposed. This deterministic perspective has been the guideline for model error treatment in sequential and variational schemes. The approach is based on a formal expression for the deterministic evolution of the model error and on the use of an approximation suitable for practical nonlinear filtering problems. The accuracy of this approximation is inherently connected to the stability property of the dynamics such as the spectrum of Lyapunov exponents. The numerical analysis is performed using low order chaotic prototypes of geophysical fluids.
First, the deterministic description of the model error evolution, incorporated into the classical extended Kalman filter equations, is used to estimate the contribution to the forecast error covariance due to model imperfections. Results reveal that substantial improvements of the filter accuracy can be gained as compared with the classical white noise assumption.
A natural step ahead implies the online estimation of the model parameters in conjunction with the system's state. The EKF in the state augmentation formulation is implemented. The dynamical evolution law for the model error allow now for the estimation of the cross covariances between errors in the state estimate and parameters, and both are efficiently evaluated online over a wide range of initial parametric error amplitude.
The extension to variational schemes is presented in the sequel. The state estimation problem requires here the model error time correlations. A straightforward way to estimate these correlations using the deterministic approach is proposed, and it is incorporated in the weak-constraint variational assimilation. Results with two simple dynamics of increased complexity suggest the potential for the successful application to more realistic situations.
Caya, A., Buehner, M. & Carrieres, T. Three-Dimensional Variational Data Assimilation in a coupled Ice-Ocean Model with Ensemble-Derived Background-Error Covariances Session 8 Friday 11  oral file 
Abstract: A three-dimensional variational data assimilation (3D-Var) system has been developed to provide analyses of the ice-ocean state for coupled ice-ocean models. The study focuses on the impact of sea-ice data assimilation on the prediction of sea-ice condition. To obtain ocean state analysis increments consistent with the sea ice analysis increments and in preparation for assimilating sea surface temperature data, an ensemble estimate of the background-error covariance matrix is used and preliminary results are shown. The ice-ocean model is driven by atmospheric wind and temperature forcings. An ensemble of forecast-analysis experiments is produced by perturbing the assimilated observations and using 20 members of the atmospheric ensemble prediction system run operationally at the Canadian Meteorological Centre. Data assimilation experiments, using various configurations of the 3D-Var, are conducted over a 4-month period during the winter of 2008 in the region of the Gulf of St. Lawrence. The analysis system assimilates RADARSAT image analyses produced by the Canadian Ice Service. The impact of additionally assimilating total ice concentration retrievals (using the NASA TEAM 2 algorithm) from passive microwave observations is also examined.
Celik, C.T. Performance of fuzzy logic in Determination of Unstable Points in Leveling Networks. Session 7   withdrawn file 
Abstract: A man-made structure or a natural phenomena like crustal movement, vertical land movements, subsidence determination, etc. is subject to deformations due to forces on them. Deformations of a structure require rigorous analysis to reflect significant movements. In general, deformation analysis relies on the reference network that is set outside the deforming area so that a comparison could be possible. However, a geodetic reference network is itself subject to movement. Therefore, the stability of the reference network needs to be conformed by some suitable methods. There exist a number of methods for performing the determination of stable and unstable points in the reference network. In this study, based on the fuzzy logic algorithm, a new method of detecting stable and unstable points in a reference network designed for detection of vertical movements is recommended. To test the proposed method a simulation was carried out and the results were presented.
Cessi, P., Wolfe, C.L. & Ludka, B.C. Eddy-balanced buoyancy gradients on boundaries and their role in the thermocline and the meridional overturning circulation Session 3 Tuesday 8  oral file 
Abstract: A model of the thermocline linearized around a specified stratification and the barotropic linear wind-driven Stommel solution is constructed. The forcings are both mechanical (the surface wind stress), and thermodynamical (the surface buoyancy boundary condition). The effects of weak diapycnal diffusivity and of eddy fluxes of buoyancy (parametrized as isopycnal diffusion) are included. These effects are especially important near the boundaries where they mediate the transport in and out of the narrow ageostrophic down/upwelling layers. The dynamics of these narrow layers can be replaced by effective boundary conditions on the geostrophically balanced flow. The effective boundary conditions allow buoyancy gradients along all solid boundaries, including the eastern one.
A special focus is on the buoyancies along the eastern and western walls, since their difference determines the meridional overturning circulation.
The advection of buoyancy by the barotropic flow is most effective on the western boundary. Because the buoyancy is coupled all along the boundaries, the barotropic flow advection affects the buoyancy everywhere, including the eastern wall. This leads to a complex interweaving of buoyancy patterns, with stacked counter rotating cells of the meridional overturning circulation. Quantitative scaling arguments are given for each of these cells, which show how the buoyancy forcing, the wind-stress and the diapycnal diffusivity, as well as the other imposed parameters, affect the strength of the overturn.
Chulli, B.Z. & Gabtni, H. Using gravity and geothermal gradients data to determine the crustal configuration and associated oil and gas accumulations in the Sahel Basin (Eastern Tunisia) Session 5 Thursday 10  poster file 
Abstract: An analysis of Bouguer gravity anomaly data and geothermal gradient data determined from bottom hole and drill stem tests temperature are used to determine the crustal configuration of the Sahel Basin in eastern Tunisia and its role in the maturation and location of the large number of oil and gas fields in the region. The regional Bouguer gravity anomaly field is dominated by gradual increase in values from the northwest to southeast and is thought to be caused by crustal thinning that is revealed by regional seismic studies. In addition, higher geothermal gradients in the same region as the Bouguer gravity anomaly maximum adds an additional constraint for the existence of crustal thinning in the region. A detailed analysis of the Bouguer gravity anomaly data was performed by both upward continuation and horizontal gradient analysis. These two techniques were combined to show that the region consists of two structure regions: 1) the NOSA-Zeramediner region in the NW part of the study area which is characterized by northwest-dipping, northeast-striking faults, thicker crust (30-31 km) and low geothermal gradients, and 2) the Mahres-Kerkennah region in the SE part of the study area which is characterized by vertical, northwest-striking faults, thinner crust (28-29 km) and a higher geothermal gradients. The correlation of mapped and geophysically-defined faults, volcanic rocks, thinned crust and high geothermal gradients with the location of known oil and gas fields indicates that the faults are a major determining factor in the location of these petroleum accumulations.
Cooke, R.M. The Risk Management Perspective in Climate Change Special Session Thursday 10  oral file 
Abstract: Economic models for costing climate change and determining the ``social cost of carbon'' are based on expected values; without proper uncertainty quantification. By comparison, Banks, insurance companies and other risk takers are, or are supposed to, operate on a different paradigm, namely capital sequestration based on value at risk. Such companies typically hold cash reserves to guarantee solvency in a 1-in-500 year loss event. The upper 5% value of climate sensitivity (the degree Celsius rise in global surface temperatures resulting from doubling atmospheric CO2) of the 22 peer reviewed studies cited in IPCC-AR4 (2007) is 7 degrees C, and some authors extract a 1% chance for values greater than 10 degrees C. At such levels, there are no arguments that society as we know it could persist. Adopting a risk management perspective would mean, (i) defensibly quantifying the uncertainty regarding future global temperatures out to a policy horizon at least equal to that currently employed by responsible risk takers (500 yrs), (ii) determining the costs of rapidly reversing the effects of CO2 emissions, should the actual global temperature values fall in the upper tail of its distribution, (iii) formulating a plan to insure the availability of such reserves in the event they are needed. The costs of such capital requirements should - but currently are not - reflected in discussions of the social cost of carbon.
Cordoba, G., Pitman, B. & Sheridan, M.F. A two-phase, depth-averaged model for geophysical mass flows in the TITAN code framework Session 1 Friday 11  oral file 
Abstract: Particle-laden surface flows are extremely dangerous natural phenomena, especially near urban and suburban areas. Tools for portraying their evolution in time, flow dynamics and run out distance are helpful in forecasting the consequences of potentially hazardous events, development of hazard and risk maps, and design of management policies in areas of potential natural disasters. One example, the TITAN2D code (Patra et al. 2005) that was developed to simulate granular flows in volcanic landscapes, is widely used in hazard assessment (see Murcia et al., 2010). Many of these models and the companion computer codes strictly apply only to dry granular flows like rock avalanches and volcanic block-and-ash flows. Because the presence of an interstitial fluid strongly modifies the dynamics of a flow, we developed a new, more general, two-phase model that is valid for flows that contain a broad range of water and solid mixtures, such as debris flows and mud flows. This model is based on balance laws for mass and momentum for each phase. The granular material is assumed to obey a Coulomb constitutive relation, and the fluid is assumed to be inviscid. The Darcy-Weisbach formulation is used to account for bed friction, and a phenomenological drag coefficient mediates the momentum exchange between phases. These equations are depth averaged, resulting in a system of 6 partial differential equations. The resulting equations correspond to the Savage and Hutter model in the limit of no fluid, and to the typical shallow water solutions (see Ortiz, 2005) for pure water. We implement the two-phase model within the familiar TITAN2D framework. The model has been tested on artificial topographic channels as well as on some volcano landscapes.
Murcia, H.F., Sheridan, M.F., Macías, J.L. and Cortés, G.P., 2010, TITAN2D simulations of pyroclastic flows at Cerro Machín Volcano, Colombia: Hazard implications, Journal of South American Earth Sciences, 29(2): 161-170.
Patra, A.K., Bauer, A.C., Nichita, C., Pitman, E. B., Sheridan, M.F, and Bursik, M.I, Rupp, B., Weber, A., Stinton, A.J., Namikawa, L.M., and Renschler, C.S., 2005, Parallel adaptive simulation of dry avalanches over natural terrain. Journal of Volcanology and Geothermal Research, 139(1-2):1-21.
Ortiz, P., Shallow water problems. in Zienkiewicz, O.C., Taylor, R.L., and Nithiarasu, P., The finite element method for fluid dynamics, 2005, Elsevier, 6th Edition, 400 pp.
Costa, A., Macedonio, G., Folch, A. & Durant, A. A model for ash aggregation in volcanic plumes Session 1 Thursday 10  poster file 
Abstract: We present a model to describe ash aggregates in a volcanic plume based on a solution of the classical Smoluchowski (1917) equation. The collision frequency function accounts for three different mechanisms: Brownian motion, ambient fluid shear and differential sedimentation. Since dry aggregation plays a minor role in presence of water, we consider sticking efficiency in wet conditions only, but the different binding effect of liquid water and ice is discerned. The model accounts for the fractal geometry of the aggregates as confirmed by laboratory experiments and observations. The proposed approach represents a first compromise between the full description of the aggregation process and the need to decrease the computational time necessary for solving the full Smoluchowski equation. We show results of a parametric study on the main model parameters and we estimate coagulation kernels and timescales of the aggregation process under simplified conditions of interest in volcanology. Moreover, the aggregation model is implemented in the FALL3D volcanic ash transport model and applied to the 18 May 1980 Mount St. Helens and the 17-18 September 1992 Crater Peak eruptions. The model provides a higher degree of agreement than previous fully-empirical aggregation models and successfully reproduces the depositional characteristics of the deposits investigated over a large range of scales, including the position and thickness of the secondary maxima. Further studies are carried out to investigate dry aggregation due to electrostatic forces.
Costa, A., Melnik, O., Sparks, S., Macedonio, G. & Gottsmann, J. Nonlinear dynamics of magma flows in conduits Session 1 Friday 11  oral file 
Abstract: Dynamics of magma flows in volcanic conduits are strongly nonlinear because of some important processes, often neglected, such as viscosity changes due to variations of temperature, crystal fraction and water content or coupling with wall-rock deformation. Crystallization of microlites and gas loss through permeable flows are time-dependent processes that can introduce strong feedback mechanisms which greatly amplify the effect on extrusion rates of small changes of chamber pressure, conduit dimensions or magma viscosity. When timescales for magma ascent are comparable to timescales for crystallization, there can be multiple steady solutions. Such nonlinear dynamics can cause large changes in dome extrusion rate and pulsatory behaviour of dome growth. Time scales of this kind of ciclicity are also affected by wall-rock elasticity. Two distinct periods can be associated to the elastic deformation of both magma chamber and dyke walls. Thermal coupling play also a crucial role in the dynamics of magma rise because can lead to significant viscosity stratifications due to temperature gradients produced across and along the conduit. In order to capture these effects mass, momentum and energy equations for magma flow inside a cylindrical conduit need top be coupled with the heat conduction equation in the surrounding host rocks with far-field conditions for rocks temperature. Simulation results show that both viscous dissipation and heat loss to the conduit walls have a pivotal role on magma dynamics allowing to distinguish three regimes: i) a conductive-heat-loss-dominated regime, ii) an intermediate regime, and iii) a viscous-heating-dominated regime. Finally, explosive volcanic eruptions have hitherto been largely interpreted in terms of multiphase flows through rigid conduits of a fixed cross-section, ranging from cylinders to parallel-sided conduits to simulate dykes. More realistic situations in which the explosive flows occur along dykes, taking into account wall rock elasticity showed major differences to results for rigid conduits, implying new processes and phenomena that affect the evolution of explosive eruptions. Dyke flows show very large deformation due to under-pressures in respect to lithostatic at the vicinity if fragmentation level where conduit significantly narrows or even closes. The role of an extensional far-field stress is also discussed.
Creyts, T.T. & Schoof, C.G. Drainage of water through subglacial water sheets Session 6 Tuesday 8  oral file 
Abstract: Dynamic lakes that cause ice surface elevation change require ample generation and supply of water from upstream sources. The hydraulic system delivers water to these lakes either through a distributed or channelized morphology. In this paper, we investigate the effects of subglacial water drainage resulting from spatially distributed water sheets. In our model, the weight of overlying ice is supported by both water pressure and various sizes of bed protrusions that penetrate the water sheet. Each of the various sizes bears a different magnitude of the overlying ice based on a linear stress recursion that balances forces at the bed. Previous results have shown that water depth can be a multi-valued function of both effective pressure (ice overburden minus water pressure) that drives sheet closure and hydraulic gradient that drives water flow (Creyts and Schoof, 2009). Curvature and structure of this multi-valued water depth function depend on the protrusion size distribution. Switches between different branches of the water depth relationship correspond to either the establishment or shut-down of a 'connected' (or efficient) drainage system. We build upon and extend previous work to show how along-path discharge affects water depth and switches from one state to another. We conclude by relating state behavior to subglacial conditions where dynamic lakes are found.
Cuffaro, M., Miglio, E. & Ruffo, P. Asymmetric Tectonics at Rift Zones: Geodynamics and Numerical Modeling Session 5 Friday 11  oral file 
Abstract: Tectonic evolution at rift zones is commonly considered symmetric along mid-ocean ridges, when modeling with relative plate motions and steady-state processes. However, the bathymetry of rift zones is generally asymmetric, being the eastern flank in average slightly shallower (100-300 m) than the western one. Also, based on surface wave tomographic models, shear wave velocities in the upper mantle indicate a difference between the western and eastern flanks of an oceanic basin. A better way to understand dynamics of the lithosphere at rift zones, and lithosphere/mantle interactions corresponds to absolute plate kinematic analyses, i.e., with respect to the mantle, modeling time-dependent tectonic processes. We performed numerical simulations of plate-driven mantle flow beneath mid-ocean ridges and we considered a time-dipendent flow induced by absolute motion of overlying rigid plates in an incompressible viscous mantle. In absolute frameworks, a net ``westward'' rotation of the lithosphere relative to the mantle can be observed, and we used velocities obtained in the hotspot reference frame, as boundary conditions. This implies that plates along a ridge, and the ridge itself, move toward the west but with different velocities, relative to the fixed mantle, and the separation between plates triggers mantle upwelling. Numerical solutions for viscosity flow beneath plates that thicken with increasing age are presented. The mantle can be modeled as a viscous fluid, and its dynamics can be described using the Stokes equations. At a first approximation the fluid is considered Newtonian. A further step in the description of the phenomena would require the inclusion of thermal effects: in this case the fluid viscosity and density have to be considered as a function of the temperature. For solving both the Stokes equations and the thermal effects, a finite element approach has been adopted. Results show an asymmetric thickening of plates along the ridge, as suggested by geological and geophysical observations, and provide useful relationships between mantle temperature and thickness of the oceanic lithosphere.
Currenti, G.M., Bonaccorso, A., Di Stefano, A., Scandura, D., Williams, C.A. & Del Negro, C. Stress changes variations around a pressured magma chamber: constraints on magma intrusion and fault slip Session 4 Monday 7  oral file 
Abstract: Pressure variations in a magma chamber may cause ground deformation and redistribution of the stress in the surrounding rocks. At Etna volcano eruptive sequences and local seismicity have revealed a significant time correspondence between volcanic unrest and rupture along fault structures. Ground deformation and stress change are often interpreted in terms of overpressure in a magma chamber without taking into account the total pressure on its walls. When studying magma reservoir failure and conditions for nearby fault ruptures, an analysis of stability is required and total pressure, rather than overpressure, has to be estimated. We performed a stress-strain analysis at Etna volcano using the Finite Element Method to investigate the perturbations to the local stress regime produced by pressure sources embedded in a inelastic domain. The main focus is to understand magma reservoir failure and characterize the conditions under which rupture on the pre-existing faults can occur. An initial state of stress is computed under the load of internal gravity body force. A gravitationally loaded model, that fails by slip on pre-existing fractures obeying a Coulomb friction law or as intact rock according to Drucker-Prager failure criterion, is considered. The investigation of elasto-plastic rheology proves useful, since elastic models are not able to explain the observed deformation without requiring unrealistically high overpressures beyond the crustal strength. Several simulations using different pressure sources, medium parameters and fault geometries were performed to demonstrate how different initial stress conditions can affect the magma chamber stability and slip pattern on the pre-existing fault structures. The magma pressure that the chamber wall can sustain without failing depends mainly on the initial stress state, which is affected by volcano topography only at shallow depth. As the magmatic pressure builds up, the simulations indicate that host rocks fail first at the crest of the magma chamber promoting earthquakes and upward magma propagation. Stress changes and slip distributions on the nearby pre-existing faults are also investigated to clarify the relationship between magmatic inflation and fault slip promoted by the induced stress changes. The proposed models place an upper limit on the pressure that a magma chamber can sustain before failing and provide a quantitative estimate of the interaction between volcano activity and tectonic processes.
D'Alpaos, A., Da Lio, C. & Marani, M. Equilibrium states and transient dynamics in tidal bio-geomorphic systems subject to climate change Session 6 Tuesday 8  oral file 
Abstract: We describe and apply a point model of the joint evolution of tidal landforms and biota to explore the equilibrium states and the transient behaviour of such systems under varying rates of sea-level change and of sediment supply. The model incorporates the dynamics of intertidal vegetation, benthic microbial assemblages, erosional, depositional, and sediment exchange processes, wind-wave dynamics. Alternative stable states and punctuated equilibria emerge, characterized by possible sudden transitions of the system state, governed by vegetation type, disturbances of the benthic biofilm, sediment availability and marine transgressions or regressions. Multiple stable states are suggested to result from the interplay of erosion, deposition and biostabilization, providing a simple explanation for the ubiquitous presence of the typical landforms observed in tidal environments worldwide. The explicit and dynamically-coupled description of biotic and abiotic processes thus emerges as a key requirement for realistic and predictive models of the evolution of a tidal system as a whole. The analysis of such coupled processes indicates that hysteretic switches between stable states arise because of differences in the threshold values of relative sea level rise inducing transitions from vegetated to unvegetated equilibria and viceversa, with implications for the preservation of tidal environments under a climate change. Finally, we explore the transient behaviour of the system forced by synthetic and observed sea-level rise forcings and identify the effects of the characteristic response time of vegetation to environmental changes on the overall system dynamics.
De Monte, S., Provenzale, A. & Sciascia, R. Plankton settling and survival Session 3 Tuesday 8  poster file 
Abstract: The balance between size-related metabolic activities and sinking is known to affect size composition of planktonic communities both on ecological and evolutionary timescales. One important driver toward the evolutionary reduction of plankton size is the fact that larger cells with negative buoyancy sink below the mixed layer. In this work, we address the role of turbulence and of its temporal fluctuations in the adaptive evolution of cell size. As pointed out by Magalef (1978), the adaptive success of planktonic cells must depend on their capability of matching mixing intensity in the water column. Turbulence in the mixed layer affects nutrient availability and light exposure, as well as the settling speed of negatively buoyant particles, therefore shaping the repartition of cells along the vertical gradients of nutrients and light. We incorporate the dependence of sinking rate on size in a model for adaptive evolution of cell size and show how different turbulent regimes (associated with different settling velocities) can have opposite effects on the adaptive value of organisms of a given size. We show that the outcome of phytoplankton competition on ecological timescales and of evolutionary dynamics on longer timescales depends on both the turbulent state of the mixed layer, affecting the dependence of the sinking rate on the cell size, and on the characteristic frequency of changes in turbulence properties. In general, larger-sized plankton is favoured for higher turbulence and for fast changes in the turbulent state. This opens a question on whether one can recognise regional-level processes, whereby the properties of small- scale turbulence (modulated by wind, for instance) cause an inversion of the evolutionary trend towards smaller sizes, and on the time scale of evolutionarily relevant variations in turbulence.
De Santis, A., Cianchini, G., Qamili, E. & Frepoli, A. The 2009 L'Aquila (Central Italy) seismic sequence as a chaotic process and implications for main shock predictability Session 4 Tuesday 8  poster file 
Abstract: In this presentation we shall demonstrate that the seismic sequence of foreshocks culminating with the recent Mw =6.3 main shock on April 6, 2009 in L'Aquila (Central Italy) evolved as a chaotic process. To do this, we apply a nonlinear retrospective prediction to L'Aquila seismic sequence and look at the behaviour in time of the error between predicted time and true occurrence of the main shock when gradually increasing parts of the sequence are considered. This is a typical nonlinear approach which is quite powerful to detect chaos in relatively short time series. The method of prediction is based on the Accelerated Strain Release (ASR) analysis that considers magnitudes and occurrence times of foreshocks. We find that i) the prediction error decays exponentially with a time constant of about 10 days and ii) at around 6 days before the main shock, ASR is quite powerful for anticipating the time of occurrence with an uncertainty of about a day. Due to their retrospective characteristics, the latter results could be affected by changes on some a-priori parameters used in the ASR technique. However, since the error behaviour in time should not depend, in principle, on the type of prediction technique, we consider these results to be strong evidence that the studied sequence of foreshocks is a chaotic process with K-entropy of about 0.1/day.
Devauchelle, O., Petroff, A.P. & Rothman, D.H. Coupling groundwater to sediment transport in ravines Session 2 Tuesday 8  poster file 
Abstract: The formation of a ravine network by seepage of an aquifer can produce a highly complex pattern, even when sediment properties and precipitation are uniform. The steephead streams of the Apalachicola Bluffs and Ravines Preserve on the Florida Panhandle are a simple example of this type of landscape formed in a homogeneous sand plateau, in the apparent absence of runoff. We investigate the mechanism of their formation by means of a two-dimensional model for the aquifer, based on both Darcy's law and the shallow-water approximation: the square of the water table elevation satisfies Poisson's equation. Once a tip is formed, water flowing through it transports sediment removed from the tip, allowing the ravine to grow forward. We assume that the growth of tips is quasi-static, that is, the value of the flow-induced shear stress is close to the threshold for grain motion on the entire stream bed. The consequences of this assumption on the streams are the following: 1) Their aspect ratio and average velocity are constant throughout the network; 2) Their width is proportional to the square root of their discharge; 3) The product of their discharge to the square of their slope is constant. Field measurements show reasonable accordance with the two first relations. As for the third one, it is a boundary condition for the Poisson equation representing the aquifer, to be satisfied along the streams. One consequence of this coupling between ground water and bedload transport is that, in the tips neighborhood, stream elevation profiles increase with the distance from the tip to the power two thirds. This prediction is consistent with observation. This result raises new questions about the structure of the network and the processes involved in its formation. In particular, preliminary numerical calculations tend to indicate that the necessity for a stream to present an upward concave profile generates new solubility conditions for the system. Equivalently, the equilibrium condition would impose a strong constraint on the length of the channels, thus impacting the network properties.
Di Martino, R.M.R., Camarda, M., Gurrieri, S. & Valenza, M. Modelling the chemical composition transients in volcanic soil gases to evaluate the source depth of the gas reservoir Session 2 Tuesday 8  poster file 
Abstract: Continuous and high frequency monitoring of the chemical composition changes of the gas emissions aimed to the volcanic activity surveillance, requires the investigation of the transport dynamics in a multiphase system. In order to study the relationships among the chemical composition of the gas emissions, the sub-surface gas flux regime, the features of the porous medium and the depth of the gas reservoir, we have designed several laboratory experiments and developed the theoretical model needed to understand both the results as well as the data collected in the field measurements. In our model, we describe the gas transport through a porous medium as the result of both diffusive and advective processes, stressing the mass transfer features during a transient state induced through the disturbance of the flux variables of the system. A cylindrical container filled with volcanic ashes of known size and shape of the grains has been used in laboratory as gas flux simulator. At the uppermost base, it was provided with both the H2 and CO2 detectors. The chemical gas disturbance was induced at the lower base of the cylinder through a gas flux ranging from $10^-3$ to $10^-1$ $cm^3s^-1cm^-2$ of H2 and CO2 in N2 matrix. This interval well includes the experimental values measured in the field. The temperature and pressure conditions were also controlled. Our results show that the first arrival times of the chemical gas disturbance depend on the height of the cylinder, the advective flux velocity and the diffusive features of each gas species. We also observe that the CO2 concentration disturbance regularly arrives later than the H2 one and the time delay between them increases as the advective flux velocity decreases. We believe that the different diffusive features of the H2 and CO2 molecules can explain such a sort of chromatographic effect. Although the advective flux is the more efficient transport process, the diffusive one operates towards the separation of the mixture in its own components. Since the arrival time of the concentration disturbance of each gas species depends also from the distance of the source, the measurements of both their time delays and the advective gas flux velocity allow the evaluation of the depth of the gas reservoir. The high frequency monitoring data acquired at Stromboli volcano have shown good agreement with laboratory results since the main changes of H2 concentration precede the CO2 ones systematically up to a few tens of hours.
van Dinther, Y., Gerya, T.V., Dalguer, L.A., Mai, M. & Morra, G. The long-term seismic cycle within geodynamic numerical simulations of a subduction zone Session 4 Tuesday 8  poster file 
Abstract: Earthquakes occur over different spatial and temporal scales. Due to a lack of direct observations, a full understanding of the long-term evolution of seismicity within subduction zones remains elusive. Without this information, it is difficult to quantify the recurrence interval of large earthquakes based on observations over our current, limited time period. Realistic numerical modeling of subduction zone physics can help to improve our understanding of the long-term seismic cycle. In this study, we aim at understanding the physics driving this cycle. By quantitatively balancing energy sources and sinks we try to estimate the types of mechanisms that are potentially driving seismic energy release. At the same time, we try to determine the causes of and space-time relationships for the activity within different clusters, where special attention is given to the potential presence of compressional events within the bending area.
We use a plane-strain finite-difference scheme with marker-in-cell technique to solve the conservation of momentum, mass, and energy for a visco-elasto-plastic rheology (code I2ELVIS). In our 4000$200 km$^2$ numerical model of ocean-continent convergence, a generic oceanic plate subducts below a continental overriding plate. In a petrologically realistic setting, localizations of plastic strain are spontaneously formed when the second invariant of the deviatoric stress tensor exceeds the Drucker-Prager yield stress. This leads to a correction of stresses to the yield stress at each marker, accompanied by a decrease in overall viscosity to mimic local weakening.
Preliminary results show the existence of spontaneous clusters of plastic strain localizations, within the thrust, outer-rise bending, and back-arc areas. Within these localizations, we observe periods of rapid deformation that are accompanied by an immediate stress drop. The activity of the clusters is closely linked to the most energetic events at the thrust interface that have a periodicity of 25-30.000 years. In periods of reduced coupling, the thrust slips more rapidly, and promotes extension within both the bending and back-arc areas. This coupling suggests that pressure unloading drives faulting within these areas. Tentatively, the main thrust event is preceded by activity within the bending area. Increasing time resolution and comparison with observational evidence, should tell us whether this could provide a numerical tool to help forecast major earthquakes.
Dioguardi, F., Dellino, P. & De Lorenzo, S. 1-D numerical simulation of the gas-particles flow in vertical conduits during large-scale experiments on the mechanics of explosive eruptions Session 1 Friday 11  oral file 
Abstract: Explosive volcanic eruptions are the most dangerous events that can affect the areas surrounding active volcanoes. It is well known that conduit exit conditions play a major role in determining the rate and style of explosive eruptions, which in turn determines the typology of impact on the areas subjected to volcanic risk. For these reasons many researchers tried to develop models for the forecasting of conduit exit conditions. Theoretical and numerical modeling are the most widespread approach in volcanology. Experimental models are not common because of the difficulties to set-up large-scale experiments. For this reason theoretical and numerical models always lack in the experimental validation. Dellino et al. (2007) developed the first large-scale experiments on the mechanics of explosive eruptions that proved to scale well to natural cases. From these experiments, Dellino et al. (2009) developed an experimental model for the forecasting of exit velocity of eruptive mixtures and the conditions of stability of the main eruptive styles. The next step is to create numerical models that reproduce the large-scale experiments and then to apply them to the real eruptions. This would be the first time that a numerical model is validated against large-scale experiments in volcanology. Here we present a steady 1-D two-phase numerical model of the gas-particles flow in the conduit, that is the first stage of the experiments. The equations of conservation of mass and momentum for gas and volcanic particles are solved via the Runge-Kutta scheme with an adaptive stepsize. All the experimental runs were simulated. The model takes in account the real shape of volcanic particles and uses the well established law of Dellino et al. (2005) for the calculation of particles terminal velocity. The pressure gradient in the conduit, which is the main driving force of the vertical two-phase flow, is obtained by the same empirical model. Finally the interphase drag force and the friction between the phases and the conduit wall are included by using the empirical laws for wall-particles and wall-fluid frictions developed in industrial engineering and classic fluid dynamics. All the experimental runs have been simulated and the model results in mixture exit velocities that are consistent with the ones measured in the experiments. Therefore the numerical model is promising and can be improved and used for the simulation of real cases.
Doronzo, D.M. & Dellino, P. Fluid dynamics of volcaniclastic turbidity currents Session 3 Tuesday 8  poster file 
Abstract: Subaqueous sediment density currents have variable flow structures, depending on particle size and volumetric concentration. They form both in ocean basins and in volcanic districts where subaqueous density currents are generated by primary pyroclastic material entering the sea. A particular type of such flows is represented by volcaniclastic turbidity currents, which form by the subaqueous remobilization of primary pyroclastic deposits. In this work, a fluid dynamic model of volcaniclastic turbidity currents is presented and applied to the currents that generated the Craco (Southern Italy) volcaniclastic deposits. The deposits consist of thick laminated to massive beds, which formed from the remobilization of thin fine-ash pyroclastic fall deposits, in an area far away from the original volcanic source. Our model is based on the sedimentological features of deposits and on particle physical characteristics. It allows the calculation of velocity, thickness and Reynolds number of the current. Furthermore, a relationship between the current thickness, laminae thickness and paleo-slope angle is obtained, which allows to evaluate the thickness of the primary pyroclastic fallout deposit. Our model can be useful for obtaining an estimation of the primary distal ash dispersal by the carachteristics of the associated subaqueous volcaniclastic deposit. Details can be found in Doronzo and Dellino (2010), in press (JVGR, doi:10.1016/j.jvolgeores.2010.01.017).
Doronzo, D.M., Dellino, P., de Tullio, M.D. & Pascazio, G. Immersed boundary simulation of pyroclastic density currents: numerical scheme and experimental validation Session 1 Thursday 10  poster file 
Abstract: Pyroclastic density currents are ground hugging, hot, gas-particle flows representing the most hazardous events of explosive volcanism. Their impact on structures is a function of dynamic pressure, which expresses the lateral load that such currents exert over buildings. In this paper we show how analog experiments can be matched with numerical simulations for capturing the essential physics of the multiphase flow. We used an immersed boundary scheme for the mesh generation, which helped in reconstructing the steep velocity and particle concentration gradients near the ground surface. Results show that the calculated values of dynamic pressure agree reasonably with the experimental measurements. These outcomes encourage future application of our method for the assessment of the impact of pyroclastic density currents at the natural scale.
Dubinkina, S., Goosse, H., Crespin, E. & Sallaz-Damaz, Y. Data assimilation using a particle filter in twin and real data experiments with a coupled climate model. Session 8 Thursday 10  poster file 
Abstract: Particle filtering is a relatively new and promising approach for data assimilation in climate models. It is a challenging task to develop a non-degenerative particle filter for a large-scale geophysical model using a relatively small number of particles. Here we implement a particle filter for data-assimilation method in a climate model, focusing on interannual to centennial timescales. Several tests are performed where we investigated various approaches to compute the cost function. These approaches are based on assigning different weights to particles and calculating the likelihood of each ensemble member from different hypotheses. In those tests, we use pseudo-observations obtained from twin experiment using the coupled climate model LOVECLIM, as well as the real data observations obtained over the last century.
Dufek, J. & Manga, M. Flow transformation in explosive volcanic eruptions: Multiphase and multiscale interactions in pyroclastic density currents Session 1 Friday 11  oral file 
Abstract: Explosive volcanic eruptions produce turbulent, multiphase flows that encompass a vast range of scales from micron-scale ash to eruptive plumes that can extend 100s of kilometers. Abrupt flow transformation (e.g. from dense to dilute pyroclastic density currents) can arise due to the energy exchange across multiple length scales and phases. Understanding these emergent phenomena is predicated on describing the phase-relative motion of many particle sizes and the multiphase transport physics of the transition from dilute to dense flows. I will discuss the use of multiphase models in addressing these different scales of fluid motion as well as how they can provide a platform to integrate microphysical, analogue experiments and observational constraints. Microphysical experiments can provide the necessary closure for statistical mechanics based models, and provide a way to examine grain-scale processes in a probabilistic manner. Such small-scale processes can dramatically alter the flow dynamics. For example, steam explosions associated with pyroclastic density currents entering the sea at Montserrat can be produced if energy exchange at the scale of the smallest particles is considered, while only considering large scale heat transfer results in passive steam production. Propagation of pyroclastic density currents and plume dispersal is also strongly controlled by the grain size distribution of particles, which can evolve as a result of both comminution and aggregation processes occurring at the grain scale. Scaled analogue experiments can provide insight into emergent flow dynamics such as secondary plumes in pyroclastic density currents, and can be compared in detail with direct numerical simulations of these turbulent phenomena. Multiphase modeling also provides a means of linking granulometry in deposits to grain size sorting by dynamic structures and flow transformations from dense to dilute flows. In particular I will focus on flow transformation related to particle-fluid thermal energy exchange, particle-boundary energy exchange and resuspension, and abrupt flow transformation due to comminution, breaks in slope and topographic barriers. Each process will be discussed in the context of a natural example using the Kos Plateau Tuff eruption, eruption of Mount St. Helens in 1980, the dome-collapse eruptions at Montserrat, and the boiling-over eruption of Tungurahua in 2006 as end-members.
Eicker, A. The determination of mass transport processes in the Earth's system from satellite gravity missions Session 5 Friday 11  oral file 
Abstract: The analysis of the Earth's time variable gravity field plays an important role in the research of the Earth's system. The satellite gravity mission GRACE (Gravity Recovery And Climate Experiment) provides, for the first time, a direct measurement of the amount of mass that is redistributed at or near the surface of the Earth. These redistribution processes include, for example, the oceanic and atmospheric circulation, water fluxes between various terrestrial water storages, melting ice, river discharge, changing sea level and convective flows in the Earth's mantle. The satellite observations provide global data sets of unprecedented homogeneity and resolution and have significantly improved our understanding of the Earth's gravity field and the processes that lead to its temporal changes.
This talk will discuss the observables of satellite geodesy and their relation to potential change and mass change. Furthermore, a variety of examples from different geophysical applications will be shown to demonstrate how the interpretation of satellite gravity data can improve our knowledge of the mass transport phenomena.
Esposti Ongaro, T., Neri, A., Clarke, A.B., Smekens, J., Voight, B. & Widiwijayanty, C. Multiphase flow dynamics and damage zonation of the May 18, 1980 lateral blast of Mount St. Helens: comparison between 3D numerical simulations and field observations Session 8 Friday 11  oral file 
Abstract: Volcanic directed blasts are among the most hazardous and fascinating of volcanic phenomena. They are characterized by an extraordinarily large rate of mass and energy discharge, extremely fast dynamics and a remarkably broad area of severe damage, despite the relatively low amount of magma involved. The rapid decompression and expansion of the eruptive mixture and the subsequent propagation of the pyroclastic density current are characterized by complex thermo-fluid dynamics, involving the propagation of a multidispersed mixture of gas and particles and interaction with the rugged, 3D topography. In this study, we adopt a 3D, transient model (Neri et al., J. Geophys. Res., 2003; Esposti Ongaro et al., Parallel Computing, 2007) based on the Eulerian multiphase transport equations to reconstruct ``a posteriori'' the dynamics of the well-documented May 18, 1980 blast at Mount St. Helens (USA). The detailed comparison of 3D model results with digitized, geo-referenced observational data demonstrates the ability of the model to reproduce the MSH blast front advancement and basic deposit characteristics and distribution, with initial conditions well constrained by available geologic data (including gas content, mass of juvenile and entrained rocks, temperature, grain size distribution, and pre-eruption pressure distribution in the cryptodome). The simulations suggest that much of the severe damage observed at MSH can be explained by high dynamic pressures in stratified gravity currents, and that the rapid decrease in damage and dynamic pressure from proximal to distal areas may be related to rugged topography beyond the North Fork Toutle River valley. Despite the considerable number of uncertain parameters in the model and in the observations, our present results demonstrate that 3D transient and multiphase flow models can accurately reproduce the main large-scale features of blast scenarios and can help in the definition and assessment of flow hazards at blast-dangerous volcanoes worldwide.
Evonuk, M. Differential flow in planets prior to core formation Session 5 Thursday 10  poster file 
Abstract: Differential rotation in the outer core plays an important role in the formation and maintenance of the Earth's magnetic field. The Earth's liquid outer core experiences a change of only 0.2 density scale heights with depth and is therefore usually modeled as a Boussinesq fluid. While this density change is small, 3D simulations do differ from constant density simulations (G. Glatzmaier, personal communication). Differential rotation in 3D Boussinesq simulations is maintained by vortex stretching of convective fluid columns due to the sloping impermeable boundaries. Meanwhile, 2D anelastic simulations of rotation in the equatorial plane show that differential rotation can form via local generation of vorticity as fluid parcels move radially, expanding or contracting with respect to the background density stratification. Constant density (Boussinesq) simulations in the same 2D geometry, of course, do not generate differential rotation. The question is then, in a 3D anelastic body, what is the relative importance of the density stratification to jet formation? Further, do multiple jet structures arise in anelastic bodies under conditions that would not form multiple jets in a Boussinesq body, i.e. in bodies without solid cores, prior to the formation of a solid inner core? This study looks at three-dimensional numerical simulations of thermal convection in inertial, rotating spherical fluid bodies with and without density stratification to explore the role of density stratification on generating and maintaining a zonal flow without a solid inner core.
Farrokhrouz, M. & Asef, M.R. Effects of Confinement on Rock Velocity Relations for different Rock Types Session 4 Tuesday 8  poster file 
Abstract: Rock velocity measurements are carried out on the cores obtained from the desired depth of target formation in laboratorial environment. In reality subsurface rocks are encountered with confining pressure of neighboring rocks. On the other hand, presented formulas to predict shear wave velocity from compressional wave velocity are usually based on laboratorial tests and under no confinement. These formulas are employed in continuous velocity calculations in oil and gas fields or medium depth constructional projects where rock samples restrained by nearby formations. Some discrepancies between real and predicted values may originate from such difference in circumstances. The effect of confinements on a series of published data gathered from different locations with various rock types was studied. The samples included main categorizations of the rocks: igneous, metamorphic and sedimentary rocks with different depths of coring. All the samples have been tested in the lab while definite confining pressures were applied to the cores. The common trends of these rocks showed that increasing the confining pressure could lead to increase of R-square between shear and compressional velocities and better estimation of Vs from Vp for clean sandstone. But the trend showed a cyclic behavior for metamorphic and igneous rocks like basalt, granite, serpentine and gneiss. In these samples, R-square value decreased and then increased while confining pressure was ascending. The trend was reverse for other sedimentary rocks, especially when impurities were present in the framework. Chert and shale samples in this study showed a trend of increase and then decrease in R-square values. In all, rock types showed different trends for R-square changes. For better estimation and conclusion, each data point velocity values was plotted vs. confining pressure and the same trend was observed. In addition, Vp/Vs ratio was also plotted as a function of confining pressure and the results conveyed similar conclusions. This may be suggested as a method for rock type definition when rock velocities are the only available data. It can be mentioned that shear wave velocity estimation only by compressional velocity could not be enough, because of changes in R-square value with various confinements. It seems that another parameter should be considered which makes prediction independent of regional pressures.
Flandoli, F., Giorgi, E., Aspinall, W.P. & Neri, A. Comparing the performance of different expert elicitation models using a crossvalidation technique Session 7 Thursday 10  oral file 
Abstract: The problem of ranking and weighting experts' performances when their judgments are being elicited for decision support is considered. Different methods exist to score experts such as the Cooke Classical model, Equal Weights, and Single Best Expert models. The purpose of our research is twofold: 1) to introduce a new model for providing optimal point-wise estimates, named the Expected Relative Frequency (ERF) model, and 2) to use a cross-validation technique to compare different scoring models' performances on real data. The new ERF model gives positive value to an expert who provides central values close to the known values. In contrast to a self-referencing conventional comparison of models, a cross-validation technique is adopted here. A single round of cross-validation involves partitioning a sample of seed items into two complementary subsets, computing expert weights from one subset (called the training set), and validating the result on the other subset (called the validation set or testing set). To control variability, average over multiple rounds of cross-validation is performed. Here, the cross-validation analysis is applied to five elicitation datasets from different expert judgment problems. For the comparison step of the cross-validation method it is necessary to introduce the concept of reward and the Cooke Classical model, the new ERF model, Equal Weights and the single Best Expert performances are compared with respect to three specific reward indices, called Calibration, Informativeness, and Expected Accuracy. These measure both an expert's ability to assess uncertainty and his accuracy in point-wise value estimation. A detailed comparison of average rewards and the probability of success of one method over another is given for the case studies. Results show that, for any specific reward index, there is only a limited probability that one method is better than another and, even in the average, no single method emerges as best for all the forms of reward. However, certain tendencies are observed: the Cooke Classical model is generally most suitable for assessing uncertainties, while the new ERF model results are preferable if the criterion is central value estimation accuracy, alone; the performance of the Equal Weights model is generally quite good but rarely, if ever, optimal, and Best Single Expert and individual single experts are ordinarily outperformed by any of the models based on pooled opinions.
Fournier, A. Combining geomagnetic data with physical models of Earth's core dynamics: An introduction to data assimilation in geomagnetism
Session 5 Friday 11  oral file 
Abstract: Data assimilation in geomagnetism designates the set of inverse methods for geomagnetic data analysis which rely on an underlying prognostic numerical model of core dynamics. Within that framework, the time-dependency of the magnetohydrodynamic state of the core needs no longer be parameterized: The model trajectory (and the secular variation it generates at the surface of the Earth) is controlled by the initial condition, and possibly some other static control parameters. The primary goal of geomagnetic data assimilation is then to combine in an optimal fashion the information contained in the database of geomagnetic observations and in the dynamical model, by adjusting the model trajectory in order to provide an adequate fit to the data.
The recent developments in that emerging field of research are motivated mostly by the increase in data quality and quantity during the last decade, owing to the ongoing era of magnetic observation of the Earth from space, and by the concurrent progress in the numerical description of core dynamics.
In this talk I will review briefly the current status of our knowledge of core dynamics, and elaborate on the reasons which motivate geomagnetic data assimilation studies, most notably a) the prospect to propagate the current quality of data backward in time to construct dynamically consistent historical core field and flow models, b) the possibility to improve the forecast of the secular variation, and c) on a more fundamental level, the will to identify unambiguously the physical mechanisms governing the secular variation and to make inferences on the interior of Earth's core, which is not directly sampled by geomagnetic observations. I will next present the different approaches to geomagnetic data assimilation which have been followed so far, and discuss in particular what those have already taught us regarding the state of the interior of Earth's core.
Gaina, C., Werner, S., Saltus, R., Medvedev, S. & group, C.-G. Structure and evolution of the Arctic in the light of new geophysical compilation and regional kinematics Session 5 Friday 11  oral file 
Abstract: New Circum-Arctic maps of magnetic and gravity anomaly have been produced by merging regional gridded data. Satellite magnetic and gravity data were used for quality control on the long wavelengths of the new compilations. The new Circum-Arctic digital compilations of magnetic, gravity and some of their derivatives have been analyzed together with other freely available regional and global data and models in order to provide a consistent view of the tectonically complex Arctic basins and surrounding continents. In particular, available tomographic models have been also scrutinised and evaluated for their potential to reveal the deeper structure of the Arctic. Tectonic boundaries (including continent-ocean boundaries and sutures) have been mapped mainly based on potential field data and their derivatives. In areas where the crustal age remains speculative we compare the crustal thickness derived from gravity inversion with other geophysical constrains. Based on our data analysis, we present a kinematic scenario as part of a larger tectonic framework, where subduction of the Pacific and South Anyui oceans led to the opening of the Amerasia Basin, motion between the North American plate relative to a fixed Eurasian/Lomonosov plate led to opening of small basins between the Lomonosov Ridge and Alpha Ridge area, and a mantle thermal anomaly (precursor of Iceland hotspot?) weakened the crust along the north/northeast part of the American craton facilitating the opening of the Canada Basin, rifting in the Baffin Bay and creating a zone of weakness in the area of future Eurekan orogeny.
Ganci, G., Vicari, A., Bonfiglio, S. & Del Negro, C. Hotsat volcano early warning system based on a combined use of SEVIRI and MODIS multispectral data. Session 9 Tuesday 8  poster file 
Abstract: Multispectral infrared observations carried out by the spacecrafts have shown that remote sensing of high-temperature volcanic features is feasible and robust enough to turn into volcano monitoring. Especially meteorological satellites have proven a powerful instrument to detect and monitor dynamic phenomena, such as volcanic processes, allowing very high temporal resolution despite of their low spatial resolution. An automated GIS-integrated system for thermal volcano monitoring, called HOTSAT, was developed at INGV-CT for analyzing both EOS-MODIS and MSG-SEVIRI satellite data. The alert system is composed by three packages: Pre-processing, Product Generation and Post-processing. Each package consists of several independent executable modules. The modules of Pre-processing are necessary for initial images geolocation and calibration, the modules of Product Generation compute higher-level products from the satellite band data as volcano hot spots and radiant flux estimation, and the ones of Post-processing project raw geometry data to a cartographic reference system and export the elaborated outputs to a Google Map platform. To locate volcano hot spots, a new contextual algorithm is introduced taking advantages from both spectral and spatial comparison methods. On a first step, the spatial standard deviation is computed on the difference between middle infrared (MIR) and thermal infrared (TIR) temperatures. These data are used to set an adaptive threshold and detect ``potential'' hot pixels. Those pixels are then further assessed as true or false hotspot detections base on statistical thresholds test derived from the MIR brightness temperatures. Following this procedure, all the computations are based on dynamic thresholds reducing the number of false alarms due to atmospheric conditions. The derivation of radiant flux is computed at all ``hot'' pixels. This is carried out using the MIR radiance method introduced for forest fires. Following this approach, the radiant flux is proportional to the calibrated radiance associated to the hot part of the pixel computed as the difference between the observed hotspot pixel radiance in the MIR channel and the background radiance that would have been observed at the same location in the absence of thermal anomalies. The HOTSAT early warning system is now suitable to be employed in an operational system of volcano monitoring. To validate and test the system some real cases on Mt Etna are presented.
Geese, A., Mandea, M. & Lesur, V. SAMS - the South African Magnetic model made of Splines Session 5 Thursday 10  poster file 
Abstract: Geomagnetic repeat surveys are conducted in many nations all over the globe in order to produce regional declination maps. To meet the scientific characteristics of magnetic fields, this requires to model potential fields on a regional scale. We revisit an approach introduced by Shure et al. (1982) relying on harmonic splines, but we limit ourselves to truncated series, i.e. we include only finite order L. This allows a direct accordance to existing global main field models. We apply the technique to a repeat station data set from southern Africa covering the years 1961-2001. In this region, the magnetic field is weak and changes rapidly. The resulting SAMS model gives a detailed insight in the exceptional behaviour of this area.
Gelman, S.E., Bergantz, G.W. & Bachmann, O. Quantifying the role of individual sources of error and uncertainty in model validation: A case study in conductive heat transfer Session 8 Thursday 10  poster file 
Abstract: The process of model validation requires a priori assessment of both error and uncertainty in model implementation (theory and numerical implementation) and constitutive and closure relations. Here we stress that validation studies in solid earth simulations have rarely included formal uncertainty analyses of transport and thermo-physical properties. Given that Earth materials are usually multicomponent, multiphase ensembles with wide ranges in uncertainty, it is difficult to assess what constitutes robust validation. In this presentation we will expand on the formal definitions of error and uncertainty for application to heat and mass transfer models as exemplified in the context of canonical heat transfer by conduction.
We define error as a deficiency not due to a lack of knowledge, typically arising for reasons of practicality and tractability. They can be acknowledged errors, such as round-off, discretization, or convergence errors, or unacknowledged, such as programming errors. Uncertainty represents a lack of knowledge, and may be classified as either aleatory or epistemic. An aleatory uncertainty is due to the natural variability in known quantities, that cannot be eliminated or entirely accounted for, and are often treated with a pdf approach. Aleatory uncertainty is not due to a lack of knowledge, but from the fundamental inability to constrain a parameter outside a window of variability. Epistemic refers to a reducible uncertainty that is due to a lack of knowledge about the system, and may be improved by additional observations or measurements. Both types of uncertainty commonly occur in solid earth applications.
Both error and uncertainty contribute to issues with validation and prediction, complicating the assessment of verisimilitude. In solid earth simulation, little attention has been paid to the interplay between sources and magnitude of error and uncertainty. In our presentation we will demonstrate how global solution uncertainty can be decomposed, or at least quantified, by identifying the greatest individual sources of inaccuracy. We exemplify this approach by considering conjugate conduction where three physical properties, density, conductivity and specific heat capacity contribute to heat transfer. We will revisit and assess the validity of the claims of Whittington et al. (2009) who measured precise thermal diffusivity for five rocks from 300 K to 1250 K, applying the average of these values to global continental crust.
Ghods, A. & Arkani-Hamed, J. Giant Impacts Control Mantle Dynamics and Cripple the Core Dynamos of Planets Session 5 Friday 11  oral file 
Abstract: A giant impact creates a huge basin on the surface and heats a planet by the shock waves travelling in the interior. We investigate the impacts that have created the largest basins existing on terrestrial planets: Caloris on Mercury, Utopia on Mars and Aitken on the Moon, all formed before 4 Ga. Due to plate tectonics on Earth and pervasive resurfacing on Venus at ~500 Ma no ancient giant basins are found on these planets. Shiva is probably the largest existing impact basin on Earth that was created at around 65 Ma, and we regard Arthemis Corona as a possible impact site on Venus. We determine the impact-induced temperature increase in the interiors of the planets using the scaling laws between an observed basin diameter and the impactor size [Holsapple, Annu. Rev. EPS, 1993], the shock wave pressure distribution model of Pierazzo et al. [Icarus, 1997] and the ``foundering'' shock heating model of Watters et al. [J. Geophys. Res., 2009]. The impact effects on the thermal evolution of a planet is investigated using 2D axi-symmetric and two-flow convection in a self-gravitating spherical shell of temperature- and pressure-dependent viscosity, temperature-dependent thermal conductivity, and pressure-dependent thermal expansion. The mantle is allowed to melt when its temperature surpasses the solidus. The impacts that created Shiva and Arthemis had minor effects on the internal dynamics of Earth and Venus. We determine the minimum size for an impactor required to influence the internal dynamics of these planets. For the small planets, the impact heating significantly enhances the temperature of the upper mantle and creates a superheated giant plume which ascends rapidly and develops a strong convection in the mantle of the sub-impact hemisphere, while the antipodal hemisphere remains almost undisturbed at least for the first 50 Myrs. The upwelling of the giant plume rapidly sweeps up the impact heated base of the mantle and replaces it with the cold surrounding material, thus hampering the effects of the impact-heated mantle on the core dynamo that might have existed prior to the impact. However, direct shock heating stratifies the core, suppresses the pre-existing thermal convection in the core, and cripples a preexisting thermally driven core dynamo. Depending on the size of the impactor and the planet, it takes 20-90 Myr for the stratified core to exhaust impact heat and resume global convection, possibly regenerating a core dynamo.
Gildor, H., Carlson, D.F., Fredj, E. & Rom-Kedar, V. Deducing an upper bound to the horizontal eddy diffusivity using a stochastic Lagrangian model Session 3 Tuesday 8  oral file 
Abstract: We present a method for estimating the upper bound of the horizontal eddy diffusivity using a non-stationary Lagrangian stochastic model. First, we identify a mixing barrier using a priori evidence (e.g., aerial photographs or satellite imagery) and using a Lagrangian diagnostic calculated from observed or modeled velocities (for instance, the relative dispersion or finite time Lyapunov exponent). Second, we add a stochastic component to the observed (or modeled) spatially non-trivial, time-dependent velocity field. The stochastic component represents sub-grid stochastic diffusion and its maximum possible magnitude is set by the eddy diffusivity. The relative dispersion of Lagrangian trajectories is computed for increasing values of the eddy diffusivity until the mixing barrier is no longer present. The value at which the mixing barrier disappears provides a dynamical estimate of the upper bound of the eddy diffusivity. The erosion of the mixing barrier is visually observed, and quantified by computing higher moments (for instance, the kurtosis) of the relative dispersion for each value of the eddy diffusivity. We demonstrate our method using a double gyre circulation model and apply it to high frequency radar observations of surface currents in the Gulf of Eilat.
Gillet, N., Schaeffer, N. & Jault, D. Rationale and geophysical evidence for quasi-geostrophic rapid dynamics within the Earth's outer core Session 5 Thursday 10  poster file 
Abstract: In this paper, we present arguments supporting the hypothesis that the flow in the Earth core for the time scales of the secular variation, is well described by a quasi-geostrophic (QG) model (that is almost invariant along the rotation axis).
A previous study showed that for axisymmetric motions, even at Elsasser numbers of order unity, the short time-scale flow is geostrophic inside the Earth core. Here, we extend this result to non-axisymmetric motions. The numerical simulations exhibit a columnar behaviour at parameters representative of the Earth core.
In addition, we present the results of several inversions of the core flow, showing that
a) for the same number of parameters, a symmetric (QG) flow model explains more of the secular variation than a tangentially geostrophic (TG) flow.
b) the flow seems to be more and more symmetric with time, possibly a consequence of the better quality of the measurements at magnetic observatories.
Goren, L., Sparks, D.W., Aharonov, E. & Toussaint, R. Modeling coupled fluid-grain deformation, with implications for landslides, fault-zones, and liquefaction Session 2 Monday 7  oral file 
Abstract: Coupled shear deformation of granular media and pore-fluids often gives rise to catastrophic failures above and below ground: landslides, soil liquifaction and earthquakes. However pysical understanding of this coupled system is surprisngly limited. We study coupled granular-deformation and fluid pore pressure changes within a shearing layer, using a new two-scale two-phase model: the grains are modeled at the grain scale using a granular dynamics algorithm, while the fluid is simulated on a slightly coarser Darcy porous-flow scale. The formulation for the fluid is based on first principles and is general enough to account for the fluid response to both elastic and inelastic finite deformation of the granular matrix. We find that the presence of fluid influences all aspects of deformation, and introduces important physical constraints on upscaling lab results to the field scale. In simplified models that assume an infinitely stiff matrix, pore pressure responds to granular deformation, but granular deformation is assumed to not be affected by pore pressure gradients, When the timescale for fluid flow (system size over flow velocity) is large compared to the timescale for deformation of the solid matrix, pore pressure responds to strain in the granular layer (porosity variations), and the magnitude of pressurization is controlled by fluid compressibility; when the fluid flow timescale is short, pore pressure responds to strain rate (rate of porosity variation). When the timescales are comparible, there is a mixed visco-elastic response. Pore pressure may rise even if the system is drained, and if the overall trend of the system is dilative, provided there are periods of abrupt or spatially localized compaction. Fully coupled simulations reveal that the relations between granular deformation, drainage conditions, and pore pressure variations that were derived under the infinite stiffness assumption remain generally valid even when this assumption is relaxed. Coupled simulations also highlight the transition in layer support (from stress-chains through the grains to fluid pressure) that occurs during pressurization. This transition is accompanied by a significant reduction of the shear stress, that potentially allows the development of runaway dynamic slip. Additional insights are offered regarding the time scales controlling dilatancy hardening and the relation between pore fluid pressurization/depressurization and stick-slip motion.
Guarracino, L. & Carrera Ramírez, J. Hydraulic and mechanical effects on tide-induced head fluctuation in coastal aquifer systems Session 2 Tuesday 8  poster file 
Abstract: The study of the dynamic relation between sea tides and coastal groundwater using analytical models has received much attention since the 1950s. Most analytical solutions only consider hydraulic effects and have been derived for aquifer systems that end at a coastline. In the present study, an exact analytical solution is derived to describe tide-induced head fluctuation in an aquifer system that extends a finite distance under the sea. The proposed model incorporates mechanical effects originated by the elastic compression and expansion of the porous formations. The aquifer system consists of a horizontal leaky confined aquifer cover by a semipermeable layer. The analytical solution is a generalization of the solutions obtained by Li and Jiao (Analytical studies of groundwater-head fluctuation in a coastal confines aquifer overlain by a semi-permeable layer with storage, Adv. in Water Resour., 24:565-573, 2001) and Li, Li and Boufadel (The enhancing effect of the elastic storage of the seabed aquitard on the tide-induced groundwater head fluctuation in confined submarine aquifer systems, J. of Hydrology, 350:83-92, 2008), that consider zero and infinite extensions of the aquifer system under the sea, respectively. An important advantage of the derived analytical solution is that it allows to separately compute the mechanical and the hydraulic components of groundwater head fluctuations. The mechanical effect is generated by the fluctuations of sea level that elastic compress and expand the offshore portion of the aquifer system. On the other hand, the hydraulic effect is due to the direct interaction between sea water and groundwater at the submarine outlet of the confined aquifer. The impact of mechanical and hydraulic effects is illustrated through a hypothetical example. The mechanical component of the total tide-induced head fluctuation is significant for aquifers that extend intermediate and large distances under the sea. This component increases the amplitude of the total head fluctuation and tends to synchronize the phase of the sea tide with the phase of the groundwater fluctuation. Then we can conclude that the observation of large fluctuations in coastal boreholes cannot be indicative of good hydraulic connection between the seawater and groundwater. Ignoring the mechanical effect will lead to significant errors in predicting tide-induced head fluctuation and in estimating hydraulic parameters from measurements of groundwater fluctuations.
Gusev, V.A. & Sobissevitch, A.L. Propagation of wideband and shock waves induced by seismic activity in the stratified atmosphere Session 3 Tuesday 8  poster file 
Abstract: Interdisciplinary study of flows of matter and energy in geospheres has become one of the most significant advances in Earth sciences. The actual contribution is the interdisciplinary study of nonlinear acoustics and physical seismology dedicated to wideband and shock wave propagation in the stratified and viscous atmosphere. Seismic activity in Earth's crust induces different process in all geospheres, and one of these processes is the generation of acoustic waves in the atmosphere. Such phenomenon influences on the atmosphere state and can become one of methods of prediction of large seismic events. The main factors to be considered are the stratification of the atmosphere, its viscous and nonlinear effects. Nature origin waves usually have wideband and even shock profiles, so the nonlinear wave equation and evolution equation for such waves are set up. There are obtained the asymptotic profiles for wideband waves, peak waves and durations for simple shock waves. The shock waves propagation through nonperturbed medium is investigated and the sound of its propagation is obtained from the exact solution of generalized Burgers' equation. The heat rate of the atmosphere due to the shock wave propagation is estimated. The efficiency of acoustic streaming induced by shock wave propagation is investigated.
Haben, S.A., Lawless, A.S. & Nichols, N.K. Conditioning of the Variational Data Assimilation problem Session 8 Thursday 10  poster file 
Abstract: Numerical weather prediction (NWP) centres use numerical models of the atmospheric flow to forecast the future weather states from an estimate of the current state of the atmosphere. Data assimilation methods combine observations with the forecast model in order to produce a `best' estimate of the current state, called the analysis. Variational data assimilation (VAR) is popularly used in operational NWP and involves minimising a weighted non-linear least-squares objective function which measures the error between the model forecast and the available observations. The solution of the minimisation is the analysis and is found using an iterative optimisation algorithm. In practice an incremental version of VAR is implemented in many operational centres. This method solves a sequence of linear approximations to the nonlinear least-squares problem. Each approximate linearised least-squares problem is solved using an `inner' gradient iteration method, such as the conjugate gradient method, and the linearization state is then updated in an `outer' iteration loop. The rate of convergence of the inner loop of the VAR iteration scheme and the sensitivity of the analysis to perturbations are proportional to the condition number, that is, the ratio of the largest to the smallest eigenvalue of the Hessian of the linear least-squares objective function. The Hessian is generally assumed to be ill-conditioned and hence operationally the system is preconditioned by transforming the state variables to new variables where the errors are assumed to be approximately uncorrelated. Experimental comparisons have demonstrated that the preconditioning improves the speed of the assimilation scheme. In this work we examine the conditioning of the variational assimilation method theoretically. We derive bounds on the condition number of the Hessian in the case of a single, periodic, spatially-distributed system parameter. We show that the conditioning is sensitive to the length-scale in the correlation structures and confirm that the condition number is generally improved by preconditioning. The theoretical bounds are used to understand how the density and accuracy of the observations affect the conditioning of the preconditioned system. Finally, we examine the conditioning of the Met Office operational VAR method experimentally. We show, using both pseudo and real observations, that our experimental results using this system are consistent with our theoretical results.
Hamiel, Y., Lyakhovsky, V. & Ben-Zion, Y. Nonlinear deformation and the elastic energy of damaged rocks Session 4 Tuesday 8  poster file 
Abstract: Crustal rocks are typically treated as linear elastic material with constant elastic moduli. This assumption is appropriate for rocks with relatively low damage, associated with low concentration of cracks and flaws, and under relatively small strains. However, detailed observational results show that rocks and other brittle materials subjected to sufficiently high loads exhibit clear deviations from linear behavior. In general, nonlinear stress-strain relationships of elastic rocks can be approximated by including higher-order terms of the strain tensor in the elastic energy expression (e.g., Murnaghan model). Such models are successful for modeling rock deformation under high confining pressure. However, values of the third (higher order) Murnaghan moduli estimated from acoustic experiments are one to two orders of magnitude above the expected values of the same moduli estimated from the stress-strain relations in quasi-static rock-mechanics experiments. The Murnaghan model also does not reproduce an abrupt change in the elastic moduli when deformation changes from compression to tension. Such behavior is observed in laboratory experiments with rocks, concrete, and composite brittle material samples. Bi-linear models with abrupt changes of the elastic moduli under stress reversal were suggested based on acoustic experiments (``clapping'' nonlinearity) and in continuum damage mechanics (unilateral damage model). We present the theoretical basis for general second-order nonlinear expression of the elastic potential. We show that a simplified nonlinear model is consistent with bi-linear elastic behavior and accounts for non-linearity for damaged solids even under small strains. We apply the simplified nonlinear model to various laboratory observations, including quasi-static modeling of composite material with different effective moduli under tension and compression; rock dilation under shear; stress- and damage-induced seismic wave anisotropy observed during cycling load of granite samples; and acoustic experiments with shifts of the resonance frequencies in rock samples. Comparison between analytical and numerical simulations and experimental results demonstrate that the suggested expression for the elastic potential for rocks accounts for both quasi-static damage accumulation and nonlinear dynamic effects.
Hayn, M., Panet, I., Holschneider, M. & Diament, M. Multiscale Feature Extraction of Potential Fields Using Poisson Wavelets
Session 5 Thursday 10  poster file 
Abstract: The Earth's geoid gives information about the inner structure and mass distribution of the Earth's interior. It is related to the geodynamical processes operating inside our planet. Analyzing the geoid is a currently often addressed scientific subject, especially since satellite data allow the global analysis of the geoid. Much effort is made on detection and interpretation of structures in the ocean. Reports about undulations of the spatial scale of about 200 km are well established and can be explained by different processes, including secondary convection patterns in the mantle. Also suggestions about undulations of larger scales, from 400 km to 1300 km, exist, but they still need to be confirmed. In our work, we address this issue. We use the method of continuous wavelet analysis in order to detect directional features in the 400-1300 km waveband in the geoid signal of the Pacific ocean. Therefore, we developed directional Poisson wavelets on the sphere and applied them to recent geoids provided from different teams derived from GRACE satellite data, as well as on satellite altimetry derived gravity anomalies and to the Gebco bathymetry. Dominant directions and scales for the different data sets are detected and compared. Dominance filtering shall avoid interpretation of wavelet transform signals that is dominated by noise. After a presentation of the directional analysis method, we discuss our results and their geodynamic implications.
Heifetz, E., Rabinovich, A., Harnik, N., Umurhan, O.M. & Lott, F. Gravity wave instability in terms of vorticity inversion and action-at-a-distance
Session 3 Tuesday 8  poster file 
Abstract: The somewhat counter-intuitive effect of how stratification destabilizes shear flows is reexamined, in what we believe to be the simplest example, in terms of action-at-a-distance interaction between ``buoyancy-vorticity gravity wave kernels''. The setup consists of an infinite uniform shear Couette flow in which the Rayleigh-Fjortoft necessary conditions for shear flow instability are not satisfied. When two stably stratified density jumps are being added, the flow however becomes unstable. At each density jump the perturbation can be decomposed into two coherent gravity waves propagating horizontally in opposite directions. We show, how the instability results from a phase locking action-at-a-distance interaction between the four waves (two at each jump), but can as well be reasonably approximated only by the interaction between the two counter-propagating waves (one at each jump). From this perspective the nature of the instability mechanism is similar to the barotropic and baroclinic ones. Next we add a small ambient stratification to examine how the critical level dynamics alters our conclusions. We find that strong vorticity anomaly is generated at the critical level due to the persistent vertical velocity induction by the edge waves at the jumps. This critical level anomaly acts in turn at-a-distance to decay the edge waves at the jumps. When the ambient stratification is increased, so that the Richardson number exceeds the value of a quarter, this destructive interaction overwhelms the constructive interaction between the edge waves and consequently the flow becomes stable. This effect is manifested when considering the different action-at-a-distance contributions to the energy flux divergence at the critical level. The edge wave interaction is found to contribute toward positive divergence, that is, toward instability, whereas the critical level-edge waves interaction contributes toward an energy flux convergence, that is, toward stability.
Hillers, G. & Ben-Zion, Y. Seasonal variations in observed noise amplitudes above 1 Hz in southern California Session 4 Tuesday 8  poster file 
Abstract: Ambient seismic 'noise' signals are ubiquitous small amplitude ground motions recorded continuously, with variable amplitudes at different frequencies and locations, across the entire observed bandwidth. Noise at periods below several tens of seconds is commonly attributed to ocean-seafloor interactions especially near coastal regions. Energies above 1 Hz are usually attributed to cultural noise with strong diurnal component and wind. Here we show that noise amplitudes at frequencies between 2 and 18 Hz exhibit strong seasonal variations in a broad southern California region. The results are based on seismic records at 30 stations between 240-245 deg longitude and 33-35 deg latitude. The data are recorded between 2002 and 2009 continuously by 3 components stations, most of which are in desert (Mojave and Anza) areas, and the strong seasonal changes are found at all examined stations and frequency bands. Focusing on continuous 6-hour nighttime segments, the 20 or 40 Hz sampled data are bandpass-filtered in 9 frequency bands between 2-18 Hz. The data are median-filtered to reduce the influence of earthquake signals, and integrated to yield half-hourly noise energy estimates. The 6-hour minimum energy value is then converted back to ground velocity and used as representative daily noise level amplitude. Notwithstanding various trends, drifts, and some erratic behavior, a common feature of the resulting time series are annual and sub-annual amplitude changes at essentially all frequencies and all stations, in both the horizontal and vertical components. The strength of amplitude variations differs across stations, and variable phase and amplitude differences occur also at individual stations, for different frequencies, at horizontal and/or vertical channels. Significantly, the amplitude variations show no correlation with distance from the coast, and some particularly clear seasonal variations are seen near topographic features. The large area coverage, spatial pattern of observed amplitudes, existence of signals in both horizontal and vertical components, and strong signals in uninhabited arid regions, imply that the signals are unlikely to originate from ocean waves or variations of ground water. The results may reflect ongoing local failures due to thermoelastic strain generated by interaction of atmospheric temperature changes with large-scale variations in topography, lithology and other surface properties.
de Hoop, M.V., Brytik, V., Cao, Q., Smith, H., Uhlmann, G., van der Hilst, R.D. & Wendt, H. Multi-scale approach to seismic inverse scattering and applications in Earth's upper mantle transition zone Session 9 Monday 7  oral file 
Abstract: We present a program for elastic wave-equation inverse scattering, based on the single scattering approximation, from two interrelated points of view, known in the seismic imaging literature as ``receiver functions'' (unknown source) and ``reverse-time migration'' or RTM (known source). We discuss: (i) the development of an anisotropic polarized-wave equation formulation, (ii) the introduction of an (anisotropic) elastic-wave RTM inverse scattering transform, and (iii) the reformulation of (ii) using mode-converted wave constituents removing the knowledge of the source by introducing the new notion of array receiver functions which generalize the notion of receiver functions. We proceed with presenting a framework for inverse scattering via imaging and partial reconstruction with finite sets of events and seismic stations using multi-scale techniques: To carry out the analysis we make use of higher-dimensional curvelets and introduce associated matrix representations for the component operators that make up the inverse scattering transform. We analyze and exploit the properties of these matrices, and arrive at an efficient method for removing the relevant normal operator. We illustrate various aspects of this research program, integrated with mineral physics and thermo-chemical convection, beneath Hawaii.
Hoteit, I., Luo, X., Pham, D. & Moroz, I. Particle Kalman Filtering: A nonlinear framework for Ensemble Kalman Filters Session 8 Thursday 10  poster file 
Abstract: Optimal nonlinear filtering consists of determining the conditional probability distribution function (pdf) of the state given previous measurements. Once the state pdf is known, one can determine different estimates of the system state, as the minimum variance estimate or the maximum a posteriori estimate. Particle filters (PF) are discrete nonlinear filters that use point-mass representation (Dirac sum) of the state pdf. In practice, these filters suffer from the degeneracy of its particles which causes very often the divergence of the filter. Another discrete solution of the optimal nonlinear filters is based on Gaussian sum representation of the state pdf. This results in a hybrid particle-Kalman filter in which the standard weight-type PF correction is complemented by a KF-type correction for each particle using the associated covariance matrix in the Gaussian sum. We refer to this filter as the particle Kalman filter (PKF). The solution of the nonlinear filtering problem is then obtained as the weighted average of an ensemble of Kalman filters operating in parallel. The Kalman-type correction reduces the risk of ensemble collapse, which enables the filter to efficiently operate with fewer particles than the PF. In this contribution, we present the PKF and discuss how the different EnKF methods can be derived as simplified variants of the PKF. We argue that the (deterministic) Square-Root EnKFs are Gaussian-based filters while the traditional (stochastic) EnKF propagates an approximation of the non-Gaussian pdf of the state. We also discuss approaches to reduce the computational burden of the PKF to make it suitable for complex geosciences applications. These are based on low-rank approximations of the Gaussian sum covariance matrices. Results of numerical experiments with the strongly nonlinear Lorenz-96 model will be presented and discussed.
Igel, H., Fichtner, A., Käser, M., Käufl, P., Pelties, C., Sigloch, K. & Wenk, S. Simulation and inversion of full waveforms for 3-D Earth structures and sources Session 4 Monday 7  oral file 
Abstract: Recently, we have seen the first applications of waveform tomography based on complete solutions to the elasto-dynamic equations in 3-D. Given the developments of computational infrastructure and resources waveform tomography may become the method of choice to determine deep Earth structure given appropriate observational coverage. The core of any inversion is the calculation of synthetic seismograms and there are as many approaches as there are numerical methodologies. Here, we discuss advances using a regular-grid spectral element approach and an algorithm based on the discontinuous Galerkin (DG) method. The latter is based on tetrahedral grids. Due to the fact that the DG approach allows the spatial fields to be discontinuous at the element boundaries the method lends itself to the simulation of rupture along of frictional boundaries (dynamic rupture). Combined with the geometrical flexibility of tetrahedral grids ruptures on arbitrarily shaped faults are possible. In addition to these technical developments we show applications of seismic modeling and inversion based on finite frequency techniques, global wave propagation, and Monte Carlo investigation of model uncertainties.
Zaliapin, I. & Ben-Zion, Y. Seismic clustering and regional physical properties: A statistical analysis Session 4 Tuesday 8  poster file 
Abstract: Discovering genuine spatio-temporal seismicity patterns characterizing local regions beyond the classical large-scale average patterns remains an extremely challenging problem. Of particular importance are patterns that can be related to key physical processes associated with specific properties of faults and the lithosphere. This study addresses two related topics: (i) Non-parametric identification of statistically significant seismic clusters, and (ii) Establishing correlations between seismic cluster patterns and regional physical properties. The first issue is solved using a cluster detection technique of Zaliapin et al. (2008), based on bimodality of the 2-D joint distribution of normalized space and time distances between earthquakes. In the second topic we study two types of patterns. First, we examine the relations between symmetry properties of spatial earthquake patterns along various faults in CA and corresponding local velocity structure images. This is done to test the hypothesis that ruptures on bimaterial faults have statistically preferred propagation directions (e.g., Ben-Zion, 2001). The results indicate strong asymmetric patterns in early-time spatially-close aftershocks along large faults with prominent bimaterial interfaces (e.g., sections of the San Andreas fault), with enhanced activity in the directions predicted for the local velocity contrasts; and absence of significant asymmetry along most other faults. Second, we compare seismic productivity within statistically significant clusters in CA with heat flow and general rock type (crystalline vs. deep sedimentary basins), which serve as proxies for the effective viscosity of the crust (Ben-Zion and Lyakhovsky, 2006). We find that (i) relatively cold regions with crystalline rock in the seismogenic zone have high aftershock productivity and low foreshock productivity, and vice versa (ii) regions with high heat flow and deep sedimentary basins have increased foreshock activity and reduced aftershock productivity. The results demonstrate that seismicity patterns do not follow universal power law statistics in all space-time domains. Assuming the observed asymmetric properties of seismicity reflect asymmetric properties of earthquake ruptures, and earthquake productivity reflects the local seismic potential of the crust, the discussed methodology and results can be used to develop refined estimates of seismic shaking hazard associated with individual fault zones and regions.
Ionescu, I.R. & Cazacu, O. Onset and dynamic shallow flow of a viscoplastic fluid. Applications to dense avalanches Session 2 Tuesday 8  poster file 
Abstract: The starting point is the 3-D setting of a sheet flow of a viscoplastic fluid. The material constitutive law may include two plasticity (flow/no flow) criteria: Von-Mises (Bingham fluid) and Drucker-Prager (Mohr-Coulomb). A stress analysis is used to deduce a Saint-Venant type asymptotic model for small thickness aspect ratio. The 2-D (asymptotic) model (constitutive law), which relate the averaged plane stress to the horizontal rate of deformation, is not anymore incompressible and represents the plane stress ``projection" of the initial 3-D viscoplastic model. The ``safety factor" (limit load) is introduced to model the link between the yield limit (material resistance) and the external forces distribution which generate the shallow flow of the viscoplastic fluid from a rest configuration. The discontinuous velocity domain splitting method, developed in [I.R. Ionescu, E. Oudet, Discontinuous velocity domain splitting method in limit load analysis, Int. J. Solids. Structures, 2010, in press], is used to evaluate the safety factor and to find the onset of an avalanche flow. A mixed finite-element and finite-volume strategy is developed. Specifically, the variational inequality for the velocity field is discretized using the finite element method while a finite volume method is adopted for the hyperbolic equation related to the thickness variable. To solve the velocity problem, a decomposition-coordination formulation coupled with the augmented lagrangian method, is adapted here for the asymptotic model. The finite volume method makes use of an upwind strategy in the choice of the flux. Several boundary value problems, modeling shallow dense avalanches, for different visoplastic laws are selected to illustrate the predictive capabilities of the model.
Isakov, V. Some Numerical Methods for Gravimetric Prospecting Session 9 Tuesday 8  poster file 
Abstract: We discuss recovery of a perturbing body from (gradient or higher order gradient) of its gravimetric potential given on a part of the external surface. For numerical reasons we mainly handle the plane case and prescribe the data on an interval of x-axis, while the unknown body is assumed to be a polygon in the lower half-plane. To understand limits of resolution we first give a detailed singular value decomposition analysis of a linear operator which continues the field to outside of a disk including unknown body. We indicate how many parameters are possible to determine depending on the distance from measurement site and size of unknown body. Then we describe two numerical algorithms using complex variables theory which are based on 1) direct (nonlinear) minimization technique and 2) Prony type method. We give many examples of numerical reconstruction for various shapes of unknown inclusion and distance to the measurement site, which show possibilities and limitations of inverse gravimetry. Finally we give a new uniqueness result for recovery of a single layer mass distribution and of lower part of an inclusion (modeling ice with snow layer on top) based on the classical Novikov's orthogonality method (and described in the author's book ``Inverse Source Problems'', AMS 1990).
Ismail-Zadeh, A., Honda, S. & Tsepelev, I. Quantitative reconstruction of lithosphere subduction using assimilation of geophysical data Session 8 Friday 11  oral file 
Abstract: This presentation intends to provide an answer: whether the evolution of the descending lithosphere can be restored quantitatively. To restore mantle structures and flow in the geological past, mathematical and numerical techniques for inverse retrospective problems should be used to constrain the initial (past) conditions for the temperature and velocity in the crust and mantle from present seismic, heat flow, geodetic and some other observations. The initial conditions so obtained can then be used to run forward models of lithosphere/mantle dynamics coupled with the models of surface processes in order to restore the evolution of the descending lithosphere. The basic principle of the inverse retrospective problem is to consider the initial condition as a control variable and to optimize the initial condition in order to minimize the discrepancy between the present observations and the model solution. We present a quasi-reversibility technique for geodynamic data assimilation and discuss applicability of this technique to restoration of a descending lithosphere in the presence of mantle phase transformations.
Isola, I. & Mazzarini, F. Model for local fluid circulation in the upper brittle crust: insights from veins' thickness distribution Session 2 Tuesday 8  poster file 
Abstract: Migration of fluids in the crust is basilar for geologic processes such as regional metamorphism and formation of hydrothermal, magmatic, volcanic systems and ore deposits. The bulk permeability of rocks is greatly enhanced by fractures, where their geometric properties (attitude, length, density and aperture) are important to define the hydraulic behaviour of the network (e.g. Bonnet et al., 2001). Veins record fluids' circulation in fractures at depth in the crust and fluid/rock interaction during regional metamorphism (Oliver, 1996) and also in geothermal systems (e.g. Bertini et al., 2006). In southern Tuscany (Italy), well-exposed Oligocene-Early Miocene sandstones hosting vein systems provide insight into the role of pore fluid and the stress state at the time of vein formation (Mazzarini et al., 2004; Mazzarini and Isola, 2007). In particular, the analysis of veins' thickness (t) allowed us to compute the average transmissivity of the veins in two different sites characterised by different veins' thickness distributions (power-law and negative exponential). The computed transmissivity for the veins is in the range $10^-3 - 10^-1 m^2s^-1$, with higher values attained by the veins having power-law thickness distribution.

Bertini, G., Casini, M., Gianelli, G., Pandeli, E., 2006, Terra Nova 18, 163-169.
Bonnet, E., Bour, O., Odling, N., Main, I., Berkowitz, B., Davy, P., Cowie, P., 2001, Reviews of Geophysics, 39, 347-383.
Mazzarini, F., Corti, G., Musumeci, G., Innocenti, F., 2004, In: Breitkreuz, C., Petford, N., (Eds.), Physical geology of High-level Magmatic Systems. Geological Society, London, Special Publication, 234, 151-161.
Mazzarini, F., Isola, I., 2007, Journal of Structural Geology, 29, 1386-1399.
Oliver, N.H.S., 1996, Journal of Metamorphic Geology, 14, 477-492.

Isola, I., Mazzarini, F. & Molli, G. Structural analysis of a percolating fracture network in karst systems: the Antro del Corchia Cave, Alpi Apuane, Italy Session 2 Tuesday 8  poster file 
Abstract: Modelling evolution of Karst aquifer in carbonate sequences requires to model flow and dissolution processes within percolation network that exploited narrow joints and bedding planes successively widened or sealed by carbonate dissolution/concretion processes. The process is very rapid with respect to geological times being generally less than 50 ky the time required for developing an integrated karst network (e.g. Bakalowicz, 2005; Siemers and Dreybrodt, 1998). The first step in modelling such a complex system is to analyse the actual fracture systems in terms of fracture attitude, distribution and conductivity. We present the first results of a structural study of a portion of the Antro del Corchia Cave in Alpi Apuane (northern Apennine, Italy) which consists of a large karst systems developed in the late Triassic-Jurassic metamorphic sequence of the Alpi Apuane (Molli and Vaselli, 2006; Piccini et al., 2008). Detailed structural analysis of brittle deformation and of its relationships with concretions along a nearly N-S trending transect in the Antro del Corchia Cave allowed us to determine which part of the actual fracture network belongs to the backbone of the percolating network. Fractures clearly sealed by concretions, fractures with clearly evidences of fluid circulation and dry fractures cutting across concretions have been classified in the cave. These observations are then compared with the geology of the area and with the observed fracture network at the surface, a level about 400 m above the analysed transect in the cave.

Bakalowicz, M., 2005, Hydrogeology Journal, 13, 148-160.
Molli, G., and Vaselli, L., 2006, Geological Society of America, Spec. Pub., 414, 79-93.
Siemers, J. and Dreybrodt, W., 1998, Water Resources Research, 34, 409-419.
Piccini, L., Zanchetta, G., Drysdale, R.N., Hellstrom, J., Isola, I., Fallick, A.E., Leone, G., Doveri, M., Mussi, M., Mantelli, F., Molli, G., Lotti, L., Roncioni, A., Regattieri, E., Meccheri, M., and Vaselli, L., 2008, International Journal of Speleology, 37, 153-172.

Kaypak, B. & Venedik, G. Three-Dimensional Seismic Velocity Structure in The Denizli Geothermal Region, Western Turkiye Session 4 Tuesday 8  poster file 
Abstract: The Denizli basin and its surroundings are one of the regions that have high seismic activity and rich geothermal fields in the west Anatolia. In the first half of a year of 2000, it was observed an increasing seismic activity around the Denizli basin without a large mainshock occurrence. A temporary seismic network consisted of 28 stations were deployed to observe the seismic activity in the region by a Turkish research institute (TÜBITAK-MAM-EMSI). In this study, three-dimensional Vp and Vp/Vs structure of the Denizli basin have been determined by using the travel times of the collected data. Firstly one-dimensional inversion schema was performed to stabilize initial velocity model and to have a more reliable hypocentral locations. Then an iterative and simultaneous three-dimensional inversion procedure was carried out to obtain 3-D high-resolution seismic velocity models. Furthermore to assess the solution quality of our inversion, we conducted a series of resolution tests.
We concluded high-resolution 3-D VP and VP/VS seismic velocity models for the upper 20 km of the crust beneath the Denizli basin and surroundings and interpret the results in the context of known geologic and tectonic units. The resulting VP models define the geometry of the basin and sediment thickness (VP $ 2 km/s) as well as regions of anomalous velocity. The VP/VS models help to constrain the location and geometry of the faults, zones of weakness and fluid saturated formations with high pore pressure zones.
Kopp, R.E. & Simons, F.J. Assessing the history of global sea level from noisy and incomplete observations of local sea level Session 7 Tuesday 8  oral file 
Abstract: The last interglacial (LIG) is an analogue for contemporary global-warming scenarios. Geological records show local sea level (LSL) in the LIG generally higher than today, but the mapping from global sea level (GSL) to LSL is modulated by the solid Earth. LSL is a function of both time since and location relative to the source of the melting events that define GSL. We determine the posterior pdf of GSL (and ice volume) through time in the LIG, from a database of sparse and scattered, noisy measurements of LSL. Noise is present in both the dependent and independent variables, sea level and age. Both posterior and prior pdfs are assumed to be Gaussian. By Gaussian process regression, we might approach this as a space-time interpolation/inversion problem, but absent large numbers of precise measurements, we need simulation to obtain the covariance of LSL and GSL. The covariance expresses the Earth system response through the sea-level equation. We use the solver developed by Mitrovica, modelling eustatic, gravitational, deformational, and rotational effects of melting ice sheets. To construct a prior, we average LSLs and GSL over multiple alternative ice sheet histories run through this forward model. The histories themselves are sampled from two underlying pdfs: one for global ice volume from isotopes and another for individual ice sheet volumes based upon random perturbations of last-glacial-maximum-to-present models. Two more Gaussian terms are added to approximate thermosteric effects and small contributions from unaccounted sources such as small mountain glaciers. To construct the posterior we start by estimating it at the points included in our database, and then interpolate to calculate values elsewhere, as in the covariance we have incorporated the full physics of the problem without constraining it to fit a particular forward model. First, we correct measurements of LSL for tectonics. Second, we employ Gaussian process regression to estimate the true sea levels. Third, we use the Metropolis-Hastings algorithm to draw a new sample of the ages, based upon the measured ages and the current estimate of the true sea levels. Repeating this Gibbs sampling sequence many times we obtain the posterior pdf for LSL and GSL in a way that satisfies the measurements. With 95% probability GSL peaked at least 6.6 m higher than today during the LIG. This highlights the vulnerability of ice sheets to even relatively low levels of sustained global warming.
Kurzon, I., Lyakhovsky, V. & Navon, O. Visco-elastic magma - Fragmentation criteria revisited Session 1 Thursday 10  poster file 
Abstract: We present a visco-elastic bubble growth model, accounting for viscous and elastic deformations and for volatile mass transfer between bubbles and melt. We define the borders between previous bubble growth models accounting for incompressible viscous melt, and our new model that also accounts for melt compressibility and elastic deformation. Elastic deformation is most prominent for magma of small vesicularities, where the growth regime depends on the shear modulus. Following an instantaneous pressure drop, elastic deformation of melt with high shear modulus is negligibly small; high gas pressure is preserved in the bubbles and their growth follows an exponential viscous growth regime. For magma of low shear modulus bubbles initially expand due to elastic deformation of the surrounding melt; gas pressure falls quickly and growth follows a square-root diffusive solution. Our model provides all the elastic components (stresses, strains and strain rates) required for defining criteria for failure and magma fragmentation. We examine two failure criteria, a stress related criteria based on the internal friction and the Mohr-Coulomb failure theory, and a strain related criteria based on fibre elongation experiments. We argue that both criteria are equivalent if we consider their dependence on the shear modulus and effect of the modulus on magma rheology.
Kurzon, I., Lyakhovsky, V., Navon, O. & Chouet, B. Amplification of pressure waves in supersaturated bubbly magma Session 1 Friday 11  oral file 
Abstract: We study the interaction of acoustic pressure waves with an expanding bubbly magma. The expansion is the result of bubble growth during or following magma decompression and leads to two competing processes that affect pressure waves. On one hand, growth in vesicularity leads to increased damping and decreased wave amplitudes, and on the other hand, a decrease in the effective bulk modulus of the bubbly mixture reduces the wave velocity and leads to wave amplification. Under certain conditions, the amplification is strong enough to account for a significant reduction in damping and may even lead to amplification. Thus, some of the energy released during growth is transformed into acoustic energy. We first examine this phenomenon analytically to identify conditions under which amplification is possible. These conditions are further examined numerically to shed light on the frequency and phase dependencies in relation to the interaction of waves and growing bubbles. Amplification is possible at low frequencies and when the growth rate of bubbles reaches an optimum value for which the wave velocity decreases sufficiently to overcome the increased damping of the vesicular material. We present two amplification phase-dependent effects: 1) a tensile-phase effect in which the inserted wave adds to the process of bubble growth, utilizing the energy of the growing bubble and therefore converting a large proportion of this energy into additional acoustic energy; and 2) a compressive-phase effect in which the pressure wave works against the growing bubbles and dissipates a large amount of the acoustic energy during the first cycle, but later gains enough energy to amplify the second cycle in comparison to the first. The two effects provide additional new mechanisms for the amplification phase seen in Long-Period (LP) and Very-Long-Period (VLP) seismic signals originating in magma-filled cracks. For example, when bubbly magma expands within a crack it builds up pressure. Upon reaching a critical pressure threshold the crack fails, allowing the release of gases through the permeable part of the magma. The resulting collapse of the permeable network allows new supply of magma and further degassing, rebuilding the pressure and initiating the next cycle. Waves emitted during the decompression stage may move freely through the expanding magma due to its lower damping and, in some cases, they may actually be amplified.
LaCasce, J.H. & Isachsen, P.E. The dynamics of the Antarctic Circumpolar Current Session 3 Tuesday 8  oral file 
Abstract: The Antarctic Circumpolar Current (ACC) is the largest current in the world, and the only truly circumpolar one. As such, our thinking about the ACC has been strongly influenced by theories of the Jet Stream in the atmosphere, which also is circumpolar. These theories are fundamentally nonlinear, as the Jet Stream structure is determined eddy fluxes of heat and momentum. There are linear theories of the ACC as well, similar to those proposed for the Gulf Stream. But these have received much less attention.
However, the linear solutions are very revealing with regards to the dynamics. The solutions depend strongly on the ``geostrophic contours'' in the Southern Ocean. These are latitude lines, if the bottom is flat, and more complex structures with stratification and topography. If the contours are unblocked by lateral barriers (continents), the ACC follows the geostrophic contours and is typically very strong. But if the contours are blocked, the dynamics are more like that in the Gulf Stream models, with western boundary currents. The solutions in the latter case bear a strikingly resemblance to the actual ACC. And the transport predicted by the model is within errors of the observed value.
This is not to say that eddies are unimportant in the Southern Ocean, as the linear models do not predict the vertical structure of the flow. But the success in reproducing the current path suggests the linear solutions can be used to understand how the ACC responds to changes in wind stress and topography.
Langmann, B., Werning, M. & Hort, M. Regional numerical modelling of volcanic ash atmospheric dispersion and deposition after the eruptions of Mt. Pinatubo and Kasatochi Session 6 Tuesday 8  poster file 
Abstract: The violent nature of volcanic eruptions involving steam can fragment magma and solid rocks surrounding the volcano into particles as small as the size of aerosols. The very fine ash particles may be carried for hundreds of miles before settling onto land or into the ocean. Fierstein and Nathenson (1992) describe several ways to determine the deposited volcanic ash fall volume and mass. These, however, only work if ash is deposited on land. They are unsuitable when volcanic ash is mainly deposited into the ocean as it was the case after the eruptions of Kasatochi in 2008 or Pinatubo in 1991. Another possibility to estimate the deposited ash volume is to use three-dimensional atmospheric models which include transport and removal processes of volcanic ash to estimate the atmospheric distribution and fall-out of volcanic ash. Here we use the three-dimensional Eulerian atmosphere-chemistry/aerosol model REMOTE (Langmann et al., 2008) to simulate the distribution of volcanic and its deposition after the eruptions of Mt. Pinatubo and Kasatochi. The aerosol dynamic and thermodynamic module M7, described in detail in Vignati et al. (2004), provides the framework for the volcanic ash size determination using log-normal size distributions. Our model simulations suggest that sedimentation was the most efficient removal process for volcanic ash mass after the eruption of Kasatochi with 70% of the total mass removed out of the atmosphere at ground level, followed by wet deposition (23 and dry deposition (7. The situation was different after the eruption of Pinatubo as a typhoon passed by and wet deposition of volcanic ash was even more important than in the case of Kasatochi. Here a series of sensitivity studies will be presented focussing on the contribution of the three ash removal processes out of the atmosphere: dry deposition, wet deposition and sedimentation.
This work is supported through the Cluster of Excellence 'CliSAP' (EXC177), University of Hamburg, funded through the German Science Foundation (DFG).
Fierstein, J. & Nathenson, M. (1992). Bull. Volcanol., 54, 156-167.
Langmann, B. et al. (2008). Atmos. Chem. Phys., 8, 1591-1607.
Vignati, E., Wilson, J. & Stier, P. (2004). J. Geophys. Res., 109, doi:10.1029/2003JD004485.
Lev, E. Extracting lava velocity and rheology from computer-vision analysis lava flow videos Session 1 Thursday 10  poster file 
Abstract: Effusive lava flows present a considerable hazard to property and infrastructure. It is therefore important to understand the processes controlling lava flow advance and channelization, of which lava rheology is among the most important factors. Unfortunately, to this time, in-situ measurements of lava velocities are sparse in both time and location. To advance our understanding of lava flow emplacement, we utilize techniques adopted from computer vision and image processing and extract surface velocities of flowing lavas. We analyze image sequences obtained from every-day videos. Thanks to the high density of the velocity field, we are able to calculate cross-channel velocity profiles, and to invert for the first-order characteristics of lava rheology. In addition, the techniques we use can be also used for studying the motion of glaciers, landslides, and pyroclastic flows, even if the only available footage is not ideal.
Longo, A., Papale, P., Cassioli, A., Saccorotti, G., Barsanti, M., Montagna, C.P., Vassalli, M. & Giudice, S. GALES: a consistent finite element numerical library for the simulation of underground magma dynamics Session 1 Thursday 10  poster file 
Abstract: The simulation of underground magma dynamics is a challenging, highly non-linear numerical problem. The finite element C++ parallel code GALES, developed to this purpose, solves mass, momentum and energy equations for multiphase multi-component homogeneous gas-liquid ($ crystals) mixtures with P-T-X-dependent properties. GALES is based on Galerkin-Least Squares Discontinuity Capturing stabilization with space-time discretization. The formulation is stable, robust and accurate, suited for computation of complex dynamics as well as for the solution of several types of conservation equations such as those governing rocks dynamics (under development) or Darcy's flows. Basis functions and the mapping between grid and reference elements depend on space and time, so that mesh deformation is implicit in the method. The use of pressure-primitive variables allows solution for compressible to incompressible flows, since the matrix of change from conservative to primitive variables is well-behaved in the incompressible limit. The stabilizing terms scale with the accuracy of the solution, approaching zero when the solution converges. The Least Squares term adds a streamline weighted transport of information that prevents pure Galerkin solution breakage, reducing and concentrating numerical oscillations in a narrow neighborhood of the internal discontinuities and boundary layers. The Discontinuity Capturing term controls these oscillations. Performed 2D numerical simulations investigate magma dynamics for Stromboli, Etna, and Campi Flegrei magmatic systems characterized by different geometrical complexities and magma composition/properties. The results show complex time-space patterns of magma convection and mixing, previously not described, and oscillatory pressure trends related to the competing effects of buoyant magma rise and dense magma sinking. Integration of the local mass contributions to the gravitational field, and solution of the Green's functions describing pressure wave propagation through the surrounding rock system, provide time-space distribution of ground deformation and gravity anomalies at the Earth's surface and the seismic signals associated with magma motion below the volcanoes. The relationships emerging between deep magma dynamics and recorded geophysical signals lead to better definition of early warning systems and forecasting of the short-term volcanic hazard.
Lopes, O.F. & de Souza, J. A Super Volcano in Southern Brazil Session 5 Thursday 10  poster file 
Abstract: A new relationship between the geological process in a normal volcanic eruption is proposed nere to explain the super volcano of the locality named Caolin Paraná. The geological history of this area is produced by a normal evolution of a thick succession of debris flow deposits of a remote and unknown volcano, the diagenetic and possible hydrothermal transformation, with the massive elimination of the matrix phase (ash) and concentration of roundly and little (with mean diameter of 20 cm) fragments of 32 different rock species, subdivided in numerous rock facies. One important cementation phase, essentially of iron oxide nature, make the whole debris flow a consolidated and resistant material, sufficiently to obstruct and interrupt a evolution of the normal volcano, may be for years or decades. This work, with many examples, actually showed as fragments with 20 or 25 thousands tones, to 10 or 20 tones (away the volcanic centre), try to establish the extension of the volcanic eruption itself (with an originally destruction of the primary volcanic structure, with no ramains now), and the environmental catastrophe that influences the earth atmosphere, possibly in a planetary scale, in the Early Phanerozoic Eon. The geodynamic model of known volcanic eruptions needs, in this context, be completely revised and seems sure that new processes need imperatively to be proposed, as answers of some questions: 1)What is the height of the volcanic column in a super eruption? 2)The energy balance exclude the pyroclastic rocks of this kind of eruption? 3)What is the area of influence (with defined ash isopac) of a super eruption? 4)What are the transport mechanisms (aerial or on/over the ground) of the giant fragments? 5)The other known examples of the same evolution, in the area, may have counterparts in other regions of the world? This seems just the beginning of a new era of study on particular kinds of volcanic evolution.
Loris, I., Simons, F.J., Voronin, S., Daubechies, I.C., Nolet, G., Judd, S.J., Fornasier, M. & Vetter, P.A. A new approach to global seismic tomography that promotes sparsity with a new three-dimensional wavelet transform in spherical geometry Session 4 Tuesday 8  poster file 
Abstract: Seismic wavespeed models of the Earth are routinely parameterized in terms of spherical harmonics, networks of tetrahedral nodes, rectangular voxels, or spherical splines. However, there as of yet few approaches to Earth model parameterization by wavelets, on the three-dimensional ball. Much is to be gained from procedures that do this efficiently, as (1) the multiresolution character of a wavelet basis allows for the models to be represented with an effective spatial resolution that varies as a function of position within the Earth. Furthermore, inversion schemes that are formulated in terms of wavelets can (2) exploit recent theoretical and numerical advances by which the most sparse solution vector, in wavelet space, is found through iterative minimization of a combination of the L2 (to fit the data) and L1 norms (to promote sparsity in wavelet space). With the continuing increase in high-quality seismic data, our focus should also be on (3) numerical efficiency and the ability to use parallel computing in constructing the model. In this presentation we propose a new wavelet basis on the sphere that has these three goals in mind. To form the numerical grid we begin with a surface tesselation known as the cubed sphere, a construction popular in fluid dynamics and computational seismology, coupled with a semi-regular radial subdivison that honors the major seismic discontinuities between the core-mantle boundary and the surface. By this mapping we transform the entire volume of the mantle into six portions. In the new variables, these chunks correspond to rectangular boxes with Cartesian coordinates. Standard algorithms can then be used to perform the wavelet transformation (or any other transformation) in each of the six bounded volumes. We choose a 4-tap discrete wavelet transform in the angular variables and use the Haar transform in the radial direction, with preconditioning and special boundary filters. We highlight the benefits of this construction and apply it to analyze the information present in several published seismic compressional-wavespeed models of the mantle, paying special attention to the statistics of wavelet coefficients across scales, and focusing on the likely gains of future inversions of finite-frequency seismic data using our new wavelet approach.
Luon, X. & Hoteit, I. Ensemble Local H-infinity filter for data assimilation Session 8 Thursday 10  poster file 
Abstract: In this work we consider the H-infinity filter and its variants for data assimilation. The H-infinity filter is derived by minimizing the loss in the worst case, a criterion different from the minimum variance used in the Kalman filter. Thus by design, the H-infinity filter is more robust than the Kalman filter, in the sense that the estimation error in the H-infinity filter has a bounded growth rate with respect to the uncertainties in estimation, while the estimation error of the Kalman filter does not possess such a guarantee. Despite the difference between the criteria in deduction, it turns out that the H-infinity filter bears a very similar structure to that of the Kalman filter. In fact, it can be shown that the Kalman filter is a special case of the H-infinity filter, while the H-infinity filter algorithm itself can be constructed based on that of the Kalman filter.
The original form of the H-infinity filter contains global constraints (in time), which may be inconvenient to solve for recursive filters. Therefore here we introduce a variant, which requires solving some local constraints instead, hence will be called the local H-infinity filter. Furthermore, analogous to the ensemble Kalman filter (EnKF), we also propose the concept of ensemble local H-infinity filter (EnLHF). We give the general form of the EnLHF, and discuss some of its special cases. We show that these special cases are similar, or even identical, to some of the EnKF methods with covariance inflation in the literature. We also use two numerical examples to illustrate the relative robustness of the EnLHF in comparison to the EnKF method (without covariance inflation), in the presence of larger than expected uncertainties in data assimilation.
Lupi, M., Geiger, S. & Graham, C. Fluid-induced seismicity in volcanic and hydrothermal systems: numerical studies of the Tjörnes Fracture Zone, Iceland Session 2 Tuesday 8  poster file 
Abstract: Pre-, co- and post-seismic fluid flow and fault behaviour have in the seismically active region of the Tjörnes Fracture Zone, North of Iceland are the subject of this study. We run fluid flow (heat and mass transport) and fault-rupture simulations of a geological model which is based on multichannel seismic reflection data. We reconstruct geological structures (i.e. faults, unconformities and lens-shape layers) of the first 10 km of the Tjörnes Fracture Zone at great detail. Numerically we resolve geological processes occurring at vastly different time-scales. The conceptual model of the toggle switch mechanism has been applied to analyse seismicity-induced fluid flow and explain key geophysical and geochemical observations in the Tjörnes Fracture Zone. The data highlight the occurrence of four distinct over-pressured areas in the crust of the Tjörnes Fracture Zone, consistently with geophysical observations. At depth critical over pressures are reached while low over pressures prevail between 7 km depth and the top of the basement. Faults often have critical over pressures and are close to failure. During failure, faults can release over pressures in less than five minutes, which is accompanied by a sudden in permeability of over seven orders of magnitude. At this time, short-lived but extreme flow rates of more than 0.01 $m s^-1$ occur. During the post-seismic phase (i.e. after failure), as the over pressure dissipates and the normal stress increases, fault permeabilities decay back to their pre-failure values in 2 to 3 years, consistently with geochemical data. Our model suggests that the permeability of seismically active hydrothermal systems is constantly evolving and can lead to short-lived catastrophic flow events. Our results show why hydrogeochemical monitoring of faults may be an important tool to improve the confidence in earthquake forecasting.
Lupi, M., Geiger, S., Thordarson, T., Carey, R.J. & Houghton, B.F. Explaining the Plinian-phreatoplinian shift during the 1875 Askja volcano eruption by coupling geological and numerical techniques Session 1 Friday 11  oral file 
Abstract: The Askja Volcano is part of the larger Dyngjufjöll complex which is located in Central Iceland. During the 1874-76 a large volcano-tectonic episode took place on the Northern rift zone and the 1875 eruption of the Askja Volcano is part of it. The explosive rhyolite eruption of the Askja caldera, one of the very few eruptions showing both phreatoplinian and Plinian styles and the only well documented by eye witnesses, is the subject of this study. The eruption began the 28th of March (9 pm) with a subplinian event (phase B) which lasted for ~1 hour. In the early morning of the 29th, after a pause of ~ 6.5 hours, the eruption continued with a phreatoplinian phase (phase C1), which lasted ~1 hour. This phase included emplacement of dilute density currents (phase C2), which became dryer with time. At 7 am the Plinian phase D commenced and lasted for about 5-6 hours. The vent position, which migrated throughout the eruption, has been constrained by previous studies. Aside from a small pond, no standing water which could have caused the phreatoplinian phase was present at that time in the caldera and therefore the source of external water driving the phreatoplinian phase is uncertain. We have undertaken a study to test if groundwater residing in the lava pile that partly fills the main Askja caldera could have been the provider of the external water involved in the phreatoplinian phase. The key questions are: Does the intra-caldera groundwater reservoir hold enough water? Can the water be transported fast enough to the vent site? The intra-caldera lava flows are strongly jointed and fractured and in order to quantify the fracture pattern and the permeability it provides we have measured and mapped out the fracture patterns of several lava flows. A discrete fracture modelling technique is used to compute the corresponding permeability and porosity. A 3D digital model of the Askja caldera that represents the pre-eruptive topography and hydrostratigraphy has been reconstructed. The calculated porosities and permeabilities provide factual data as model input parameters. Results show that the eruption resulted in ultra-fast radial flow of groundwater towards the conduit and provides sufficient water flux to drive the phreatoplinian phase of the 28-29 March Askja eruption. Furthermore, the model also produces rapid draw-down of the intra-caldera groundwater which explains the drying-out of the phreatoplinian phase depicted by the dilute density current deposits.
Lyakhovsky, V. Damage rheology model: from local quasi-static to non-local dynamic formulation Session 4 Monday 7  oral file 
Abstract: Many damage models attempt to quantify material deformation beyond the reversible elastic regime using continuum and discrete approaches, scalar and tensorial damage state variables, linear and nonlinear elasticity, and other categories of ingredients and formulations. In previous studies, we developed a nonlinear visco-elastic continuum damage rheology for evolving elastic properties of rocks sustaining irreversible brittle deformation. The model accounts for three general aspects of brittle rock deformation: (1) Dependency of the effective elastic moduli on the existing crack density (damage) and loading conditions; (2) Evolution of the crack density and elastic moduli with the ongoing deformation; (3) Macroscopic brittle instability at a critical level of damage. The previous formulation assumes local relation between the energy function and material damage, leading to complete localization during failure, and may be classified as ``quasi-static'' since it does not allow a detailed analysis of dynamic events associated with macroscopic instability. We present a new non-local and dynamic formulation of the damage rheology model that eliminates these major shortcomings. The model accounts for strong micro-crack interaction in a highly damaged region prior to the total macroscopic failure by incorporating the spatial derivative of the damage as an additional state variable. This leads to a definition of an internal length scale that does not exist in elasticity and in local damage models. The new model predicts a finite width of a newly-created highly-damaged zone. Accounting for the dynamics require two additional modifications: (1) Incorporating the strain rate tensor in the state variables provides a thermodynamic formulation for the Kelvin-Voigt visco-elastic model and produces some damping during the dynamic processes; (2) We suggest that the loss of convexity signifies a phase transition from a highly damaged solid to a granular material. In the vicinity of the critical level of damage, we expect a formation of a mixture zone of partially solid and granular media similar to the mushy region discussed in a context of the Stefan problem of classical thermodynamics. We discuss briefly key relevant aspects of granular mechanics and suggest a simplified approach to characterize the transition between damaged solid and granular material. We demonstrate that the new damage formulation allows a self-consistent analysis of the dynamic rupture process.
Maimon, O., Lyakhovsky, V., Melnik, O. & Navon, O. Factors controlling propagation of a dyke filled with gas-saturated magma Session 1 Thursday 10  poster file 
Abstract: We present mathematical and numerical model of the propagation of a dyke filled with volatile-saturated magma and a gas cap at its upper part. The magma flows due to its buoyancy and the overpressure at the dyke origin. The model solves for the elastic stresses in the host rocks, the shape of the dyke, the rate at which the tip of the dyke propagates, and for the two-phase (gas and melt) channel flow. The numerical code allows studying various regimes of dyke propagation and the importance of the controlling parameters. Following a sudden increase in the pressure at the dyke origin, the pressure wave propagates until it reaches the magma front. Then the dyke front propagates at a rate that depends on the sub-critical crack propagation of the tip, the gas mass in the cap cavity, the gas flux from the bubbly magma to the cap and on other physical parameters of the system. Pressure decrease in the magma during the flow leads to nucleation and growth of gas bubbles and gas filtration. When pressure variations are small, we neglect gas transport, assuming a constant gas mass in the cap. In this case the propagation of the dyke is self-similar. The rate of propagation is controlled by the rate of fracturing at the tip, and the rate of magma flow in the dyke. In the case of relatively low magma viscosity and high energy needed to fracture the host rock, the rate of the dyke propagation is controlled by the rate of fracturing at the dyke tip (fracture-controlled regime). The dependence of the dyke propagation rate on the amount on the gas in the cap is very weak. When the energy to fracture the host rock is small, the rate of dyke propagation is controlled by the rate of the magma flow (magma-controlled regime). In this case the propagation rate is very sensitive to the gas mass in the cap cavity. The gas-filled cap cavity opens the dyke in front of the magma and leads to higher propagation rates than predicted by models ignoring the existence of a gas cap. The maximal rate of propagation is obtained at the transition between these regimes. The propagation of the dyke with a small amount of gas in the gas cap, the rate of propagation is fracture-controlled. At later stages, the amount of gas in the gas cap increases, the crack growth is faster than magma channel flow and the dyke propagation is in the transitional or magma-controlled regime. If the gas mass is high enough, the cap may separate and form a distinct unconnected balloon that propagates on its own.
Mair, K., Jettestuen, E. & Abe, S. Characterising contact force networks in 3D sheared granular fault gouge Session 4 Tuesday 8  poster file 
Abstract: Faults often exhibit accumulations of granular debris, ground up to create a layer of rock flour or fault gouge separating the rigid fault walls. Numerical simulations and laboratory experiments of sheared granular materials, suggest that applied loads are preferentially transmitted across such systems by transient force networks that carry enhanced forces. The characterisation of such features is important since their nature and persistence almost certainly influence the macroscopic mechanical stability of these systems and potentially that of natural faults.
3D numerical simulations of granular shear are a valuable investigation tool since they allow us to track individual particle motions, contact forces and their evolution during applied shear that are difficult to view directly in laboratory experiments or natural fault zones. In characterising contact force distributions, it is important to use global structure measures that allow meaningful comparisons of granular systems having e.g. different grain size distributions, as may be expected at different stages of a fault’s evolution. We therefore use a series of simple measures to characterise the structure, such as distributions and correlations of contact forces that can be mapped onto a force network percolation problem as recently proposed by Ostojic and coworkers for 2D granular systems. This allows the use of measures from percolation theory to both define and characterise the force networks. We demonstrate the application of this method to 3D simulations of a sheared granular material. Importantly, we then compare our measure of the contact force structure with macroscopic frictional behaviour measured at the boundaries of our model to determine the influence of the force networks on macroscopic mechanical stability.
Makedonska, N., Goren, L., Sparks, D.W. & Aharonov, E. Friction versus dilation revisited: insights from theoretical and numerical models.
Session 4 Monday 7  oral file 
Abstract: The intimate relation between apparent friction of shearing granular layers and their dilation was already discussed by Mead in 1925. Motivated by the importance of this connection to the frictional strength of geological faults and to earthquake generation, many laboratory and numerical experiments on sheared granular layers investigated the relation between the apparent friction and dilation. Although the nature of the connection is not very well established, apparent friction is often cited to be a linear combination of two contributions: 1. The surface friction coefficient of the grains and 2. the dilation rate (which is equivalent to the ratio of the work done in dilating the system to the work required for shearing). The contribution of the dilation to apparent friction arises since dilation is required to allow grain rearrangement during shear, yet dilation requires input of work against the normal force. We revisit this connection using theoretical treatment of two-dimensional sheared uniform granular layers and complementary Discrete Element simulations, both for gouge layers and rough surfaces without gouge. The theoretical calculation shows that fluctuations in both apparent friction and dilation rate, that occur during a particular type of grain-scale shear motion, follow a relationship that is non-linear, although in practice appears close to linear. Results show that dilation, and hence apparent friction, are connected to shear localization. In numerical simulations of granular layers (without grain breaking or chemical processes) shear localization occurs but does not persist; instead the systems fluctuate between a state of distributed shear and dilation and localized motion on short-lived shear planes, with overall compaction. The transitions between these two types of shear involve changes in the dilation rate-friction relationship. Large systems of uniform and nonuniform grains show significant scatter about the linear relationship, which are connected to variations in the extent of shear localization.
Manga, M. & Dufek, J. In-situ production of ash and clast breakup in pyroclastic density currents: model results and a preliminary field validation Session 1 Thursday 10  poster file 
Abstract: Abrasion and comminution of clasts during the propagation of pyroclastic density currents (PDCs) have long been recognized as a potential source for the enhanced production of volcanic ash. However, the importance of comminution on both ash production and flow characteristics has not been quantified. The amount of ash produced during transport can potentially affect runout distance, deposit sorting, the volume of ash introduced in the upper atmosphere, and internal pore pressure. We develop a numerical multiphase flow model to evaluate the amount of ash produced by particle-particle interactions during transport. The subgrid-scale physics in this model is based on laboratory experiments of particle collisions. We find that for typical pyroclastic density currents, 10-20% of the initial clasts comminute into ash with the percentage increasing as a function of initial flow energy. Most of the ash is produced in the high-energy regions near the vent, although flow acceleration on steep slopes can produce ash far from the vent. On level terrain, this ash generates gravity currents that detach from the main flow and can more than double the effective runout distance. Currents that descend steep slopes produce the majority of their ash in the flow head, and flow snouts likely develop sub-angular to rounded pumice during this process. To test our model predictions we examined the roundness of pumice clasts deposited in pyroclastic density currents at Mount St Helens, Washington and Mount Lassen, California. Because clast comminution should progressively round clasts during transport, we expect a relationship between mass loss and roundness. We first determined experimentally a relationship between mass loss and particle roundness for pumice from both volcanoes. With this empirical calibration, we can relate roundness of clasts in actual deposits to a mass loss and hence to the amount of ash produced. Given the uncertainties in this calibration and the parameters that are used in the numerical simulations, we find good agreement between expected and measured roundness at both volcanoes.
Marzocchi, W. Some physical insights on earthquake statistics derived from a physics-based earthquake seismicity simulator Session 4 Tuesday 8  poster file 
Abstract: Earthquake occurrence stems from a complex interaction of processes that are still partially unknown. This lack of knowledge is revealed by the different statistical distributions (sometimes antithetical) that have been so far proposed, and by the different beliefs about the role of some key components as the tectonic setting, fault recurrence, seismic clusters, and fault interaction. In this presentation, we explore these issues through a numerical model based on a realistic interacting fault system. The model has three main features: 1) it may account for a realistic 3-D fault distribution; 2) it imposes, a priori, a simple behavior of earthquake occurrences on each fault (purely recurrent or purely Poisson); 3) it allows faults to interact through a co- and post-seismic Coulomb Failure Function (CFF) model. The model has been applied (in collaboration with Italian and US colleagues) to areas with different tectonic style, Central Italy and California. Here, we report the results of these applications in order to understand which physical components are more important to describe the statistics of earthquake occurrences. Specifically, we focus the attention on the time domain, comparing the time features of real seismic catalogs with the synthetic catalogs obtained by the model. The results highlight many interesting issues, such as the role of recurrent fault behavior versus earthquake time clusters, insights on long-term modulation of earthquake occurrence rates, and possible intrinsic limitations of earthquake forecasting based on CFF.
Marzocchi, W. Forecasting large earthquakes and eruptions: is it a scientific issue? Special Session Thursday 10  oral file 
Abstract: The forecast of (large) events is of primary importance for two very different purposes. The first one is practical since reliable and skillful forecasts are the basic components for establishing sound risk mitigation actions. The second one is more philosophical, and it is related to the fact that forecasting is a cornerstone of scientific knowledge. At this purpose, the World's Largest General Scientific Society AAAS reports ``the growing ability of scientists to make accurate predictions about natural phenomena provides convincing evidence that we are really gaining in our understanding of how world works'' (AAAS 1989, 26).
The basic problem in forecasting turns to be if we are able to check reliability and skill of forecasts. In general, this is not a problem if large amount of data are available. For instance, in seismology and volcanology we can usually check the reliability and skill of forecasting models for small to moderate events. The largest earthquakes and eruptions - that are obviously the most important to forecast - are rare and it is impossible to get a reasonable amount of events in a reasonable period of time to conduct satisfactory testing experiments. Seismology and volcanology face the lack of large events for testing purposes in two different ways: in one case, it is assumed that the largest events are a sample of the distribution of the small-to-moderate events (no ``black swan'' exists). This means that the statistical distribution of small-to-moderate events is extrapolated for the largest ones. The second strategy is to assume that largest events have some peculiarities that make them distinguishable from smaller events (black swans do exist!). This means that the above-mentioned extrapolation is no longer valid and a specific statistical distribution is needed. The problem is that we do not have enough data to build the model and, even worse, we do not have enough data to check their reliability and skill. This is the domain where the expert opinion plays a basic role. Scientists tend to abandon the scientific domain based on hypothesis and testing, in order to build consensus models based on expert opinion or, equivalently, on the so-called ``best available science''. In this presentation I will discuss these issues, showing recent examples derived both from volcanology and seismology, and trying to give some insights on the basic question ``Can our models only predict the irrelevant?''
Mazurova, E.M. Computation of Height Anomaly through the Fast Wavelet Transform
Session 9 Tuesday 8  poster file 
Abstract: Now, the central place in geodesy has been taken by satellite methods of co-ordinate determination using the signals of the GPS and GLONASS satellite navigation systems. Geodetic heights, however, computed from satellite measurements and normal heights determined by geometric leveling, exist independently from each other. The relationship between geodetic heights in a co-ordinate space and normal heights is executed by height anomalies which should be known with high accuracy. The problem of accurate determination of the gravity field transformants is caused by the fact that the classical methods of solving problems of physical geodesy are based on deriving and resolving some integrated equation. Continuous errorless values of gravity anomalies across the whole surface of the Earth are needed for that. But it is practically impossible to make a continuous gravimetric leveling, and the classical approach is not capable to produce an accurate answer without gravimetric leveling being completed. In practice, the initial information is discrete, with a lot of measuring errors, and known not across the whole surface of the Earth. Notice that Neumann's integral from which it is possible to calculate height anomaly is a convolution integral. Therefore, it is possible to use linear transformations for its calculation, such as the Fast Fourier Transform and wavelet-transform. In computing height anomaly, M.S.Molodensky's combined method is applied; with it, numerical integration is conducted only within a certain ``limited area''. The influence of the distant field zone is taken into account by expending the gravity anomaly into series by spherical functions. The present paper discusses the algorithm for calculating the height anomaly in the ``limited area'' in flat approximation. Neumann's integral calculations have been done on the basis of the fast discrete Wavelet-Transform. The results of the research into compressing the Neumann integral kernel that can be performed in calculating this integral by the Wavelet-Transform method are presented. An efficiency comparison of applying the Fast Fourier Transform and Fast Wavelet-Transform to computing the Earth's gravity field transformants has been made.
Mazzieri, I., Smerzini, C., Stupazzini, M. & Rapetti, F. Mortar Spectral Element Method applied to the study of seismic response of 2D alluvial deposits Session 9 Monday 7  oral file 
Abstract: Mortar element techniques, are discretization methods for partial differential equations, which use separate finite or spectral element discretization on nonoverlapping subdomains. In general the different meshes on the subdomains do not match on the interface, so the exact continuity condition at the skeleton of the decomposition is replaced with a weak one, that can be written in terms of a Lagrange multiplier and the jump of the traces. Then suitable quadrature formulas are used to evaluate the associated integrals, involving discrete functions defined on different non-matching grids or of different polynomial degree. Such strategies improve the geometrical flexibility of the spectral discretization, allow for mesh adaptivity, and help to circumvent the loss of accuracy near singularities. In the present work, we consider the Mortar Spectral Element Method coupled to second order time integration scheme to simulate the propagation of elastic waves in complex unbounded media. We consider geometrically non-conforming domain partition where meshes in to the subdomains are independently generated from the neighboring ones and associated with different spectral approximation degrees. We thus subdivide the domain into $n$ non-overlapping subregions $k$, and in each $k$ the grid is assumed to be conforming. The resulting skeleton, i.e., the intersection of subdomain boundaries, is then decomposed into mortar elements that are $(d-1)$-dimensional geometrical entities for a $d$-dimensional problem. The coupling between possible discontinuities at the interfaces of the subdomains, either in h (mesh size) or in N (spectral degree), is handled by mortars. In this contribution we present a first set of benchmark assessing the convergence, accuracy and flexibility of the previous method, implemented in the well known spectral elements based code GeoELSE. Furthermore we illustrate how this recent enhancement of GeoELSE can be effectively used also for the numerical analysis of soft soil alluvial deposit lying on the top of stiff sediments or rock. An applicative example will be provided focusing on the 2D seismic response of the Gubbio alluvial basin (central Italy). This numerical simulation includes the combined effect of strong lateral variations of soil properties and topographic amplification. Some hints on the work in progress to effectively handle 3D problems with MEM are also given.
Melnik, O.E., Gorokhova, N.V., Plechov, P.Y., Blundy, J., Rust, A. & Muir, D. Crystal growth in ascending magmas: from individual crystals to crystal populations Session 1 Friday 11  oral file 
Abstract: Ascent of magmas during extrusive volcanic eruptions involves large changes in pressure-temperature conditions which result in non-equilibrium degassing-induced crystallization. Erupted rocks record magma ascent conditions through zonation in individual crystals and through size distributions in populations of crystals. We have developed a mathematical model of diffusive growth of an individual plagioclase crystal from an undercooled magma. The model considers diffusion of 3 components: anorthite (An), albite (Ab) and the rest of the melt. Diffusive fluxes contain cross-terms that account for the relative motion of the components. On the boundary of the growing crystal, mass conservation of An and Ab is satisfied. The growth rate of a crystal strongly depends on undercooling. Liquidus temperature depends on the An content and pressure. Crystal composition is a non-linear function of the An content of the melt. These features make the multi-component model of crystal growth highly non-linear. Numerical simulation of crystal growth in response to a linear drop in temperature with time shows complicated non-monotonic patterns. Crystal growth leads to depletion of a boundary layer in plagioclase components, meanwhile diffusion leads to re-equilibration of the melt with growing crystal. As a result, the crystal grows in pulses with sharp changes in composition due to changes in the melt composition in the boundary layer. As magma ascends, the competition between growth of preexisting crystals and nucleation of the new crystals determines the crystal size distribution (CSD). If nucleation dominates, the CSD will be shifted to smaller crystal sizes. In steady state, CSD evolution is described by a linear hyperbolic equation with coefficients that depend on ascent velocity, crystal growth and nucleation rates. If these coefficients are known as a function of temperature, pressure and crystal content, it is possible to reconstruct variations of these parameters along the conduit based on measures of the CSD of eruptive products. Obtained values of parameters from the CSD of volcanic dome samples are in good agreement with predictions of mechanical models of conduit flows.
Merino, E. Self-Accelerating Dolomite-For-Calcite Replacement: Dynamics of Burial Dolomitization and Associated Ore Mineralization Session 2 Monday 7  oral file 
Abstract: Attempts at understanding burial dolomitization became a comedy of errors after Weyl (1959, 1960), the first geochemical modeler, recommended ignoring petrography and assumed that replacement happens by dissolution-precipitation. This has led to models that are ad-hoc, contain wrong reactions, ignore feedbacks and petrographic insights, and avoid checking predictions against reality. Dolomitization is a beautiful case of geodynamics. The new model correctly predicts a suite of geological, textural, microstructural, rheological, paragenetic, geochemical and morphological properties of dolomites. The model assumes that the dolomitizing brine must be Mg-rich and - counterintuitively - dolomite-undersaturated. Deep brines from the Mississippi Salt Dome basin do satisfy that requirement, and happen also to be rich in Pb, Zn, Sr, Ba, and SO4, hinting at why dolostones host Mississippi-Valley-type ores. The fact that the dolomite-for-calcite replacement is self-accelerating via Ca2+ has astonishing consequences. (1) It makes dolomitization self-organized in time (in the form of many growth `pulses' per `slice') and space (in sequential limestone `slices', each several meters thick, with dissolution porosity accumulated at the entry into each slice, thus spatially periodic - as observed). (2) Because dolostones are strain-rate-softening and replacement happens by induced-stress-driven pressure solution, the self-accelerating replacive growth must pass continuously to displacive dolomite growth - before shutting itself down abruptly. (3) The same Ca2+ involved in the self-accelerating replacement also drives increasing insertion of submicroscopic calcite slivers into the displacive dolomite crystals, deforming them. The displacive growth causes bilateral veins to form normal to the crystal growth direction. The occurrence in most dolostones of zebra veins consisting of displacive, saddle dolomite crystals and a replacive/displacive contact which is seamless stunningly confirms the new dynamics. (4) The high aqueous Ca2+ left when dolomite growth shuts down drives calcite growth and dedolomitization, and growth of other Ca-bearing minerals such as anhydrite and fluorite, and of sphalerite and galena - as observed. (5) Dolomite growing in late pulses of one slice should replace ores formed in earlier pulses, concentrating the ores downflow - as observed.
de' Michieli Vitturi, M., Clarke, A.B., Neri, A. & Voight, B. Modelling extrusion cycles of dome-forming eruptions Session 1 Thursday 10  poster file 
Abstract: We investigated the transient dynamics of magma ascent along a dome-forming conduit coupled with the extrusion of a solid plug. The numerical model DOMEFLOW, developed by the authors for this work, has been extended in order to simulate stick-slip extrusion cycles. DOMEFLOW is a transient 1.5D isothermal two-phase flow model of magma ascent through an axisymmetric conduit of variable radius, which accounts for gas exsolution, bubble growth, crystallization induced by degassing, permeable gas loss through overlying magma and through conduit walls, as well as viscosity changes due to crystallization and degassing. During ascent, magma pressure decreases and water vapor degasses from the melt as the melt simultaneously crystallizes, causing changes in mixture density and viscosity, and eventually the formation of a solid plug sealing the conduit. The modified model describes the displacement of the vent plug by considering frictional forces between the walls and the plug, producing as a result stick-slip extrusion cycles. In particular, when the magma pressure at the base of the plug exceeds the sum of the magmastatic pressure and the static friction, the plug begins to slip along the conduit walls and the pressure drops, eventually leading to a new stick-phase and a new cycle. Results are compared to well-documented cyclic phases of the ongoing eruption of the Soufrière Hills volcano, Montserrat, in order to demonstrate the appropriateness of the formulation. The model is then used to understand basic controls on cycle period and its sensitivity to magma ascent rate for a given set of magma characteristics and imposed friction coefficients.
Morra, G. Particle simulations with Stokeslet and Stresslet kernels Session 9 Monday 7  oral file 
Abstract: Smoothed particle methods in applied mechanics and in hydrodynamics (SPAM and SPH) are commonly employed in fluid and solid mechanics. In this work I explore how a similar formulation but with a different particles kernel enables solving problems such as RT-instability and convection in non-homogenous media. The innovation of this work is based on the observation that BEM is in fact a particle approach, in which the particle kernels are described by the anisotropic field due to the integral of a Green function on a panel. The new kernels are based on Stokeslets and Stresslets, and their integrals on panels (common name for Boundary Integral Elements in BIEM simulations). Classically SPAM and SPH cannot be applied to problems that treat sharp surfaces or interfaces and to materials under tension. I introduce Stokeslets and Stresslets as efficient kernels for pure ``particle'' or ``particle in panels'' simulations to show how they can simulate such sharp interfaces. This approach presents fundamental advantages compared to SPAM and SPH, because of the combination of different kernels in one simulation, which allows solving non-homogenous Stokes flow. I show how this approach allows advances in modelling sharp surfaces and interfaces, and suggest how it can be applied to fundamental problems in geophysical fluid-dynamics.
Morra, G., Seton, M. & Muller, D. Hierarchical self-organisation of tectonic plates
Session 5 Thursday 10  poster file 
Abstract: It is well known that the earth's surface is divided into plates of different sizes. It has been recently pointed out that the plate sizes display a hierarchical organisation. Using recently released reconstructed plate boundary polygons for the past 140 Myrs, we analyse the behaviour of their hierarchical structure. We find that (i) the distribution of the largest and the smallest plates are always decoupled in the last 45Myrs, confirming a match of the critical exponent of the small plates with a fragmentation mechanism; (ii) the distribution of the largest plates during the last 140Myrs is also a power law, involving the 6-7 largest plates; (iii) the fluctuations of the power law exponent for the largest plate sizes varies in a timeframe of tens of millions of years, reaching a maximum of almost one at 60-50 Ma, and a minimum of almost zero at 110-100 Ma; (iv) in the period 100-80 Ma the growth is the fastest ever observed; (v) in the last 100 Ma the fractal exponent appears to behave as a fast excitation followed by a slow relaxation. The coincidence between the initiation of spreading at mid-ocean ridges and jumps in ridge position during excitation, and plate rejoining or accretion during relaxation suggests that the plate tectonic system is a global self-organising system. The coupling between plate forces and mantle convection is analysed in detail to find the dominant physical process.
Morra, G. & Yuen, D.A. Pulsations of a plume rising through a non-monothonic mantle
Session 9 Tuesday 8  poster file 
Abstract: Recent observations in seismic tomography and mineral physics of the spin transition at mid mantle depth suggest a non-monothonic viscosity radial patterns in the lower mantle. We investigate its effect to the upwelling of a mantle plume, explicitly modelling its rise. Models are performed to study the short or long term effect of the non-monotonicity to the fluxes in the plume tail. The numerical method is based on a Multipole-accelerated Boundary Element (FMM-BEM) implementation, for solving the momentum equation in spherical coordinates. We test several radial viscosity patterns, starting from the well-known hypothesis of a viscosity hill in the middle of the lower mantle and expanding it to more complex scenarios. In particular we test the combined contribution of viscous and density radial profiles. The constraints on the profile come only from first principle and laboratory mineral physics data. We test both the effect of the spin transition and/or of the chemical layering. We identify critical viscosity ratios and differential densities that enable the system to pulse, inducing a punctual manifestation at the earth surface. Such sharp transition between linear and pulsating plume pattern is relevant for identifying the origin of the punctuation of the islands in long oceanic hot-spot seamounts.
Naylor, M., Mudd, S.M. & Yoo, K. Markov Chain Monte Carlo (MCMC) Inversion of Hillslope Elevation and Soil Thickness Data for the Baselevel History Session 7 Thursday 10  poster file 
Abstract: Hillslope topography and soil thickness respond to changes in river incision or deposition. For example, accelerated river incision leads to a wave of steepening and soil thinning that begins at the channel and moves upslope [1]. Because of the coupled response of topography, soil thickness and channel incision or deposition rates, it may be possible to use hillslope properties to reconstruct the erosional or depositional history of channels. A prerequisite for such inversion of hillslope properties to reconstruct historical landscape dynamics is a method that allows one to quantify both the most likely channel history as well as the uncertainties in changing channel erosion or deposition rates through time. Here we present robust methods ideally suited for this purpose: Monte Carlo Markov Chain (MCMC) methods. Specifically, MCMC methods [2] involve (i) taking some assumed base level history, (ii) perturbing that history (iii) running a forward model to estimate new hillslope profiles (iv) choosing whether to accept the new history by a Metropolis-Hastings Algorithm (v) storing the favoured history and repeating. In this way we iterate towards the most likely channel history whilst exploring parameter space in such a way that confidence intervals can quantified. Here we demonstrate how this approach returns not only a best estimate history, but also credibility intervals which reflect the progressive loss of information with time. These techniques are generic and should be employed more generally within geomorphology.
[1] Mudd, S.M., and D.J. Furbish (2007), Responses of soil mantled hillslopes to transient channel incision rates, Journal of Geophysical Research-Earth Surface, 112, F03S18
[2] K. Gallagher, K. Charvin, S. Nielsen, M. Sambridge and J. Stephenson (2009) Markov chain Monte Carlo (MCMC) sampling methods to determine optimal models, model resolution and model choice for Earth Science problems Marine and Petroleum Geology 26(4)
Naylor, M., Touati, S., John, G., Main, I.G. & Bell, A. Interpretation of statistical signals in earthquake data Session 4 Tuesday 8  poster file 
Abstract: Linking physical process with empirical statistical trends associated with extreme events is challenging. Non-uniqueness in models, biasing in data selection/processing and the apparent emergence of structure in complex stochastic systems all act to limit our ability to agree on a common null hypothesis, let alone quantify inherent predictability within complex earth systems. For example, here we use earthquake statistics to explore evidence for Characteristic earthquakes [1] and universality in earthquake interevent time distributions [2,3] to demonstrate how the statistics of extreme distributions [4] and selection bias can give rise to emergent trends that have previously been used to support physical models rather than being identified as statistical artefacts.

[1] Naylor, M., J. Greenhough, J. McCloskey, A.F. Bell, and I.G. Main (2009), Statistical evaluation of characteristic earthquakes in the frequency-magnitude distributions of Sumatra and other subduction zone regions, Geophys. Res. Lett., 36
[2] Touati, S., Naylor, M. & Main, I.G., (2009) Origin and nonuniversality of the earthquake interevent time distribution, Phys. Rev. Lett .102, 168501
[3] S. Touati, M. Naylor, I.G. Main and M. Christie (Submitted) Masking of earthquake triggering behaviour by a high spontaneous rate and implications for ETAS inversions
[4] Naylor, M., Main, I.G. & Touati, S. (2009) Quantifying uncertainty in mean earthquake interevent times for a finite sample, J. Geophys. Res, 114, B01316.

Neri, A., Andronico, D., Aspinall, W.P., Baxter, P.J., Bertagnini, A., Cioni, R., Esposti Ongaro, T., Flandoli, F., Giorgi, E., Macedonio, G., Mulas, M., Papale, P., Pistolesi, M., Rosi, M. & Todesco, M. Assessing pyroclastic density current hazard from sub-Plinian eruptions at Vesuvius (Italy) and associated uncertainties Session 7 Thursday 10  poster file 
Abstract: Pyroclastic density currents (PDC) represent a potential enormous risk at Vesuvius. The chances of survival in locations affected by these flows are so low that these areas need to be evacuated in advance of any eruption that produces such flows. The aim of this work is to produce probabilistic hazard maps of PDC invasion and runout on the basis of field reconstructions and numerical simulations, supplemented by expert judgement. Fieldwork studies have produced detailed stratigraphic analyses of deposits produced by sub-Plinian events at Vesuvius, such as the AD 1631, 472 and Greenish eruptions, indicating the areas that were invaded and the maximum distances reached by the flows. In complementary work, transient and 3D numerical simulations of column collapse events and associated PDC have been carried out for partial and near-total collapse of a sub-Plinian volcanic column. Results simulate the complex dynamics of the near-vent collapse and indicate the main effects of volcano topography on propagation of flows. These different strands of empirical and modelling evidence were critically analysed, appraised and combined using two expert judgement pooling algorithms: the Cooke ``Classical model", and a newly designed point-wise estimation model named ERF (see Flandoli et al., Proceedings CMG 2010). These were used to account for uncertainties, and quantify probabilistic hazard maps. The assessment included consideration of several sources of uncertainty including vent location, eruption intensity (i.e. mass flow-rate), flow generation mechanism (partial column collapse vs low-fountaining), and the properties of the eruptive mixture and flow. Epistemic limitations of the stratigraphic studies and epistemic and aleatory limitations of numerical modelling were also accounted for in the elicitation process. The findings reveal, for any vent location area, the remarkable spreading of the probability of flow invasion and maximum runout due to the numerous sources of uncertainty investigated. Moreover, the uncertainty in vent location within the caldera strongly determines the relative probability of pyroclastic flow invasion of different sectors around the volcano, as well as corresponding runouts. Although it is likely, and desirable, that future studies will significantly reduce epistemic sources of uncertainty, the present results address the key problem of how to handle uncertainty, and raise the question of how to communicate it to decision makers.
Nermoen, A., Galland, O., Jettestuen, E., Fristad, K., Podladchikov, Y.Y., Svensen, H. & Malthe-Sørenssen, A. Experimental modeling of piercement structures Session 1   withdrawn file 
Abstract: Piercement structures such as mud volcanoes, hydrothermal vent complexes, pockmarks and kimberlite pipes are observed in a range of geological settings. These structures are formed during the release of over-pressurized fluids in the upper part of the crust and consist of circular pipes with highly brecciated rocks. The goal of this work is to predict the conditions at which piercement structures form in general. We have conducted sand box experiments by injecting compressed air through an inlet at the base of a bed of glass beads. We varied systematically the depth $h$ and width $w$ of the inlet. During the experiments, we manually increased the inlet gas velocity until fluidization occurred. The fluidized zone consisted of a diverging cone-like structure where the grains are ejected up through the centre and flow down along the margin, producing inward dipping strata similar to those observed in nature. In each experiment, we measured the critical velocity $v_f$ at which fluidization initiates. By using dimensional analysis, we show that $v_f$ normalized by the intrinsic Darcy velocity is correlated to the ratio $h/w$. We derived an analytical model to estimate the critical velocity $v_f$. The model consists of a force balance between the buoyant weight of the fluidized zone and the seepage forces due to the fluid flow through the grains. The total seepage force is obtained by integrating the pressure gradient, that was calculated analytically by solving the Laplace equation with similar boundary conditions to those of the experiment. The analytic model reproduces the observed correlation between $v_f$ and $h/w$, but it slightly under-estimates the fluidization velocity. This result suggests that the seepage force is the main force triggering fluidization. However secondary forces, such as internal friction, may play a minor role. From the solution of the Laplace equation combined with the experimental results, we defined the Venting number, which is a criterion for the onset of fluidization. From the Venting number we estimated the critical pressure of fluidization in a variety of geological environments. Notably so, a comparison of the analytical pressure estimate and measurements of the fluid-pressure profile for the Lusi mud volcano (Indonesia) shows a good consistency. This suggests that the derived model for the Venting number is a useful tool to predict fluidization in a variety of geological settings.
Olhede, S.C. & Simons, F.J. How to measure the strength of the lithosphere without using the admittance or coherence between gravity and topography Session 5 Thursday 10  poster file 
Abstract: The lithosphere is modeled using a differential equation characterized by a set of parameters, at least one of which, under the assumption of elastic behavior, is generally thought of as a proxy for its strength: the flexural rigidity (D), or, by extension, the elastic thickness. This lithospheric system then takes an input: topographic loading by mountain building and other processes, and maps it into an output: the gravity anomaly and the final, measurable, topography. The input is not measurable but some of its properties can be characterized. The outputs are measurable but the relation between them is obfuscated by their stochastic nature and the presence of unmodeled components. Estimating D, most usually in the spectral domain, generally involves constructing summaries of gravity and topography. Both admittance and coherence are popular; both are ratios of the cross-spectral density of gravity and topography to the power spectral densities of either, the whole sometimes squared. Despite the fact that neither admittance nor coherence are Gaussian, estimating D usually comes down to the least-squares fitting of a parameterized curve, where Gaussian behavior is tacitly assumed. In this two-step procedure, admittance or coherence are first estimated, and subsequently inverted for the strength parameters. Rarely, if ever, are lithospheric models found that satisfy both coherence and admittance to within their true error. Why don't they? Poorly characterized errors of admittance and coherence are not the only problems with this procedure. There is also the implicit annihilation of information during the construction of these statistics (coarsely sampled, sometimes squared, ratios, measures of the data as they are) themselves. Then there is the fact that we do not want to know coherence and admittance at all - we want to know properties of the lithosphere! In this presentation, we intend to abandon coherence and admittance studies for good, by proposing an entirely different method of estimating flexural rigidity, which returns it and its confidence interval, as well as a host of tests for the suitability of the assumptions made along the way, and the possible presence of correlated loads and anisotropy in the response. The crux of the method is that it employs a maximum-likelihood formulation that remains very grounded in the data themselves, and which is formulated in terms of variables that do have a Gaussian distribution.
Padilla, L. & Vallis, G.K. Joint effects of natural variability and forcing uncertainty on observational estimates of climate sensitivity Session 6 Tuesday 8  poster file 
Abstract: We explore the joint effects of natural variability and forcing uncertainty on the transient climate sensitivity (i.e, the short-term response of the surface temperature to a doubling of CO2) of the Earth. We assimilate observations of temperature in the 20th century and various estimates of the radiative forcing into a simple energy balance model using a ``sigma point'' Kalman filter. Use of the filter, along with estimates of natural variability and forcing uncertainty, allows us to estimate both the climate sensitivity and its uncertainty. Currently, the forcing uncertainty provides the largest source of errors in our estimates. If our knowledge of the radiative forcing were to improve we could expect to provide better estimates of climate sensitivity sometime in the future, then constrained by natural variability. We quantify this statement, providing an estimate of when in the future we will have a good estimate of climate sensitivity for given levels of forcing uncertainty and natural variability.
Pais, M.A. & Schaeffer, N. Extending the Quasi-Geostrophic flow model for core flow inversions. Session 5 Thursday 10  poster file 
Abstract: Tangentially Geostrophic (TG) or Quasi-Geostrophic (QG) flow model do not allow the flow at the core surface to cross the equator. In addition, neither model is valid near the equator. The goal of this study is to try to overcome these limitations.
Inversion of the geomagnetic secular variation for core surface flow models with a gradually softer TG constraint, makes it possible to identify regions on the core surface where deviations to the TG condition first occur. These are localized near the equator as expected from the anhilation of the Coriolis forces there and the breakup of TG or QG force balances that may prevail in most of the core surface or core volume regions, respectively. Deviations to TG are mostly equatorially symmetric, suggesting an anti-symmetric flow component emerging in equatorial regions.
We perform a numerical experiment whereby an antisymmetric forcing, localized at the top and bottom of a spherical core, induces a flow with an axial component that changes sign at the equator and is close to a linear dependence on z.
Based on these observations, we extend the QG model with a new velocity component that depends linearly on z. We deduce the kinematic conditions to be obeyed by this flow at the core surface, and use them to invert for an extended QG-flow model. The regularizing quadratic forms we use are also inspired by the numerical simulations that show formation of zonal jets with steep latitudinal gradients.
The effect of the Inner-Core on the resulting flow does clearly show up, without having been imposed. Furthermore, the flow is mostly symmetric in a region spanning from the cylinder tangent to the inner core to latitude 20 to $30^rc$ from the equator. Length of Day predictions seem also to be improved.
Panet, I., Holschneider, M., Diament, M., Kuroishi, Y. & Mikhailov, V. Multi-scale modeling and analysis of the Earth's gravity field using Poisson wavelets Session 5 Friday 11  oral file 
Abstract: The Earth gravity field reflects the mass distribution and mass displacements inside the Earth system, from the surface and fluid envelops to the core. With the advent of satellite gravity and the launch of the CHAMP, GRACE and GOCE missions, measurements of unprecedented quality are available and provide a new view of the static and time-varying gravity field at low and medium resolutions. Combined with the surface gravity data within high resolution models, they allow to study geodynamic processes and their interactions in a wide range of spatial and temporal scales.
To derive such models and analyze them, Poisson multipole wavelets (Holschneider et al., 2003) appear particularly promising. Because of their localization properties both in space and frequency, they are well-suited to derive regional gravity models from heterogeneous datasets. The analysis of the obtained models allows to unfold the different components of the gravity field in different scales and locations and thus, to highlight specific geodynamic signals of interest.
We present a few examples of applications of this approach to model and analyze the static and time-varying gravity field, for studying the Earth's structure and the mass redistributions after large earthquakes.
Parker, W.S. Confirmation and Testing of Scientific Models, Revisited Session 7 Tuesday 8  oral file 
Abstract: In a well-known paper in Science magazine, Oreskes et al. (1994) counseled against the oft-used terminology of model ``verification" and ``validation". The concept of confirmation, however, remains on the table both in their discussion and in many subsequent treatments of model evaluation.

I argue that usually what we can sensibly aim to confirm are not scientific models themselves, but rather their adequacy for particular purposes. I illustrate how testing for adequacy differs from testing for the truth of modeling assumptions. Explicitly tying model evaluation to purposes in this way, I suggest, has the advantage of discouraging the mistaken assumption that successes had by a model in one context constitute good evidence of its general adequacy.
Nevertheless, a concept other than confirmation might be even more useful when it comes to the practice of model evaluation. On a confirmation-driven approach, attention is not necessarily focused on testing models in the most informative ways; we may favor tests of convenience, even if they have little chance of revealing that our model is inadequate for the purpose of interest (if in fact it is inadequate). An alternative is to approach model evaluation as an activity aimed at severely testing, rather than confirming, the hypothesis that the model is adequate for the purposes at hand. A ``severe test" of some hypothesis H is a procedure that is likely to indicate that H is false, if and only if H is in fact false (Mayo 1996).
After some remarks on how severe testing for adequacy-for-purpose might proceed, I argue that even if severe testing is often not achievable, the concept of severe testing would be valuable as a regulative ideal. A focus on severe testing might encourage better use of the limited resources available for model evaluation and, in the face of predictive successes, would direct attention to important questions about the testing procedure that are often overlooked in model evaluation at present. I briefly illustrate how this might play out in the context of climate model evaluation.

Pasquero, C. & Tziperman, E. Statistical parameterization of heterogeneous oceanic convection Session 3 Tuesday 8  poster file 
Abstract: A statistical convective adjustment scheme is proposed that attempts to account for the effects of mesoscale and sub-mesoscale variability of temperature and salinity typically observed in the oceanic convective regions. Temperature and salinity in each model grid-box are defined in terms of their mean, variance and mutual correlations. Sub-grid scale instabilities lead to partial mixing between different layers in the water column. This allows for a smooth transition between the only two states (convection on and convection off) allowed in standard convective adjustment schemes. The advantage of the statistical parameterization is that possible instabilities associated to the sharp transition between the two states, that are known to occasionally affect the large scale model solution, are eliminated. The procedure also predicts the generation of correlations between temperature and salinity and the presence of convectively induced up-gradient fluxes that have been obtained in numerical simulations of heterogeneous convection and that cannot be represented by standard convective adjustment schemes.
Patra, A.K., Dalbey, K., Stefanescu, E.R., Pitman, E.B., Calder, E., Bursik, M.I. & Sheridan, M.F. Using Computer Models and Uncertainty Quantification to Construct Hazard Maps Session 7 Thursday 10  oral file 
Abstract: The ultimate benefit of computer modeling of hazards arises in the ability to forecast the effect of a catastrophic event and to prepare maps of the probability of such hazard incorporating these predictions. The simplest way to use such a model to forecast hazard would be to use a Monte-Carlo like sampling strategy over the range of possible inputs and then compute the probability of different outcomes. However, in a setting where simulations can cost hours on modern supercomputers the use of such methodology is prohibitively expensive though occasionally still attempted. We present here an alternate approach using a systematically chosen set of simulations to construct a statistical model and then sampling to construct the hazard map. To overcome the cost of inverting the extremely large covariance matrices to assimilate data from new simulations we use a simple localization procedure which also enables the use of parallel computing. The computational time required is thus reduced from one that was estimated to take O(100) days to 1 day on 1024 processors.
Peltier, W.R. & Griffiths, S.D. Internal waves in the atmosphere and oceans and explicit models of the internal tide Session 3 Tuesday 8  oral file 
Abstract: Topographically forced internal waves are ubiquitous features of both the atmospheric and oceanic general circulation. In the atmosphere the breaking of such waves above high elevation topography provides the explanation of the intense downslope windstorms that occur frequently in the lee of topographic extrema. The best known examples of such events are those which occur frequently over the Rocky mountains to the west of Boulder Colorado. In the coastal oceans the same phenomenon is also a recurrent feature in circumstances in which tidally forced flow over a sill, say at the mouth of an estuary, also generates intense internal wave activity and wave breaking. The best known example in the coastal ocean is from Knight Inlet on the coast of the Canadian province of British Columbia. There, as a consequence of the existence of the overlying free surface, the phenomenon is often accompanied by the generation of intense upstream propagating solitary waves. On a global scale in the oceans the same process as operates in the coastal environment is responsible for the the generation of the so called ``internal tide'' that arises due to the interaction of the barotropic tide with bottom topography. Since the internal tide is responsible for a substantial fraction of tidal dissipation and is therefor heavily implicated in the evolution of the orbit of the Moon, it is a phenomenon of planetary scale relevance. Detailed analyses have been performed of the impact of the modest changes in the geometry of the ocean basins associated with the Late Quaternary ice-ages on the partition of tidal dissipation between bottom friction and the internal tide. It is demonstrated that intense resonances develop at high latitudes in both hemispheres that are expected to have extremely important climatological consequences. These are primarily associated with the excitation of gravest mode Kelvin waves at the M2 frequency.
Peruzzo, E., Barsanti, M., Flandoli, F. & Papale, P. The stochastic quantization method and its application to the numerical simulation of volcanic conduit dynamics under random conditions Session 7 Thursday 10  poster file 
Abstract: The use of numerical codes to solve theoretical models of volcanic systems allows a thorough description of the physical and chemical processes occurring at volcanoes. Due to intrinsic uncertainties inherent in the definition of specific volcanic systems and according to the definition of volcanic hazard, the probabilities of possible volcanic scenarios should be estimated as a function of probability distributions describing the input quantities. However, since the computational costs of the numerical codes involved are often very high, it is fundamental to reduce the number of simulations required to reasonably cover the spectrum of possible conditions expected. We accomplish this task through implementation of the Stochastic Quantization (SQ) method [1]. This technique allows optimal approximation of a continuous probability distribution with a discrete one. Our strategy can be applied to general systems with uncertain inputs. We show the results of a benchmark test based on a one-dimensional steady model of multiphase magma flow in volcanic conduits: the output distribution estimated by SQ is very close to a Monte Carlo (MC) estimate of the true one based on $10^3$ simulations. In general, the SQ method turns out to provide substantially better results than the MC method with the same number of simulations. Our analysis shows that probabilistic MC scenarios with order $10^3$ simulations are well reproduced by order 10 SQ simulations.
[1] S. Graf, H. Luschgy, Foundations of quantization for probability distributions, Springer, Berlin-Heidelberg, 2000.
Petroff, A.P., Devauchelle, O. & Rothman, D.H. The bifurcation of channel heads cut by springs Session 2 Tuesday 8  poster file 
Abstract: When groundwater emerges from a spring with sufficient intensity to remove sediment, it channelizes the landscape. Proposed examples of such ``seepage channels" include centimeter-scale rills on beaches and levees, hundred-meter-scale valleys on Earth, and kilometer-scale valleys on Mars. Past work has identify the growth and formation of channel heads as responsible for the formation of a network of streams. Here we extend this work by describing a mechanism by which new channel heads form. Based on observations of a kilometer-scale network of seepage channels, we develop a model that describes how the focusing of groundwater into the head causes it either to grow forward or to bifurcate. Rain falling around the network seeps into the ground where variations in the height of the water table cause it to flow into the channels. The resulting subsurface flow is determined by the Poisson equation. Because the size of a stream is very small compared to the spacing between streams, we describe seepage channels as a one-dimensional network growing in a Poisson field. If a random perturbation causes a channel head to bifurcate when it drains water flowing mainly from a single direction, one of the newly bifurcated heads will quickly out-compete the other and a single head will grow forward. In this case, the channel head is stable. Conversely, if a channel head bifurcates when it drains water flowing from all directions, the heads continue to grow in different drainage basins. In this case, the channel head is unstable. We develop a model in which this different behavior is exhibited by a pitchfork bifurcation when the channel reaches a critical length. This bifurcation is a consequence of a difference between the growth of one-dimensional curves in a locally driven, Poisson field versus an infinite Laplace field fed from far away: while a bifurcation in an infinite Laplace field is quickly damped by competition, a Poisson field allows bifurcated heads to remain stable by growing in different directions. This model predicts the critical length of a channel, which we compare to observed bifurcations in the field.
Phipps Morgan, J., Hasenclever, J. & Shi, C. Better Strategies for Finite Element Solutions of Variable Viscosity Stokes Flow Session 9 Monday 7  oral file 
Abstract: The numerical solution of variable viscosity Stokes Flow has received relatively little attention by the modeling community until recently. Here we discuss several known computational hurdles to progress, and suggest strategies that offer promise in overcoming them. We will also discuss possibly even more effective MATLAB-based flow-algorithms that can better utilize GPGPUs within a distributed parallel-computing environment.
The choices for solving the discrete pressure equation arising from Stokes flow typically involve several tradeoffs between speed and storage requirements. Our preferred recipe currently has the following ingredients: (1) Use accurate and efficient LBB-stable elements that lead to good convergence in iterative solution schemes (Silvester and Thatcher, 1986; Pelletier et al., 1989). We currently recommend Cornell Macroelements, as they are LBB-stable, are as easy to code as Taylor-Hood elements and have the same convergence properties as them, yet are much more locally incompressible and have no known incompressibility-related artifacts. In contrast, Taylor-Hood elements are 'too squishy' and also have a potential numerical artifact related to capturing an accurate hydrostatic pressure solution (Pelletier et al., 1989). (2) Use an inexact CG pressure preconditioner to drastically reduce the necessary tolerance for the inverse-velocity calculation within it. This leads to a work-reduction of $2-3-fold. (3) Use a coarse-grid preconditioner as part of a Dryja-Widlund Additive Schwartz multilevel preconditioner. The key point is to never compute any additional global fine-grid pressure-matrix operations except for the single operation needed to calculate each new CG search direction (i.e. avoid smoothing operations). We are currently exploring an apparently promising $>2$ level cascade-type improvement of this general preconditioning strategy.
To make more effective use of GPGPU-enhanced computing we are exploring the use of tile-based DD-like preconditioners and smoothers, and block-adaptations of element-by-element matrix-vector multiplies, which are straightforward to implement in a MATLAB+Jacket programming environment. As a general observation, we find that the relative ease of list-making, book-keeping, vectorizing, and loop-unrolling in MATLAB makes it an extremely useful environment for development of easy-to-enhance parallel finite-element codes for HPC-based research in geodynamics.
Pilkey, O.H. & Pilkey-Jarvis, L. Predictive Modeling of Processes on the Surface of the Earth Doesn't Work Special Session Thursday 10  oral file 
Abstract: We argue that quantitative modeling of earth surface processes, at accuracies useful for engineering purposes, is not possible. Many such models (such as climate models) involve human behavior, always an impossible parameter. Simplification of parameters is often extreme such as the frequency, nature and magnitude of future storms and floods in sediment transport models. The parameters in quantitative earth surface process models are assumed to be the most important ones but when ``unusual'' events occur other parameters are involved, not accounted for in the model. Frequently quantitative models are verified by hind casting but ordering complexity makes this impossible. Ordering complexity refers to the interactions among the numerous components of a complex system that occur in unpredictable and unexpected order and sequences. Successful reproduction of an event in a complex natural system is no guarantee that the model will successfully predict the next event.
Pitman, E., Patra, A.K., Calder, E., Bayarri, M., Berger, J.O., Dalbey, K., Spiller, E. & Wolpert, R. Quantifying Input Uncertainty in Models of Volcanic Flows Session 7 Thursday 10  oral file 
Abstract: Mathematical simulations of volcanic mass flows supplement field deposit studies, providing data from which to make predictions of volcanic hazards. Several parameters and functions must be specified upon input, in order to start the simulation. These inputs include the total volume of the mass flow, the initial direction of the flow, friction and dissipation factors, and the terrain over which the mass moves. These inputs are often poorly characterized, or are known to be inaccurate. These inaccuracies affect the output results of the simulations. Here we provide a methodology to account for uncertain inputs and the effect of uncertainty on outputs, in order to provide a quantitative measure of hazard at locations downstream from the volcano. This approach draws on insight from the geo-science community, methods of high performance computing, and novel use of Bayesian statistics and mathematical modeling.
Polacci, M. Implications of rock texture characterization on the modelling of volcanic processes Session 1 Thursday 10  poster file 
Abstract: Volcanic processes are complex phenomena that require a multidisciplinary approach to understand their behaviour and to monitor and forecast volcanic activity. It has long been proved that the nature and physical characteristics of volcanic rocks provide fundamental information about the history and evolution of the volcanic system that has generated them. In the following we report examples of how studies conducted on rocks erupted from Italian volcanoes and other volcanic areas in the world have been integrated with numerical models to understand the physics of volcanoes and the physico-chemical properties of their magmas. High-resolution 2D (via scanning electron microscopy) and, more recently, 3D (via X-ray computed microtomography) imaging on volcanic rocks has led us to 1) to investigate the degassing behaviour of basaltic and trachytic magmas and to provide implications for the eruptive dynamics and style of the related volcanic systems. Permeability measurements and lattice Boltzmann simulations of fluid flow have been also used as tools to constrain magma degassing; 2) to find a relationship between results from geophysical measurements and the distribution of gas vesicles in basaltic scoria; 3) to define the rheological properties of magmas of mafic to felsic composition and use the results in the modelling of magma ascent dynamics during Strombolian, fire-fountain and effusive activity as well as during large Plinian eruptions. The results of these works have strongly pointed out that parameters as vesicularity, vesicle number density, vesicle interconnectivity, permeability, crystallinity, crystal shape etc. need to be included in numerical and theoretical models in order to be able to reproduce accurate conduit flow conditions and forecast reliable eruptive scenarios.
Power, W.L. & Tolkova, E. Obtaining natural oscillatory modes of bays and harbors via Empirical Orthogonal Function analysis of tsunami wave fields Session 8   withdrawn file 
Abstract: To a tsunami wave, bays and harbours represent oscillatory systems, whose dynamical properties determine the response to tsunami and consequently the potential hazard. It has been observed that tsunami records at tide gauges are often dominated by frequencies present in the local seiche records, and not the tsunami frequency in the open ocean. Thus the tsunami energy is distributed among local oscillatory modes. In such cases, near-shore observations, in particular gauge records, do not provide the comprehensive picture of the wave behavior in a harbour. For instance, the largest waves or currents may occur in the center of the harbour or at its entrance. To comprehend the associated danger, it is desirable to know the resonance modes of a bay or harbour and the likelihood of their excitation.
The direct way to obtain resonance modes of a water reservoir is solving the linearized shallow-water wave equation for its eigenfunctions, under the homogeneous boundary conditions (zero normal flow on the land-water boundary, and purely normal flow into an open ocean). One of the difficulties of posing such a problem for a harbour is specifying its boundary with the ocean. For harbours with a narrow entrance, it suffices to assume a normal flow across the harbour mouth. However, for too many other bays and harbours with wider entrance, the oscillator boundary with the open ocean is as much an unknown as the oscillatory modes themselves.
Using the Empirical Orthogonal Functions (EOFs) formalism, the resonance properties of a bay can be deduced from ``observations" (modeled tsunami wave fields) instead of being derived from dynamical equations. Analysis of the wave fields in a resonator can also answer a question of how likely, and/or in which combinations, different modes are exited, due to feasible natural events. The weak spot of this approach is that in general, EOFs capture features of variability in a particular wave propagating in a body of water, which may not necessary represent the normal modes of the water body. This study, which encompasses both real-valued and complex-valued EOF analysis, focuses on techniques for obtaining and identifying ``true" oscillatory modes. It has been applied to deducing resonance properties of a number of bays and harbours, and interpreting a tsunami wave evolution in these harbours.
Ptitsyna, N., Levitin, A., Dremukhina, L. & Gromova, L. Modeling external magnetic field dynamics during extreme events Session 5 Thursday 10  poster file 
Abstract: The use of many geomagnetic activity models for examining extreme events can easily lead to incorrect estimates. This is due to the fact that most of these models have been established on the basis of statistical processing of space-age data sets that hardly contained any extreme values of the solar wind velocity and density and interplanetary magnetic field, because the data for extreme magnetic storms is rather scarce. For example, only one super-intense magnetic storm has been recorded (geomagnetic index Dst=-640 nT, March 13, 1989) during the space-age. One of such extreme geomagnetic event is the famous Carrington storm on 2 September 1859. For modeling the magnetic field recorded at the Colaba Observatory during this storm we used a new self-consistent version of a time-dependent magnetospheric paraboloid model PM [Feldstein et al., J. Geophys. Res., 2005, A11, A11214 10.1029/2004JA010584]. The model is valid for storms of any intensity and quantifies the contribution of current on the magnetopause, ring current, and magnetotail current to Dst variation. The result of our modeling suggests rather good agreement with the magnetic storm dynamics, which include a very rapid and big decrease in the horizontal intensity of the geomagnetic field during the main phase, followed by a sharp recovery. The main features of the storm are related to temporal dynamics of the magnetospheric tail current system. It is characterized by the rapid Earthward-directed propagation of the front edge of the plasma sheet in the main phase of the storm to a distance of 2-3 RE from the Earth's center and the subsequent rapid return back to a distance of $sim 7-8$ RE in the magnetospheric tail.
Qamili, E., De Santis, A. & Cianchini, G. Shannon information of the geomagnetic field for the past 7000 years and implications on the present field understanding Session 5 Thursday 10  poster file 
Abstract: The present behaviour of the geomagnetic field as expressed by IGRF (International Geomagnetic Reference Field) model deserves particular attention when compared with that shown over the past few thousands of years by several paleomagnetic/archeomagnetic models, such as CALS3k, CALS3k$$cst, ARCH3k, ARCH3k$$cst, SED3k and CALS7k. In order to exploit the dynamical properties of the geomagnetic field in its whole complexity, we apply some Information Theory concepts to the above models, in particular analysing the behaviour over time of Shannon Information and Kolmogorov Entropy. The results show how the present geomagnetic field is rather distinct, at least for the past 400 years. This is consistent with a significant global critical state started at around 1750, and still present, characterized by significant low geomagnetic dipole energy and Shannon information and high K-entropy. The details of how these characteristics may develop are not clear, since the present state could move toward an excursion or a geomagnetic polarity reversal, although we cannot exclude the possibility that the ``critical'' behaviour will become again more ``normal'' and the present geomagnetic field will arrest its decay.
Reeves, D. & Rothman, D. The role of spatial heterogeneities in the apparent age dependence of dissolution and precipitation rate constants in porous media Session 2 Monday 7  oral file 
Abstract: As fluid flows through porous rock, dissolution and precipitation reactions exchange material between the fluid and the reactive surface of the medium. The effective reaction rate at which this exchange occurs is of fundamental importance to many applications, such as carbon sequestration, but remains poorly understood. Laboratory and field measurements of rate constants for common materials yield quantitatively different values. Furthermore, rate constants have been observed to decrease with time according to an inverse power law. We propose that transport and material properties at the pore network scale, such as spatial heterogeneities in the pore size distribution, are critical to the quantitative understanding of reactivity and aging. Using a combination of length and time scale analysis and pore network simulations, we examine situations both near and far from equilibrium. In the latter case, we build on the reactive infiltration literature to show that rapid reactions cause changes in the network properties and flow patterns, thereby affecting the reaction rate itself and introducing feedback. Near equilibrium, fluctuations may lead to a slow evolution of the pore size distribution and spatial heterogeneities. We find that spatial heterogeneities play a critical role in controlling transport of solute through the medium, and can modify the measured reactivity even at apparent concentrations within the so-called dissolution plateau.
Rothman, D.H. Singular Blow-up in the End-Permian Carbon Cycle Session 6 Tuesday 8  oral file 
Abstract: About 252 million years ago, as the Permian period gave way to the Triassic, roughly 90% of all living species disappeared from the fossil record. This event, the most severe extinction in Earth history, was accompanied by a rapid ($sim 10,000$ year) change in the carbon isotopic composition of seawater. By transforming these chemical changes to physical fluxes, we show that the isotopic event is consistent with an incipient singularity in the growth of the oceans' reservoir of dissolved inorganic carbon. The singular influx of CO2 indicates a fundamental nonlinearity in the carbon cycle that, in principle, allows prediction of the extinction event about 100 Kyr in advance. Its identification also suggests that any hypothesis for the extinction's cause should predict such a blow-up. One popular hypothesis suggests that the extinction is caused by Siberian trap volcanism. If Siberian volcanism is indeed the cause, then our findings suggest that it must have led to the progressively rapid combustion, over a period of 100 Kyr, of at least 10 Gt of isotopically light organic carbon (e.g., coal). As an alternative, we identify a biological mechanism that naturally leads to a singular flux of CO2 and discuss its relevance to observations.
Rückelt, J. & Slawig, T. Sensitivity analysis and parameter estimation for a marine biogeochemical model. Session 8 Thursday 10  poster file 
Abstract: Results of parameter optimization and uncertainty analysis for a one dimensional marine biogeochemical model of NPZD type are presented. The model, developed by Schartau and Oschlies, simulates the distribution of nitrogen, phytoplankton, zooplankton and detritus in a water column and is driven by an ocean circulation model. It is tested against BATS (Bermuda Atlantic Time-Series Study) data. We study the model numerically, using sqp solvers with exact gradients provided by automatic differentiation. Contrary to published results, we find only little sensitivity of the optimal parameters to observational error alone. We compare the results for optimized initial values with those obtained by so called spinup. We assess the sensitivity of optimal parametersets to temperature, diffusivity etc. i.e. the forcing.
Samuel, H. & Evonuk, M. Fast and Accurate Modeling of Advection in Geophysical Flows with Dynamic Implicit Surfaces Session 9 Tuesday 8  poster file 
Abstract: Advection is one of the major processes that commonly acts on various scales in nature (core formation, mantle convective stirring, multi-phase flows in magma chambers, subduction dynamics, salt diapirism ...). While this process can be modeled numerically by solving conservation equations, various geodynamic scenarios involve advection of quantities with sharp discontinuities. Unfortunately, in these cases modeling numerically pure advection becomes very challenging, in particular because sharp discontinuities lead to numerical instabilities, which prevent the local use of high order numerical schemes.
Several approaches have been used in computational geodynamics in order to overcome this difficulty, with variable amounts of success. Despite the use of correcting filters or non-oscillatory, shock-preserving schemes, Eulerian (fixed grid) techniques generally suffer from artificial numerical diffusion. Lagrangian approaches (dynamic grids or particles) tend to be more popular in computational geodynamics because they are not prone to excessive numerical diffusion. However, these approaches are generally computationally expensive, especially in 3D, and can suffer from spurious statistical noise.
As an alternative to these aforementioned approaches, we have tested the Level Set and Particle Level set methods against the robust and popular tracer-in-cell method on well-known 2D thermochemical benchmarks and typical 3D convective flows. Although the Level Set method has some advantages over the tracer-in-cell method it generally fails to resolve fine (sub-grid) scale features and requires high resolution in order to conserve mass accurately.
However, the use of Lagrangian tracers in the particle level set method yields much more accurate solutions of purely advective transport. In every case we ran we found that the Particle Level Set method accuracy equals or is better than other Eulerian and Lagrangian methods, and leads to significantly smaller computational cost, in particular in three-dimensional flows, where the reduction of computational time for modeling advection processes is most needed.
Scandura, D., Currenti, G., Bonaccorso, A., Di Stefano, A. & Del Negro, C. Finite-Element modeling of static stress changes induced by recharging and intrusive phases at Etna volcano Session 4 Tuesday 8  poster file 
Abstract: The last thirty years have been a very active time in the eruptive and seismic history of Mt Etna. Significant correlation between the occurrence of large earthquakes along the main volcano-tectonic structures and volcanic unrest has observed throughout the analysis of eruptive sequences and local seismicity. The most outstanding tectonic features at Mt Etna are clearly recognizable on the eastern flank of the volcano, where the clearest morphological evidence of active faulting exists. The volcano is characterized by a structural setting resulting from a complex regional tectonic with a compressive regime along a near N-S trend and an extensional regime tending approximately E-W observable along the eastern coast. Inside the volcano edifice the eastern flank shows an eastward sliding, which is delimited in the northern border by the Pernicana Fault System (PFS). PFS is one of the most active tectonic structures in the Etna area, characterized both by aseismic continuous ``slow'' movements associated with the eastern flank sliding and by shallow earthquakes. We examine the possible relations between PFS ruptures and volcanic unrest during the most energetic earthquakes (about M4) occurred in the last three decades. The effect of magmatic intrusions and inflation processes on the PFS are studied through the estimate of stress redistribution using deformation models based on the Finite Element Method (FEM) and constrained by high-quality geodetic data. The results and the principal conclusions suggested by the modeling reveal that the loading and rupture of PFS occur mostly as a response to accommodate the stress changes induced by pre-eruptive injection of the magma when this penetrates inside the northern sector of the volcano. The results from FEM have shown that the effect of topography, crust rigidity layering and the presence of fault heterogeneities are significant and cannot be neglected in quantitative analyses of faulting or earthquakes.
Schachtschneider, R., Hayn, M. & Holschneider, M. Anisotropy of the geomagnetic field detected with directional Poisson wavelets on the sphere Session 5 Thursday 10  poster file 
Abstract: The calculation of geomagnetic field models is a fundamental task in geomagnetism. Up to now models have been calculated with isotropic priors about the magnetic field. As the geomagnetic field is not isotropic that approach might not be optimal. The main goal of this project is to analyze the recently obtained high quality models of the geomagnetic field based on CHAMP satellite data with the help of the newly developed directional Poisson wavelets on the sphere. This will allow us to asses and quantify the position and scale dependent locally anisotropic correlation structure of the crustal field. In a first step we apply the wavelet transform on the radial component of the geomagnetic field, obtaining information about dominant directions and scales of the locally anisotropic features. By constructing a full covariance matrix this information can be used to increase the accuracy of the field inversions if this kind of non-isotropic information are taken into account as prior information for the Bayesian inversion process.
Schaeffer, N. Towards a faster Spherical Harmonic transform Session 9 Monday 7  oral file 
Abstract: Several fast spherical harmonics algorithm have been proposed in the past decades. However, their use has not spread, showing that they still suffer of limitations.
First, their speed advantage only appears for quite large cutoffs degree (around 512). They often lack a vector transform, have no flexible interface, and some of them have stability problems.
On the other hand, the Gauss-Legendre algorithm requires less grid points for an exact quadrature and is computationnaly very efficient if one precompute everything.
The new algorithm and implementation presented here is aimed at numerical simulations. It is not a fast algorithm as it has almost the same complexity as the Gauss-Legendre one. However, it is faster, and mostly efficient for anisotropic truncation, that is relevant for fast rotating fluids. It makes use of a regular grid, while using almost the same amount of point for an exact quadrature as the Gauss-Legendre algorithm.
We compare its speed with a highly opitmized Gauss-Legendre transform, and make other suggestions for improving the speed of numerical simulations in spherical shells. We also propose a highly optimized reference implementation, with both scalar and vector transform, and a flexible interface.
Scholz, C.H. Large Earthquake Triggering, Clustering, and the Synchronization of Faults Session 4 Monday 7  oral file 
Abstract: Large earthquakes are sometimes observed to trigger other large earthquakes on nearby faults. The magnitudes of the calculated Coulomb stress transfers presumed to cause the triggering are $10^-2- 10^-3$ of the earthquake stress drops and the triggering delay times are similarly small with respect to the natural recurrence time of the earthquakes. This requires that both faults be simulaneously very close to the ends of their seismic cycles. Paleoseismological data show that for the same regions prior large earthquakes have occurred on nearby faults in clusters in space and time separated by long quiescent periods. Both observations suggest that synchronization is occurring between faults. Theory and observations indicate that synchronization can occur between nearby faults with positive stress coupling and intrinsic slip velocities within an entrainment threshold. In the south Iceland seismic zone, the central Nevada seismic belt, and the eastern California shear zone, several synchronous clusters, that apparently act independently, can be recognized. This behavior is the equivalent, in three dimensions, of the phase locking of fault elements that results in the seismic cycle of individual faults being dominated by large ``characteristic'' earthquakes, and for synchronization of fault segments along a given fault. Rupture patterns of repeated individual earthquakes or earthquake clusters are not identical in either the two or three dimensional cases. The state of this system, which exhibits strong indications of synchrony without exact repetition, may be called fuzzy synchrony.
Schoof, C. Drainage network formation under glaciers
Session 6 Tuesday 8  oral file 
Abstract: Ice flow velocities in glaciers and ice sheets are known to depend sensitively on conditions at the base of the ice. In particular, switches in subglacial drainage network configurations can have a major impact on glacier and ice sheet sliding, as subglacial water at high pressure tends to favour rapid sliding. Abrupt changes in the style of subglacial drainage have therefore often been invoked as an explanation for the onset of fast ice flow in ice sheets in the form of ice streams (narrow bands of ice that move at speeds orders of magnitude larger than what can be explained by shearing in the ice), and for surges in glaciers (relaxation oscillations in ice velocity and ice thickness). To date, spatially extended models that could test some of the hypotheses about drainage networks have been conspicuously lacking. Here we review the state of the art in subglacial drainage, and present a new model that captures the dynamics of water-filled cavities and subglacial channels in a unified mathematical description, which is based on a network of discrete drainage elements. We show that this model can capture the formation and disappearance of an arborescent, channelized network under changes in water input and sliding velocity. We also show that this mechanism is insufficient to explain surging behaviour, and show how feedbacks with water storage are essential in causing the abrupt drainage transitions observed in some glacier surges.
Sciascia, R., Pasquero, C. & Provenzale, A. Settling plankton settling Session 6 Tuesday 8  oral file 
Abstract: Plankton live in the upper ocean, where light allows photosynthesis. They often are denser than water and sink to deeper levels. It has been hypothesized that ocean turbulence favors plankton suspension but studies show contrasting results, indicating that settling can be even faster than in still water. Here, we present theoretical and numerical studies that indicate under what conditions turbulence induces a significant suspension of heavy particles and allows for survival of plankton in the euphotic layer. The study is also important for determining detrital settling velocity in the ocean water column, a key parameter in carbon cycle modeling.
Scotti, A., Galimberti, R. & Ruffo, P. Darcy based models for oil generation and expulsion from source rock. Session 9 Tuesday 8  poster file 
Abstract: Oil industry is resorting more and more often to mathematical and numerical modeling to locate new oil and gas reservoirs, and to exploit the existing ones at their best. The numerical simulation of generation and migration processes can support geologic analysis in tracking hydrocarbons ``from source to trap'' and therefore reduce the risk in oil exploration. This work deals with the mathematical modeling and numerical simulation of the generation of hydrocarbons, and their expulsion from the source rock, with the aim of predicting the amount and quality of oil and gas generated in a region during million years. The chemical reactions of generation are strongly coupled with selective retention processes which act as filters altering the chemical composition of the oil, and with the fluid flow in the rock. We model the flow of water and oil in the rock as a two phase multicomponent Darcy flow. The numerical approximation of the problem relies on a time splitting approach to decouple the equations. Mixed hybridized finite elements are employed for the Darcy problem for an accurate velocity field in presence of strong inhomogeneities in the permeability. The degenerate parabolic equation for the oil phase saturation is solved with an advection-diffusion-reaction splitting, employing the Godunov method for the advective part and expanded mixed finite elements for the degenerate diffusion part. Integration schemes that preserve that are conservative and preserve the positivity of the solution are proposed for the integration of chemical reactions. Moreover, the retention processes introduce a discontinuity in the system that requires a special treatment from the analytic and numerical point of view. The result of monodimensional and bidimensional simulations show that generation, retention and fluids flow strongly influence each other and in particular, the inclusion of retention processes delays the expulsion and enriches the expelled hydrocarbons in gas and saturates, in agreement with typical experimental data collected at the wells.
Selva, J., Marzocchi, W., Papale, P., Civetta, L., Del Pezzo, E. & Sandri, L. Community-based short-term eruption forecasting at Campi Flegrei Session 7 Thursday 10  oral file 
Abstract: A key element in emergency preparedness is to define in advance tools to assist decision makers and emergency management groups during crises. Such tools must account for all of the expertise and available scientific knowledge. During a pre-eruptive phase, the key for sound short-term eruption forecasting is the analysis of the monitoring signals. This involves the capability (i) to recognize anomalous signals, relating single or combined anomalies to physical processes, and to assigning them probability values, and (ii) to quickly provide an answer to the observed phenomena even when unexpected. Here we present a $>$ 4 years long process devoted to define the pre-eruptive Event Tree (ET) for Campi Flegrei. A community of about 40 experts in volcanology and volcano monitoring participating to two Italian Projects on Campi Flegrei funded by the Italian Civil Protection, has been constituted and trained with periodic meetings on the statistical methods for the ET definition, based on the model BET$$EF (Marzocchi et al., 2008). Model calibration has been carried out through public elicitation sessions, preceded and followed by devoted meetings and web forum discussions on monitoring parameters, their accuracy and relevance, and their potential meanings. The calibrated ET allows anomalies in the monitored parameters to be recognized and interpreted, assigning probability values to each set of data. This process de-personalizes the difficult task of interpreting multi-parametric sets of data during on-going emergencies, and provides a view of the observed variations that accounts for the averaged, weighted opinion of the scientific community. We first discuss all mathematical aspects of the probabilistic model used to translate available data, expert opinions, and monitoring measures into a probability assessment. Then, we discuss the results of this 4-years process. Notably, the results draws a clear and rational picture of what the community, as a whole, considers well established and then useful for practical applications. The elicitation process enables us to filter out possible biases introduced by personal and/or specialist beliefs that the community (still) considers not enough constrained or too speculative. Such a common picture can be also particularly useful to guide future implementations of the monitoring network, as well as research investments aimed at substantially improving the capability to forecast the short-term volcanic hazard.
Shabanzadeh, H., Bagheri, G., Ghasemi, M. & Nekoofar, M. Investigation of the results of seismography to determine shear wave velocity and soil shear strength parameters Session 8 Thursday 10  poster file 
Abstract: Correlation between the results of field tests such as impact and standard penetration (SPT), and dynamic parameters of soil layers such as shear wave velocity(Vs) has been verified and several empirical relations were presented in different texts. In the absence of field geophysical methods or laboratory tests, these relations should be used with caution and require individual experience. Besides, most of these relations are valid only when the SPT value (Nspt) is less than 50 impacts. In this paper, correlation between soil shear strength parameters (based on Mohr-Coulomb failure criterion) and soil shear wave velocity obtained from cross hole seismography method on soil layers by boreholes of maximum 30 meters in depth (in accordance to Iranian Code of Design 2800), are investigated. These parameters are based on direct shear tests on the samples taken from boreholes drilled in Zanjan Petrochemical Plants site investigations.
Simons, F.J. Slepian functions and their use in geophysical signal estimation and spectral analysis Session 9 Tuesday 8  poster file 
Abstract: All geology is local, yet many modern techniques to analyze geophysical signals, from satellite gravity and the geoid, the lithospheric magnetic field, to the making of seismological earth models, heavily rely on a global function basis --- the spherical harmonics. When the available data are incomplete and/or noisy, or when the desired models or the features to be extracted are geographically localized (earthquakes, ocean topography, regional-scale anomalies, etc), these have easily demonstrable shortcomings. Enter Slepian functions: firmly rooted in the theory of potential fields, they are designed to combine the best of both the spatial and spectral viewpoints by achieving optimal spatio-spectral localization. A basis of Slepian functions for geophysical signals is sparse, and when used as windows to taper data for spectral analysis, they allow for the robust extraction of localized power.
It is a well-known fact that mathematical functions that are timelimited (or spacelimited) cannot be simultaneously bandlimited (in frequency). Yet the finite precision of measurement and computation unavoidably bandlimits our observation and modeling scientific data, and we often only have access to, or are only interested in, a study area that is temporally or spatially bounded. In the geosciences we may be interested in spectrally modeling a time series defined only on a certain interval, or we may want to characterize a specific geographical area observed using an effectively bandlimited measurement device. It is clear that analyzing and representing scientific data of this kind will be facilitated if a basis of functions can be found that are spatiospectrally concentrated, i.e. localized in both domains at the same time. Here, we give a theoretical overview of one particular approach to this concentration problem, as originally proposed for time series by Slepian and coworkers, in the 1960s. We show how this framework leads to practical algorithms and statistically performant methods for the analysis of signals and their power spectra in one and two dimensions, and on the surface of a sphere.
Among the applications that I will discuss are the analysis of the time-variable gravity field for the recovery of coseismic gravity perturbations, the sparse analysis and representation of the lithospheric magnetic field, the recovery of the power spectral density of the cosmic microwave background radiation, and the filtering of the geoid for geodynamic applications.
Skiba, Y.N. & Parra-Guevara, D. Pollution level assessment and control of emission rates Session 6 Tuesday 8  poster file 
Abstract: The main objective of the mathematical modeling in environment protection is the prediction of concentrations of pollutants and development of the methods protecting the people from dangerous pollution levels. In this work, we suggest methods for solving important problems related to point sources (industrial pollution, oil spill problem, aquatic system cleaning, etc.) and linearly distributed (automobile) sources. The adjoint approach utilizing direct and adjoint pollution concentration estimates in selected zones is used. The direct estimates require solving the pollution transport model and enabling a comprehensive analysis of ecological situation in the whole area. By contrast, the adjoint estimates use solutions to an adjoint model and depend explicitly on the number, positions and emission rates of the sources and on the initial distribution of pollutants. Besides, the adjoint model solution serves in such estimates as influence (weight) function providing valuable information on the role of each source and initial data in polluting a selected zone. These estimates are effective in studying the model response to variations in the emission rates and initial conditions, and in developing control strategies. The control over emission rates of industrial plants develops quantitative criteria which permit to avoid violations of existing sanitary norm by means of margin reductions of the original industrial emission rates. Such criteria are designed taking into account the complexity inherent to the processes of dispersion and transformation of pollutants in the atmosphere, the number of point sources to control, their locations in a region and the corresponding ecological laws. With adjoint estimates we have developed both non-optimal and optimal control strategies. The method is illustrated by a simple example with only two plants in operation. Another example shows that the same adjoint approach can be used to build an efficient control strategy for aquatic system cleaning. A method for detecting the enterprises, which violate the emission rates prescribed by a control, is given. It serves to detect the infractions and apply sanctions. A method of determining an optimal position of a new factory in the region is also described. The method is illustrated by a simple example when only two types of climatic wind, eleven enterprises and three ecologically most important zones are taken into account.
Smith, L.A. Extracting Insight from Predictions of the Irrelevant: Can the Diversity in Our Models Inform Our Uncertainty of the Future? Special Session Thursday 10  oral file
Abstract: The open question of whether or not ``physically reasonable'' solutions of the Navier-Stokes equations even exist is one of great mathematical interest and nontrivial monetary reward. For a mathematician such a question is of interest for its own sake. For the geophysicist (or policy maker) interested in the Earth System, however, the interesting questions focus on our ability to adequately simulate partial differential equations thought to describe the climate system. Perfection is a non-starter here; even before we compile our models, we know they are not ``valid'' representations of the Earth System, as we have built blatant untruths into each and every one of them. So the question is not whether they are perfect but whether they are useful. Can they predict only the irrelevant? Given that we know our models are empirically inadequate, precise probability forecasts based on model-output should not be interpreted naively as decision-relevant probabilities; how then might our simulations provide insight?
Section 11 of Tukey's 1962 paper, ``The Future of Data Analysis'', titled Facing Uncertainty, was often quoted by Albert Tarantola. Although written for 1960's data analysts, it is of equal value to geophysical modellers of this century: ``The most important maxim for data analysis to heed, and one which many statisticians seem to have shunned, is this: `Far better an approximate answer to the right question, which is often vague, than an exact answer to the wrong question, which can always be made precise.''' Nonlinearity exposes the limitations of least squares methods, as Tarantola stressed in ``Inverse Problem Theory''. State uncertainty in nonlinear systems undermines the use of least squares methods in a manner not unlike the way structural model inadequacy undermines the use of model-based probability distributions. What might it mean to put a prior on an empirically vacuous model-parameter? Could it be rational to base policy on the posterior probability distribution of a model-variable, knowing that the model (class) considered was empirically inadequate? Is the theory of probability, as such, irrelevant for decision-support in extrapolation problems like climate change?
``Kitchen sink'' models aim to provide the ``best available'' answer to an intractable problem, by including every process that will fit into a computational package and still allow that package to execute before the grant ends. Without question one gets seemingly precise answers in terms of probability distributions: an exact answer to the wrong question. But what is it that corresponds to Tukey's ``approximate answer to the right question'' in an extrapolation problem where model complexity is constrained by computational power (as in climate modelling)? We will examine the various flavours of ``vague'' on offer. What does one judge an ``approximate answer'' when the aim is insight yet the accessible model class in inadequate? Is a credible model of a less Earth-like planet (say, a planet with Earth-like topography and no oceans) more informative than a model that sorta looks like the Earth as we know it, but with physical processes (clouds) and boundary conditions (the Andes) significantly different from well-understood physics and the known topography?
Questioning the policy relevance of model-based ``implied-probabilities'' leads to the maintenance of parallel pure and applied research programs over the timescales on which climate policy will demand scientific support (five years and fifty). This also raises the possibility of moving away from probabilities, perhaps towards non-probabilistic odds. How wide is the gulf between programs advancing scientific understanding and those constructing a basis for evidence-based policy making? If applying cost-benefit analysis directs us to the wrong questions, might designing climate experiments to inform a risk management framework answer the right question, approximately? One cannot take today's climate models literally given the range of systematic errors even in global mean temperature, a range far exceeding the observed centennial increase. In terms of insight, however, today's climate models play a critical, positive, informative role. Not one of them suggests emergent feedbacks that cast doubt on the fundamental insight drawn from the hierarchy of previous state-of-the-art models: that increasing greenhouse gases will result in a warmer planet. As long as our ensemble forecasts from different state-of-the-art models disagree (in distribution), then we know that the details of model structure matter. Nevertheless I will argue that we can draw policy relevant insights from today's state-of-the-art models and take action rationally, even in the case that they ``predict the irrelevant''.
Sokol, Z. & Zacharov, P. Assimilation of extrapolated radar reflectivity into a NWP model and its impact on forecasts of convective precipitation Session 8 Thursday 10  poster file 
Abstract: A method assimilating observed and extrapolated radar reflectivity is presented. The method is focused on a nowcasting of heavy convective precipitation. The data are assimilated into the COSMO NWP model, which is integrated with a horizontal resolution of 2.8 km, and the extrapolation method is based on the COTREC algorithm. The aim of the combination of observed and extrapolated data is to utilize the fact that extrapolation techniques usually yield good forecasts for short lead times. The accuracy of hourly precipitation forecasts will be evaluated and compared with forecasts using only observed radar reflectivity. Current results show that the extrapolated data improve precipitation forecasts for the second and third lead hours. Convective precipitation is difficult to forecast in a deterministic way. Therefore the simulated forecasts will be used to estimate uncertainty of the forecast and its dependence on the lead time and horizontal scale.
Solodoch, A., Boos, W., Kuang, Z. & Tziperman, E. Excitation of slow MJO-like Kelvin waves in the equatorial atmosphere by Yanai wave-group via WISHE-induced convection Session 6 Tuesday 8  oral file 
Abstract: The intraseasonal Madden-Julian oscillation (MJO) involves a slow eastward-propagating signal in the tropical atmosphere which significantly influences climate yet is not well understood despite significant theoretical and observational progress. We study the atmosphere's response to nonlinear ``Wind Induced Surface Heat Exchange'' (WISHE) forcing in the tropics using a simple shallow water atmospheric model. The model produces an interestingly rich interannual behavior including a slow, eastward propagating equatorial westerly multiscale signal, not consistent with any free linear waves, and with MJO-like characteristics. It is shown that the slow signal is due to a Kelvin wave forced by WISHE due to the meridional wind induced by a Yanai wave group. The forced Kelvin wave has a velocity similar to the group velocity of the Yanai waves, allowing the two to interact nonlinearly via the WISHE term while slowly propagating eastward. These results may have implications for observed tropical WISHE-related atmospheric intraseasonal phenomena.
Song, H., Hoteit, I., Cornuelle, B.D. & Subramanian, A.C. An Adaptive Approach to Mitigate Background Covariance Limitations in the Ensemble Kalman Filter Session 8 Thursday 10  poster file 
Abstract: A new approach is proposed to address the background covariance limitations arising from under-sampled ensembles and unaccounted model errors in the ensemble Kalman filter (EnKF). The method enhances the representativeness of the EnKF ensemble by augmenting it with new members chosen adaptively to add missing information that prevents the EnKF from fully fitting the data to the ensemble. The vectors to be added are obtained by back-projecting the residuals of the observation misfits from the EnKF analysis step onto the state space. The back-projection is done using an optimal interpolation (OI) scheme based on an estimated covariance of the subspace missing from the ensemble. In the experiments reported here, the OI uses a pre-selected stationary background covariance matrix, as in the hybrid EnKF/3DVAR approach, but the resulting correction is included as a new ensemble member instead of being added to all existing ensemble members.
The adaptive approach is tested with the Lorenz-96 model. The hybrid EnKF/3DVAR is used as a benchmark to evaluate the performance of the adaptive approach. Assimilation experiments suggest that the new adaptive scheme significantly improves the EnKF behavior when it suffers from small size ensembles and neglected model errors. It was further found to be competitive with the hybrid EnKF/3DVAR approach, depending on ensemble size and data coverage.
Stewart, L.M., Dance, S.L. & Nichols, N.N. Estimation and modelling of observation error correlations for numerical weather prediction Session 8 Friday 11  oral file 
Abstract: Variational data assimilation algorithms are designed using the ideas of Bayesian statistics to compute the most probable state of the atmosphere taking into account the uncertainty in the available observational and forecast data. Thus for good results, accurate specification of the forecast and observation error distributions is vital. Remote sensing data often have correlated errors, arising from observation pre-processing and mismatches between observation and model resolutions. However, the correlations are typically ignored in operational numerical weather prediction. This approximation is made since the details of the correlation structure are often unknown. It also allows a simplification of the calculations and a reduction in computational cost. Problems are avoided by some combination of data thinning, bias correction, and increasing the size of the observation error variance used in the assimilation. Unfortunately, these measures do not fully exploit the observations and significant information is lost in the assimilation.
Our recent work using the IASI instrument with the UK Met Office assimilation system has shown that it may be possible to estimate observation error correlations using data assimilation output diagnostics. In this paper, we present some simplified model studies examining how well we can expect such an approach to approximate the true observation error covariances. The results are extremely encouraging and show that, if care is taken to include sufficient statistics, diagnostic approaches to estimation of error correlations can produce a good estimate of the true observation error correlations. We also investigate computationally efficient approximations of correlation matrices that may be feasible for operational applications and their effect on analysis accuracy. The results demonstrate that it is often better to model error correlation structure incorrectly than not at all, and that including error correlation structure in data assimilation algorithms is both feasible and beneficial.
Storto, A., Masina, S., Dobricic, S. & Di Pietro, P. Modelling global ocean background-error covariances via ensemble simulations for use in reanalysis system Session 8 Thursday 10  poster file 
Abstract: Variational assimilation techniques for use within global ocean reanalysis systems require the estimation of background-error covariances. In order to simulate the true error evolution in the assimilation and forecast step of our global ocean analysis system, we performed ensemble variational assimilation by adding independent random Gaussian perturbations to the entire observation vector - formed of in-situ and space-borne observations - multiplied by the observational error covariance matrix. Additionally, we perturbed surface forcing as well (wind stress and sea-surface temperature) and introduced an auto-correlated multi-variate perturbation in the non-advective terms of the state parameters tendencies along with the ocean model integrations. The ensemble size consisted of 6 members and covered the period from 1993 to 2005. This strategy allowed us to depict the horizontal and vertical distribution of model error variances and horizontal correlation length-scale, for temperature and salinity apart. The background-error covariances structure showed evidence of i) the need of separately accounting for temperature and salinity error horizontal correlations, as the former primarily respond to the mesoscale dynamics while the latter to the fresh water budget (e.g. increase of error variances are found in correspondence of river mouths, etc.); ii) the need of accounting for seasonal variations within the background-errror vertical covariances; iii) the need of allowing anisotropic correlation length-scale in the Tropical belt, whereas the zonal horizontal correlation scales are found constantly longer. We used these results to further model the background-error covariances within our deterministic analysis system: vertical background-error covariances are modelled from monthly bivariate EOFs of temperature and salinity computed from the ensemble error dataset at native ensemble simulation resolution. Horizontal correlations are obtained from a four-iteration application of a recursive filter, whose coefficients allow for locally smoothed, seasonal and vertical variations of the correlation radius, further to increased zonal correlation length-scales in the Tropical region. Preliminary results prove that the analysis system is benefited by the new background-error formulation in terms of verification skill scores.
Subramanian, A.C., Hoteit, I., Neef, L. & Song, H. Implementation of the nonlinear filtering problem and balanced dynamics Session 8 Thursday 10  poster file 
Abstract: We investigate the role of the linear analysis step of the ensemble Kalman filters (EnKF) in exciting spurious gravity waves in atmospheric motion models. This is achieved through the comparison of the behaviors of an EnKF and a fully nonlinear particle-based filter (PF) with a simple model of balanced dynamics. The filters have very similar forecast step but their analysis steps are different. More specifically, the analysis step of the PF generalizes the optimality of the EnKF analysis to non-Gaussian distributions. The model admits a chaotic vortical mode coupled to a comparatively fast gravity wave mode. It can be initialized such that it evolves on a so-called slow manifold, where the fast motion is suppressed. It can also be initialized such that the fast varying variables are diagnosed from the slow varying variables. This is called the slaved mode of the model. To determine how well the nonlinear analysis step preserves dynamical balance in the solution, identical twin assimilation experiments are performed, wherein the true state is balanced, but the observational errors project onto all degrees of freedom, including the fast modes. EnKF and PF capture the variables in slow manifold well since, once the variables are attracted towards the slow manifold, they stay there. PF captures slaved modes better implying nonlinear jumps in dependent variables are captured better. Results also include testing other PFs on the model to see which PF performs best in the fully nonlinear Lorenz86 model.
Swaters, G.E. Cross-equatorial flow of grounded abyssal water masses Session 3   withdrawn file 
Abstract: Oceanographic data clearly show and OCGM simulations unambiguously illustrate that the abyssal circulation is quite happy to cross the equator. Away from the equator, but still removed from their source regions (where equatorward motion is initiated via the Sverdrup vorticity balance associated with Stommel-Arons Theory), these deep meridional flows are geostrophic and appear to be topographically steered. Theoretical analysis of equator-crossing abyssal flows is challenging not only because geostrophy breaks down but these flows apparently do not conserve potential vorticity. The two existing dynamical paradigms for the underlying dynamical balances that permit cross-equatorial abyssal flow can be summarized as eddy-forming (nonlinearity balances the pressure gradient) and viscous (bottom friction balances the pressure gradient). However, the low current speeds and small Ekman numbers associated with these currents suggest neither of these processes is phenomenologically relevant. In this talk we will show that the slow steady inviscid cross-equatorial flow of a grounded abyssal water mass within a meridional trench (i.e., a meridional channel with sloping boundaries) can be described by the linear non-hydrostatic reduced-gravity equations in which the horizontal or tangential components of Earth's spin vector are retained. Such a model predicts that away from the equator, the abyssal water mass (assumed to lying directly on top of the sloping channel boundary) is topographically-steered and in geostrophic-hydrostatic balance. As the equator is approached, the abyssal layer height decreases (in accordance with the conservation of planetary geostrophic vorticity) until very near the equator where the horizontal components of the Coriolis acceleration become important and initiate down slope and zonal motion in which the abyssal water mass migrates to the other side of the channel and continues flowing meridionally eventually re-establishing a topographically-steered geostrophic-hydrostatic balance in the other hemisphere.
Szymczak, P. & Ladd, T. Dissolution in porous media and fractures: initial vs moving front instability Session 2 Monday 7  oral file 
Abstract: A reactive fluid flowing through porous or fractured rock and dissolving the rock matrix will trigger an instability, leading to spontaneous formation of pronounced channels or wormholes. We present a linear-stability analysis of this system and show that there are two different instabilities. One is associated with an initial, uniform porosity state whereas the other describes a steadily propagating one-dimensional dissolution front. We discuss the origin of both instabilities and the physical conditions under which they can be observed. In particular, it is argued that the ``initial state" instability is relevant to the dissolution of fracture in carbonate rocks, giving rise to the formation of limestone caves.
Talagrand, O. Data Assimilation in Geophysics. An Update Session 8 Friday 11  oral file 
Abstract: Assimilation of observations, which originated from the need of defining initial conditions for numerical weather forecasts, and then progressively extended to many diverse applications, has become a basic component of numerical modeling of geophysical systems. Assimilation of observations is best expressed as a problem in Bayesian estimation. Namely, determine the probability distribution for the state of the system under observation, conditioned to the available information. That information essentially consists of the observations proper, and of the physical laws that govern the system, available in practice in the form of discretized numerical model. For a number of reasons, starting with the very large dimension of the corresponding state space (up to 108 for numerical weather prediction models), Bayesian estimation is not possible in practice, and one has to restrict to more modest goals. The main two present types of assimilation algorithms, Variational Assimilation (VA) on the one hand, and Ensemble Kalman Filter (EnKF) on the other, are empirical extensions, to weakly nonlinear situations, of Best Linear Unbiased Estimation (BLUE), which achieves Bayesian estimation in linear and Gaussian cases. Variational Assimilation globally adjusts a model solution to observations distributed over a period of time. Ensemble Kalman Filter evolves in time an ensemble of estimates of the system state, and updates those estimates as new observations become available. The update is performed according to the theory of BLUE. The performances of both types of algorithms is presented, and their respective advantages and disadvantages is discussed. Particle Filters, like EnKF, evolve an ensemble of state estimates, but use an exactly Bayesian updating scheme. However, their cost has so far been too high for use for large dimension systems. Their performance, and possible future evolution, are discussed.
Tanimoto, T. Mechanical Interactions among the Atmosphere, the Oceans and the Solid Earth Session 3   withdrawn file 
Abstract: Seismologists have learned in the past decade that there are important seismic noise in seismograms, generated by the atmosphere and the oceans. For example, the noise cross-correlation approach (Campillo and Paul, 2003) has been so successful that the technique has now become a standard approach for Earth-structure study. This approach is just one example of techniques to take advantage of seismic noise, and there could be many other approaches if we understood the mechanisms of interactions among the atmosphere, the oceans and the solid earth.

The purpose of this paper is to summarize the nature of seismic noise for frequencies between 1 mHz and 1 Hz and discuss outstanding issues for seismology, in particular for earth-structure retrieval. We will focus on three different frequency bands:
(1)2 mHz and lower frequencies: For this frequency band, seismic noise shows rapid amplitude decay with frequency because the source of noise is density change in the atmosphere which is turbulent in the frequency band. This phenomenon was examined from 1970s to early 90s but its investigation seems to have slowed since then. We now have much larger amount of high-quality barometer and seismic data and thus it is possible to gain more insight for the nature of seismic noise for this frequency band.
(2)The Hum (3-15 mHz): Enigmatic background free oscillations of the Earth were noted in 1998 (Suda and Nawa, 1998; Tanimoto et al., 1998). A long list of investigation since then has tracked its source to the long-period (~a few hundred seconds) ocean waves (Rhee and Romanowicz, 2004; Tanimoto, 2005; Webb, 2007), although the detailed processes are not clear yet. Most signals are spheroidal modes (Rayleigh waves) but existence of toroidal mode signals was pointed out recently (Kurrle and Widmer-Schnigdrig, 2008). Mechanisms of excitation are still vague for spheroidal modes and those for toroidal modes are not at all clear.
(3)Microseisms (0.05 - 0.5 Hz): This is the most energetic band among the entire seismic band (1 mHz - 10 Hz) and the successful cross-correlation approach mostly takes advantage of this band. The mechanism of excitation seems to be the one proposed by Longuet-Higgins (1950) for Rayleigh wave signals. But there are occasional strong Love wave signals for which the mechanism of excitation is not clear yet.
For each of the three frequency bands, I will provide brief summary on our understandings and potential practical uses.

Tarquini, S., Favalli, M. & Fornaciai, A. DEM-due uncertainties in gravity-driven mass movements on Earth surface Session 7 Thursday 10  poster file 
Abstract: The gravitational force, on Earth, drives a large variety of mass movements bounded by the ground surface. As a first approximation, these movements are always driven along the steepest descent path of the surface. As a consequence, the morphology of the bedrock is crucial in determining the path of surficial flows, and the characterization of its mathematical representation (Digital Elevation Models - DEMs) is of primary importance in environmental studies. We present two novel maps characterizing in an effective way the interplay between a gravitational mass flow and the surface over which the mass moves. In the first map (type 1 map), each point express the susceptibility of mass movements starting from this point to stay channellized or to spread over the surface. In the second map (type 2 map), each point express how the resulting mass movements remain stable as the starting point is shifted by a given distance. We take lava flows as mass movement test case and the flanks of the Mount Etna volcano (Sicily, South Italy) as surface test case. We perform our analysis by elaborating a large database of lava flow simulations carried out by using the DOWNFLOW code, which has been extensively used for lava flow forecasting at several basaltic volcanoes (Mount Etna, Mount Nyiragongo and Mount Cameroon). These simulations are triggered from about 70,000 possible future vents located at the nodes of a 80 m-spaced grid covering a wide area of the volcano, from the summit craters down to the lower flanks. These maps are of use in hazard and risk assessment because at the onset of a new eruption, as far as the new forming vent is set on the maps, they allow the characterization of two important properties of the expected lava flow. Type 1 maps tell, at each point, if the flow venting at that point or reaching that point is expected to spread or not over topography, helping in designing barriers for lava flow diversion. Type 2 maps provide the DEM-due uncertainty on lava flow path as a function of the precision in vent position, and is very useful to design a new hazard map. Maps like the one presented here can be derived for different kinds of gravity-driven mass movements over the Earth surface in other areas.
Todesco, M. & Rinaldi, A.P. The rocks and the fluids. The complex sound of the hydrothermal activity Session 2 Monday 7  oral file 
Abstract: Volcanic gases and aqueous fluids permeate the rocks in volcanic areas and feed a wide range of surface phenomena, such as hot springs, fumaroles, and diffuse degassing. Measurements of fluid temperature and composition provide information on the source of discharged fluids, and are commonly carried out in volcanic areas to gather information on the state of magmatic degassing at depth. The presence of hydrothermal circulation, however, generates different kinds of signals, including changes in gravity, electrical resistivity, or ground deformation. These signals form as temperature, pressure, and phase distribution change through time and space, possibly due to a variable activity of different fluid sources. At the same time, the circulating fluids are affected by the nature and properties of the subsurface rocks. The extent and velocity of fluid propagation, and the magnitude of surface phenomena are governed by the rock permeability and by the presence and geometry of fractures. The hydrothermal activity and the related signals, result from the complex interactions between one or more fluid sources and the porous medium through which the fluids propagate. Numerical modeling of multi-phase and multi-component hydrothermal circulation is used to simulate the signals generated by the hydrothermal activity under various choices of rock properties and boundary conditions. Result highlight the importance of heterogeneous rock permeability and illustrate the different sensitivities of the considered observable parameters to changes of the source and of the medium properties. Modeling results may provide a useful framework for the interpretation of monitoring data collected in volcanic areas.
Tosi, N., Yuen, D. & Cadek, O. Dynamical consequences in the lower mantle with the post-perovskite phase change and strongly depth-dependent thermodynamic and transport properties Session 9 Tuesday 8  poster file 
Abstract: We have carried out numerical simulations of large aspect ratio 2-D mantle convection with the deep phase change from perovskite (pv) to post-perovskite (ppv). Using the extended Boussinesq approximation for a fluid with temperature- and pressure-dependent viscosity, we have investigated the effects of various ppv phase parameters on the convective planform, heat transport and mean temperature and viscosity profiles. Since ppv is expected to have a relatively weak rheology with respect to pv and a large thermal conductivity, we have assumed that the transition from pv to ppv is accompanied by both a reduction in viscosity by 1 to 2 orders of magnitude and by an increase in thermal conductivity by a factor of 2. Furthermore, we have analyzed the combined effects of a strongly decreasing thermal expansivity in pv and steeply increasing thermal conductivity according to recent evidence from high-pressure experiments and first-principle calculations. As long as the thermal expansivity and conductivity are constant, ppv exerts a small but measurable effect on mantle convection: it destabilizes the D" layer, causes focusing of the heat flux peaks and an increase of the average mantle temperature and of the temporal and spatial frequency of upwellings. When the latest depth-dependent thermal expansivity and conductivity models are introduced, the effects of ppv are dramatic. On the one hand, without ppv, we obtain a very sluggish convective regime characterized by a relatively cool mantle dominated by large downwellings that tend to stagnate beneath the transition zone. With ppv, on the other hand, we observe an extremely significant increase of the average mantle temperature due to the formation of large sized and vigorous upwellings that in some cases tend to cluster, thus forming superplumes. If a very large thermal conductivity at the core-mantle boundary is assumed ($ksimeq 20$ W/mK) we obtain a quasi steady-state regime characterized by large and stable plumes with long lifetimes. The combination of strongly depth-dependent expansivity and conductivity is a viable mechanism for the formation of long-wavelength, long-lived thermal anomalies in the deep mantle, even if a low-viscosity ppv atop the core-mantle boundary is included.
Toussaint, R., Niebling, M., Flekkoy, E.G., Schmittbuhl, J., Vinningland, J., Maloy, K., Schelstraete, B., Johnsen, O., Lindner, A., Clement, E., Chevalier, C. & Koval, G. Sedimentation and aerofracture: Sedimentation instability and fracturing/channeling transitions in mixed grain/fluid flows - impact of viscosity and compressibility Session 2 Monday 7  oral file 
Abstract: We present systems where we inject a fluid at high pressure in a porous material saturated with the same fluid, or where gravity overpressurizes the fluid in interticial pore space during granular deformation. This fluid is either a highly compressible gas (air), or an almost incompressible and viscous fluid (oil or water). We compare both situations. These porous materials are designed as analogs to real rocks, but their cohesion are tuned so that the hydrofracture process can be followed optically in the lab. The fluid is injected on the side of the material. At sufficiently high overpressures, the mobilization of grains is observed, and the formation of hydrofracture fingering patterns and transition to thin branching fractures is followed and analyzed quantitatively. We study the pattern formation in the decompaction process starting from free boundaries, and the formation of channels and fractures around injection chambers or gravitationally pressurized chambers. Thresholds of pressures are determined for the formation of these preferential paths. The geometry of these channels, their fractal dimension and other characteristics, are extracted from experiments and simulations. The two situations with air or oil are compared together. Many similarities are observed about shape selections and dynamics, when time is rescaled with the viscosity of the fluid. We also investigate systematically the change of intersticial fluid, in compressibility and in viscosity, over a sedimentation Rayleigh/Taylor instability. The adaptation of hybrid granular fluid modelling algorithms to tackle the different types of flows will be presented. In practice, these problems are relevant for important aspects in the formation of increased permeability networks as seen in nature and industry: e.g., in active hydrofracture in boreholes, piping/internal erosion in soils and dams, sand production in oil or water wells, and wormholes in oil sands. It is also important to understand the formation of channels, and the behavior of confined gouges when overpressured fluids are mobilized in seismic sources. Indeed, the formation of preferential paths in this situation can severely affect the fluid and heat transport properties in this situations, and thus affect the pore pressurization effects. Eventually, the impact of lubrication effects in avalanches with long runoff, is also urging for models, as tested on the simple flows presented here.
Tyasto, M., Danilova, O., Ptitsyna, N., Sdobnov, V. & Villoresi, G. Evaluation the uncertainty of the Earth's Magnetospheric Magnetic Field Models Session 7 Thursday 10  poster file 
Abstract: As the magnetospheric dynamics depends on both the solar wind conditions and the previous history of the magnetosphere, and the currently available models use only the present values as input, they are not very reliable in predicting the magnetospheric state. Moreover, the magnetic field (MF) at any point of near-Earth space is the sum of fields from various current sources as currents of internal sources, on the magnetospheric boundary, ring currents, currents in high-latitude region and in the magnetospheric tail. However, the contribution of each component cannot be accurately estimated, because this task requires solving an inverse problem, which has multiple solutions. Therefore, all models describing magnetospheric MF by means of these current systems are approximate models. Even the level of the uncertainties cannot be correctly defined, since the real contribution from any of the above sources is not known. We proposed a method for validating magnetospheric MF models by means of cosmic ray data. Cosmic ray intensity distribution on the Earth surface contains information on the space distribution of MF encountered in the magnetosphere through which charged particles propagate. We have analyzed feasibility and limitations of cosmic ray data to be a tool for estimation of magnetospheric MF models. Our approach is based on the fact, that variations of cosmic rays are related to changes in geomagnetic cutoff rigidities. We have compared the cutoff rigidity changes obtained by the trajectory tracing method in the model MF with those obtained on the base of experimental cosmic ray data. The obtained results have shown that cosmic ray data can be successfully used for validation of models in presenting the dynamics of magnetospheric MF for mid latitudes.
Ueno, G., Higuchi, T., Kagimoto, T. & Hirose, N. Maximum likelihood estimation of error covariances in ensemble-based filters and its application to a coupled atmosphere-ocean model Session 8 Thursday 10  poster file 
Abstract: We propose a method for estimating optimal error covariances in the context of sequential assimilation, including the case where both the system equation and the observation equation are nonlinear. When the system equation is nonlinear, ensemble-based filtering methods such as the ensemble Kalman filter (EnKF) are widely used to deal directly with the nonlinearity. The present approach for covariance optimization is a maximum likelihood estimation carried out by approximating the likelihood with the ensemble mean. Specifically, the likelihood is approximated as the sample mean of the likelihood of each member of the ensemble. To evaluate the sampling error of the proposed ensemble-approximated likelihood, we construct a method for examining the statistical significance using the bootstrap method without extra ensemble computation. We apply the proposed methods to an EnKF experiment where TOPEX/POSEIDON altimetry observations are assimilated into an intermediate coupled model, which is nonlinear, and estimate the optimal parameters that specify the covariances of the system noise and observation noise.
Valentini, L. Numerical Study of the role of Korteweg Stress in Magma Dynamics Session 2 Thursday 8  poster file 
Abstract: The Korteweg Stress theory predicts the existence of capillary stresses arising at the interface between two miscible liquids characterised by a sharp compositional gradient. Such stresses may result in a series of effects that include interface relaxation and breakup, in analogy to surface tension-induced instabilities in immiscible liquids. In the last decade both computer models and fluid dynamical experiments have been devised, in various areas of research, in order to investigate Korteweg stresses and provide both a deeper theoretical insight and unequivocal evidence of their effects. However, only recently attention has been drawn on the possible role of Korteweg Stress in geophysical systems. In this study, a series of numerical experiments have been performed with the aim of ascertaining whether magma dynamics may be affected by the action of Korteweg stresses. In particular, attention was focussed on the simulation of the development of micro-scale emulsion textures in systems comprising two miscible magmas. Such textures have been observed, in the literature, in both natural and experimental samples. The results show that, depending on the physico-chemical properties of the interacting magmas (e.g. viscosity, diffusivity), surface minimization can be induced by processes of interface relaxation, breakup and coalescence, potentially leading to the formation of emulsion-like textures. Furthermore, a parametric study is performed with the aim of investigating the extent of Korteweg stress-driven instabilities in the presence of other processes related to buoyancy, momentum diffusivity, etc. New dimensionless numbers are defined, which relate Korteweg stresses to the other forces acting on the system, and may be used to infer at which scales and for which combinations of melts Korteweg stresses may affect the evolution of a magmatic system.
Vallis, G.K. The maintenance of stratification in the ocean and atmosphere: from conveyor belts to geostrophic turbulence. Session 3 Tuesday 8  oral file 
Abstract: We will discuss the maintenance of stratification in the ocean, comparing and contrasting with the maintenance of stratification in the atmosphere. In the upper ocean stratification is maintained by a combination of wind-driving (producing the ventilated thermocline) and an advective-diffusive processes, producing an internal thermocline. These processes are now reasonably well understood, and are also have no real atmospheric analog. The deep ocean is less well understood, and its stratification seems to be maintained, at least in part, by the action of baroclinic eddies originating in the Southern Ocean. We will present an analytic theory for this along with some numerical simulations, and discuss whether the role of the eddies is comparable to the role of baroclinic eddies in the atmosphere.
Vanadit-Ellis, W. US Army Centrifuge: Modelling Effects and Consequences of Earth Dynamics Session 7   withdrawn file 
Abstract: Large models can be built and tested at normal gravity (1-g), field tests can be conducted, but these are expensive, difficult to perform and most importantly they are often inaccurate. Predictions from 1-g model are also difficult if the materials are time, stress, and/or strain state dependent. However, if a small-scale model could be subjected to an acceleration field of magnitude many times earth's gravity, self-weight stresses and other gravity dependent processes can be simulated. Acceleration above normal gravity is achieved through the use of a centrifuge. While gravity exerts a dominant effect in all geotechnical situations other forces and loading can be imposed by mechanical devices on board a centrifuge. Dynamic problems such as earthquake can be simulated using a servo-hydraulic shaker. The paper will discuss different approaches of centrifuge modelling for quantifying the effects and consequences caused by earthquakes, floods and contaminant migration.
Vassalli, M., O'Brien, G.S., Bean, C.J., Lokmer, I. & Saccorotti, G. Importance of structural and rheological complexity on ground deformation inversion: a numerical study Session 4 Tuesday 8  poster file 
Abstract: Understanding the source of ground deformation at active volcanoes is of primary importance for a better assessment of the volcanic hazards. The inversion of recorded data is a powerful tool to infer features of the source responsible for such deformation. It is a computationally complicated and costly procedure due to the complexity of the rock domain and the number of unknown variables such as the source position, shape and extent. As a consequence, inversion methods are normally based on models that simplify the source geometry and/or rock properties. We have examined the degree to which the complexities of the rock matrix may be neglected when interpreting data associated with awakening episodes. The ground deformation associated with a volume change at depth has been simulated using a numerical code based on a discrete elastic lattice method (O'Brien and Bean, 2004) which accounts for the presence of complex topography, rock heterogeneities and fractures/faults. The deformation due to the same source but different rock models, from homogeneous to heterogeneous with fault discontinuities and topography have been computed. These synthetic data have been then inverted through the commonly used model by Davis (1986) that accounts for a point pressurized cavity of ellipsoidal shape, arbitrarily oriented in an homogeneous half-space. The results show that the inversion of ground deformation data is very sensitive and that the assumption of simplified models (i.e. ignoring rock heterogeneity) can lead to a wrong interpretation of the ground deformation itself, in particular in term of the source mechanism and magnitude. The present study highlights that the incorporation of realistic rock models is crucial for determining the correct source of volcanic ground deformation. A possible solution could be the use of Green's function calculated for the specific site under-investigation but it involves a large computational cost.
Vicari, A., Bilotta, G., Bonfiglio, S., Cappello, A., Ganci, G., Herault, A., Rustico, E. & Del Negro, C. GPU-based models to perform numerical simulations of lava-flow dynamics Session 1 Thursday 10  poster file 
Abstract: Lava flows represent a problem particularly challenging for physically based modeling because of the mechanical and thermal features of lava change over time. In order to generate complex trajectories due to the interactions between lava flows and the underlying topography, we need to model the main mechanical features of lava and the way they evolve over time depending on temperature. That challenge has inspired the INGV-CT to develop two different models to study lava-flow dynamics. The first model, named MAGFLOW, represents the top of the evolution of cell-based models for lava-flow simulations. It is based on a cellular automata structure in which the evolution function is a steady state solution of the Navier-Stokes equations and heat transfer (due to radiative losses) and solidification effects are modeled via a temperature dependent viscosity. The model was validated comparing simulated lava flows with the real ones occurred at Mt Etna during the 2001, 2004, 2006, and 2008 eruptions, showing that the code works properly fitting well-constrained eruption data sets. The second numerical model, based on Smoothed Particle Hydrodynamics (SPH) approach, was recently developed by INGV-CT. The SPH method allows us to tackle the different aspects of direct simulation of a lava flow on a real topography, including the numerical resolution of the Navier-Stokes equations coupled with the energy equation with a non-linear rheology (four different rheologies have been implemented: Newton, Bingham, power-law and Herschel-Bulkley) and taking into consideration phase changes. The model allows us to obtain the distribution of solid parts and velocity, temperature and viscosity values across the whole lava flow. For both models the computational cost is however considerable, especially for SPH model, when we need to be able to handle long simulations with millions of particles. To this purpose, we implemented both models on modern graphic processor unit (GPU) by using a NVIDIA CUDA architecture, obtaining a speed-ups respectively of over 40x for MAGFLOW code and 120x for SPH code.
Vilotte, J., Festa, G. & Raous, M. Dynamic rupture along interfaces including damage and friction: initiation, propagation and radiation Session 4 Monday 7  oral file 
Abstract: Seismic rupture of large earthquakes initiates and propagates along pre-existing faults that have a complex internal structure (fault zone). The mechanical and numerical modeling of dynamic earthquake rupture is of great importance for seismic engineering and seismic risk management.
The propagation and the radiation of a seismic rupture has long been considered in seismology as a friction dominated process, and formulated as a propagating shear crack problem under the assumption of a Barenblatt-type surface energy. Slip weakening and rate-and-state friction laws are today commonly used. Both contain intrinsic length scales leading to a finite fracture energy, and a characteristic cohesive length-scale, insuring finite stress and slip velocity. As the quality and the density of the seismological and geodetic observations are continuously improving, new mechanical and numerical models are developed to study realistic dynamic rupture propagation in heterogeneous faults and extract new information from the recorded radiated waves and energy.
We first provide a short description of recent development of numerical methods based on non smooth contact mechanics and high-order variational approximation together with interface constitutive laws that include the interactions between volumetric damage breakdown processes and frictional dissipation.
Important scaling issues remain in the modeling of rupture dynamics and radiation. We first investigate the radiation and directivity of simple rupture dynamical models along complex fault geometries, including rupture initiation, rupture propagation and rupture arrest phases, comparing classical slip-weakening and interface laws including damage and friction. We then investigate effects related to off-fault dissipation in terms of volumetric dynamical damage, and its potential implication for the scaling of the radiative energy, for different fault geometries.
Finally, we discuss some open issues regarding the multi-scale modelling of earthquake rupture dynamics, and high frequency generation, when keeping a formulation based on surface energy.
Viskari, T., Pekka, K., Tatu, A. & Heikki, J. Simultaneous retrieval of aerosol variables through data-assimilation Session 8 Thursday 10  poster file 
Abstract: Detailed knowledge of the atmospheric aerosol size distribution and chemical composition is required for better understanding of the effects of aerosols on climate, air quality and health. In situ measurements provide indirect information of aerosols, such as electric mobility. Mathematical inversion techniques are successful in aerosol size distribution retrieval, but the results are limited to the parameters in the particle size range resolved by the instrument. Combining different measurement types or measurements distributed in time has proven difficult for inversion techniques. Simultaneous retrievel of chemical composition and size distribution through inversion remains a significant scientific question.
An alternative retrieval approach based on data assimilation is introduced here. A sectional evolution model of aerosol microphysics and an array of measurement data is intended to be assimilated for the retrieval of aerosol size and composition parameters. In contrast to purely mathematical inversion techniques, data assimilation enables the use of aerosol dynamics in the retrieval process, thus effectively constraining the set of possible solutions in a physically consistent way. Essentially, simultaneous use of multiple instruments, temporally distributed measurement data and an evolution model as a strong constraint for the retrieval solution provides enough information for fixing all degrees of freedom of the evolution model, even some of those not directly observed.
The modelling component of the proposed research is based on the University of Helsinki Multi-component Aerosol model. The measured variables include the particle size distribution in different particle sizes, particle mass distribution, chemical composition, light absorption, ambient vapour concentrations and meteorological parameters.
The method has been tested with synthethic size distributions and idealized conditions, with the method succesfully determining the correct initial state even when the initial estimation error was significant. Additionally, assumptions required for the application of data assimilation were found to be valid for the aerosol dynamical models and observations except for certain conditions concerning the size distribution of particles smaller than $10 nm$. In these cases the evolution of the number concentration error covariance becomes problematic.
The next step of the research is to apply the method with real observations. This is currently underway.
Weatherley, D.K. On the role of boundary conditions in particle-based numerical simulations of brittle failure Session 4 Tuesday 8  poster file 
Abstract: Brittle failure is a complex phenomenon involving nonlinear interactions within brittle media and trans-scale coalescence of cracks. The macroscopic constitutive behaviour of brittle media is a combination of both the intact strength of the constituent minerals and the frictional strength of pre-existing discontinuities within the medium. Simulating brittle failure using continuum-based numerical methods is challenging both from a mathematical standpoint, to define appropriate nonlinear constitutive laws, and a computational standpoint, requiring adaptive re-meshing to properly characterise the propagation of cracks within the numerical domain in many instances.
A popular alternative numerical method is the Discrete Element Method (DEM). The DEM represents a brittle medium as an assembly of spherical particles connected via brittle-elastic springs (called bonds). At each timestep, the net force acting on each particle is computed and Newton's Laws are integrated to update particle positions and velocities. External forces are applied to the assembly via elastic interactions between boundary particles and movable walls. Although conceptually simple, the DEM has proven popular for numerical studies on brittle failure since it does not require complex constitutive laws and expensive re-meshing algorithms to mimic aspects of the phenomenology of brittle failure.
In a series of high-resolution DEM simulations involving uniaxial compressive failure of cylindrical samples, the patterns of broken bonds have been closely examined. Using typical methods to apply stress to the samples, broken bonds are distributed uniformly throughout the sample and only during late stages of failure do wide shear bands develop which cause the sample to fragment. To simulate localised crack coalescence and formation of realistic fracture surfaces, a new type of boundary condition has been developed. Rather than monotonically increasing the stress applied to the sample, the loading rate is reduced as bonds break and the compressive strength of the sample begins to degrade. The extreme sensitivity of brittle failure patterns in DEM simulations to the boundary conditions is a surprising result with broad implications for the application of the DEM to understand brittle failure of heterogeneous media.
Wicht, J. Towards Realistic Planetary Dynamo Simulations Session 5 Friday 11  oral file 
Abstract: The last years have witnessed an impressive growth in the number and quality of numerical dynamo simulations. These models successfully describe many aspects of the geomagnetic field and also set out to explain the distinctly different fields of other planets. The success is somewhat surprising since numerical limitation force dynamo modelers to run their computations at unrealistic parameters. In particular the Ekman number, a measure for the relative importance of viscous diffusion to the Coriolis force, is many orders of magnitude too large. After giving a brief introduction into the basics of modern dynamo simulations we discuss the fundamental dynamo regimes and address the question how well the modern models reproduce the geomagnetic field. First-level properties like the dipole dominance, realistic magnetic field strength, convective flow vigor, and an Earth-like reversal behavior are already captured by larger Ekman number simulations. However, low Ekman numbers are required for modeling torsional oscillations which are thought to be an important part of the decadal geomagnetic field variations. Moreover, only low Ekman number models seem to retain the huge dipole dominance of the geomagnetic field once the Rayleigh number has been increased to values where field reversals happen. These cases also seem to resemble the low-latitude field found at Earth's core-mantle boundary more closely than larger Ekman numbers cases.
Yuen, D.A., Wright, G.B., Sanchez, D.A. & Barnett, Jr., G.A. The Coming Role of GPU in Computational Geosciences
Session 9 Monday 7  oral file 
Abstract: For the last couple of years GPU has been catching the attention of many interested parties in computational world, because of its alluring speed and inherent potential growth over CPU. With the arrival of CUDA in 2007, which has really facilitated the GPU programming of scientific software, GPU has really expanded into diverse scientific areas, ranging from chemistry to astrophysics. Geosciences proved to be no exception to this trend and we are seeing already some incursions by GPU into seismic wave propagation, flow in fractured porous media and tsunami wave studies. Today we will discuss our experience with GPU in thermal convection. We will review our venture into using a translator of MATLAB codes into CUDA via JACKETS, which is a product of Accelereyes. This effort has allowed us to study 3D thermal convection at the Infinite Prandtl limit at a resolution of 400x400x200 finite-difference points with a second-order accurate in space and third-order in time. We have also implemented a CUDA code from scratch of a 2-D thermal convection code with the same spatial-temporal accuracy as the 3D code. This code can attain a speed on a single Tesla C-1060 GPU on the order of O(0.1) microsec/ time step grid point at the limit of asymptotically large number of grid points. We will present results at large Rayleigh numbers, in excess of a hundred million and compare them with CPU results.
Zurita-Gotor, P. & Vallis, G.K. Circulation sensitivity to heating in a simple model of baroclinic turbulence Session 3 Tuesday 8  oral file 
Abstract: In this presentation we discuss the sensitivity of the equilibrium state of an idealized two-level primitive equation model on the form and strength of the heating. We compare two different versions of the model, forced using two different heating formulations: a standard model forced by Newtonian cooling and a ``prescribed heating model'', in which the net heating reaching the lower level is prescribed. To interpret the results, it is useful to distinguish between two facets of the heating: the differential heating, defined as the meridional heat transport from low to high latitudes in the model, and the vertical destabilization, or heat transport from the lower to the upper troposphere. While the differential heating is always internally determined (in both heating formulations), the vertical destabilization is internally determined in the Newtonian cooling model but prescribed from the outset in the model forced from below.
An important constraint for understanding the results is that, when the eddies are quasi-adiabatic, the eddy mixing slope should scale as the isentropic slope. This implies that the isentropic slope is ultimately constrained by the structure of the heating, and in particular by the ratio between vertical destabilization and differential heating. Consistent with this, the criticality always decreases (i.e., the isentropic slope always flattens) with increasing differential heating in the model forced from below. In contrast, the criticality may either increase or decrease with increasing differential heating in the Newtonian forcing model, depending on whether the vertical destabilization -internally determined now- increases faster or slower than the differential heating. When they change proportionally, the criticality is relatively insensitive to the external heating.

Created by JabRef on 15/04/2010.