RESEARCH INTERESTS
Alik Ismail-Zadeh
What
do I do in science?
A
principal objective of my research is to understand interplay between crustal
and mantle dynamics and surface processes by means of theoretical and
experimental studies. My approach is three-pronged. 1: I study geophysical
phenomena on the basis of available observations. 2: I develop mathematical and
numerical methodologies to model geodynamic processes (such as deformations and
stresses in and around the descending lithosphere) or fluid flow (such as lava
flow); thermal convection in the mantle; sedimentary basin evolution and salt diapirism; and gravitational, thermal and buckling
instabilities of geological structures. 3: I test model predictions against
observations and explain specific features of the processes and phenomena.
Why
should the general public be interested in what I do?
When
nature unleashes itself (resulting in earthquakes, volcanic eruptions, lava
flows, landslides or tsunami), it strikes where and when it decides. Extreme
events strike without warning and cause devastation in loss of human life and
environmental damage. There is no doubt that man will never be able to prevent
these occurrences entirely. However, geoscientists can gain a better
understanding of the complex mechanisms behind the geohazards, dynamics of the
lithosphere and its surface manifestations. Convolving the knowledge of
physical processes with the knowledge of vulnerability and exposure we can
construct a map of risks alerting policymakers and society on disaster risk
reduction.
Why
does it interest me?
As a geophysicist, I interest it to
understand a complex dynamic behavior of the Earth. As an applied
mathematician, I interest it to bring formalism to Earth sciences. As a
habitant of our planet, I interest it to know how to cope with disasters.
CURRENT to PAST RESEARCH
Natural hazards, climate change and drivers of disasters
Many nations face challenges in assessing,
understanding, and responding to the time-dependent nature of disaster risk. Changes
in the intensity of occurrences of extreme events coupled with changes in vulnerability and
exposure alter the impacts of natural hazards on society in mostly negative ways. An
interrelationship between natural hazard (NH), climate change (CC), vulnerability (V),
exposure (E), and decisionmaking (DM) is considered. While NHs trigger disasters and CC is
likely to intensify occurrences of disasters, V and E present major drivers of disasters.
Informed DM on disaster risk reduction should be based on scientific evidence from NH
and CC, knowledge of V and E, and relevant options for actions on preventive disaster
measures as a part of preparedness and public awareness. Published in
Natural Hazards, 111,
2147–2154, 2022.
Lava dome morphology and viscosity inferred from data-driven numerical
modeling of dome growth at Volcan de Colima, Mexico during 2007-2009
Magma extrusion, lava dome growth, collapse of
domes, and associated pyroclastic flow hazards are among important volcanological studies.
We analyze the influence of the magma viscosity and discharge rates on the lava dome morphology
at Volcan de Colima in Mexico during a long dome-building episode lasting from early 2007 to fall
2009 without explosive dome destruction. Camera images of the lava dome growth together with
recorded volumes of the erupted lava have been used to constrain numerical modeling and
hence to match the history of the dome growth by nudging model forecasts to observations.
Our viscosity model incorporates crystal growth kinetics and depends on the characteristic time
of crystal content growth (or CCGT) and the crystal-free magma viscosity. Initially, we analyze
how this viscosity, CCGT, and the rate of lava extrusion influence the morphology of the growing
dome. Several model scenarios of lava dome growth are then considered depending on the crater
geometry, the conduit location, the effective viscosity of dome carapace, and the extrusion rates.
These rates are determined either empirically by optimizing the fit between the morphological
shape of modeled domes and that of the observed dome or from the recorded lava dome volumes.
The maximum height of the modeled lava dome and its horizontal extent are in a good agreement
with observations in the case of the empirically-derived extrusion rates. It is shown that the
topography of the crater at Volcan de Colima is likely to be inclined toward the west. The viscosity
of the modeled lava dome (~10^12 Pa s) is in a good agreement with the effective viscosity
estimated experimentally from lavas of Volcan de Colima. Due to the interplay between the lava
extrusion and the gravity forces, the dome reaches a height threshold, and after that a horizontal
gravity spreading starts to play an essential role in the lava dome evolution. The model forecasts
that the dome carapace of higher viscosity (~10^14 Pa s) influences the dome growth and its
morphology during long dome-building episodes by retarding horizontal advancement and
developing steep-sided eastern edge of the dome at the volcano. The developed model can be
used in assessments of future effusive eruptions and lava dome growth at Volcan de Colima
or elsewhere. History matching modeling of lava dome growth sheds a light on dynamic processes
inside the dome and may assist in assessing stress state in the dome carapace and in forecasting
the dome failures. The results are published in
Frontiers in Earth Science, 9,
735914, 2021.
Inferring lava viscosity from lava dome morphology
The viscosity of magma nonlinearly depends
on the volume fraction of crystals and temperature. We estimate magma viscosity based on
a comparison of observed and simulated morphological forms of lava domes.
We consider a two-dimensional axisymmetric model of magma extrusion on the surface
and lava dome evolution, and assume that the lava viscosity depends only on the
volume fraction of crystals. The crystallization is associated with a growth of the liquidus
temperature due to the volatile loss from the magma, and it is determined by the characteristic
time of crystal content growth (CCGT) and the discharge rate. Lava domes are modeled using
a finite-volume method implemented in Ansys Fluent software for various CCGTs and volcanic
vent sizes. For a selected eruption duration a set of morphological shapes of domes (shapes
of the interface between lava dome and air) is obtained. Lava dome shapes modeled this way
are compared with the observed shape of the lava dome (synthesized in the study by a random
modification of one of the calculated shapes). To estimate magma viscosity, the deviation
between the observed dome shape and the simulated dome shapes is assessed by three
functionals: the symmetric difference, the peak signal-to-noise ratio, and the structural similarity
index measure. These functionals are often used in the computer vision and in image processing
(AI methods). Although each functional allows to determine the best fit between the modeled
and observed shapes of lava dome, the functional based on the structural similarity index
measure performs it better. The viscosity of the observed dome can be then approximated by
the viscosity of the modeled dome, which shape fits best the shape of the observed dome.
This approach can be extended to three-dimensional case studies to restore the conditions
of natural lava dome growth. The results are published in
J. Volcanol. Seismol., 15(3),
1159–168, 2021.
Lava dome morphology inferred from numerical modelling
Lava domes form when highly viscous
magmas erupt on the surface. Several types of lava dome morphology can be
distinguished depending on the flow rate and the rheology of magma. We develop
a two-dimensional axisymmetric model of magma extrusion on the surface and lava dome
evolution and analyse the dome morphology using a finite-volume method implemented
in Ansys Fluent software. The magma/lava viscosity depends on the volume fraction of
crystals and temperature. We show that the morphology of domes is influenced by two
parameters: the characteristic time of crystal content growth (CCGT) and the discharge
rate (DR). At smaller values of the CCGTs, that is, at rapid lava crystallisation,
obelisk-shaped structures develop at low DRs and pancake-shaped structures at high DRs;
at longer CCGTs, lava domes feature lobe-shaped to pancake-shaped structures.
A thick carapace of about 70% crystal content evolves at smaller CCGTs. We demonstrate
that cooling does not play the essential role during a lava dome emplacement, because
the thermal thickness of the evolving carapace remains small in comparison with the
dome’s height. A transition from the endogenic to exogenic regime of the lava dome
growth occurs after a rapid increase in the DR. A strain-rate-dependent lava viscosity
leads to a more confined dome, but the influence of this viscosity on the dome
morphology is not well pronounced. The model results can be used in assessments
of the rates of magma extrusion, the lava viscosity, and the morphology of active lava
domes. The results are published in
Geophys. J. Int., 223,
1597–1609,
2020.
Examples of morphological shapes of lava domes: on the left are images of domes
of the Soufriere Hills volcanic complex (according to Watts et al., 2002); on the right are
lava domes according to the studied models (colors indicate the distribution of viscosity
in lava domes; the most viscous part of lava in red, the least viscous part in blue).
Nonlinear dynamics of crustal blocks and faults and earthquake occurrences
in the Transcaucasian region
The Caucasus is a part of ongoing collision
of the Arabian and Eurasian plates, where moderate to strong earthquakes caused
significant losses of lives and livelihood in the past. To better understand seismic hazard
in the region, we develop a model of block-and-fault dynamics for Transcaucasia, the
largest part of the Caucasus to the south of the Greater Caucasus Mountains, to simulate
regional earthquakes. The model structure is developed by employing the results of the
morphostructural analysis to delineate crustal blocks and the geodetic observations on
crustal movements in the region. The model incorporates a nonlinear rate-dependent slip
of the faults separating the blocks. A set of numerical experiments has been performed
to address the following questions: (i) where strong earthquakes occur and what their
reoccurrence time is; (ii) how rigid crustal blocks react to the Arabian plate push and
to movements of the ductile part of the crust in Transcaucasia; and (iii) whether the
fault slip rates and the block displacements in the model correlate with observed
GPS-velocities. The model results confirm that the contemporary crustal dynamics
and seismicity pattern in Transcaucasia are determined by the north-northeastern
motion of the Arabian plate relative to Eurasia and by the movement of the ductile crust
underlying the rigid crustal blocks. Variations in the rheological properties of the fault zones
and/or of the ductile crust influence the pattern of seismicity. The number and maximum
magnitude of synthetic earthquakes change with the variations in the movements of the
crustal blocks and in the rheological properties of the lower crust and the fault zones.
The model results can be used in comprehensive seismic hazard assessment of the
Caucasus region based on instrumentally observed, historical and synthetic seismicity.
The results are published in Physics
of the Earth and Planetary Interiors, 297, 106320, 2019.
Distribution of synthetic seismicity (with earthquake magnitudes greater than 6)
imposed on the map of observed seismicity: blue stars mark earthquakes of magnitude Ì7+,
blue circles Ì6.5+ events, and light blue dots Ì6+ events.
Forging a paradigm shift in disaster science
Despite major advancements in
knowledge on disaster risks and disasters caused by natural hazards, the number
and severity of disasters are increasing. Convolving natural, engineering, social
and behavioral sciences and practices with policymaking should significantly reduce
disaster risks caused by natural hazards. To this end, a fundamental change in
scientific approaches to disaster risk reduction is needed by shifting the current
emphasis on individual hazard and risk assessment dominant in the geoscientific
community to a transdisciplinary system analysis with action-oriented research on
disaster risk reduction co-produced with multiple stakeholders, including
policymakers. This paradigm shift will allow for acquisition of policy-relevant
knowledge and its immediate application to evidence-based policy and decision
making for disaster risk reduction. The need for the paradigm shift is more critical now
than ever before because of the increasing vulnerability and exposure of society to
disaster risk and the need for cross-cutting actions in policy and practice related to
climate change and sustainability. The results are published in Natural Hazards, 86, 969–988, 2017.
Reconstruction of thermal
and dynamic characteristics of volcanic lava from surface thermal measurements
We study a model of volcanic lava
flow to determine its thermal and dynamic characteristics from thermal
measurements of the lava at its surface. Mathematically this problem is reduced
to solving an inverse boundary problem. Namely, using known conditions at one
part of the model boundary we determine the missing condition at the remaining
part of the boundary. We develop a numerical approach to the mathematical
problem in the case of steady-state flow. Assuming that the temperature and the
heat flow are prescribed at the upper surface of the model domain, we determine
the flow characteristics in the entire model domain using a variational
(adjoint) method. We have performed computations of model examples and showed
that in the case of smooth input data the lava temperature and the flow
velocity can be reconstructed with a high accuracy. As expected, the noise
imposed on the smooth input data results in a less accurate solution, but still
acceptable below some noise level. Also we analyse the influence of
optimization methods on the solution convergence rate. The proposed method for
reconstruction of physical parameters of volcanic lava can also be applied to
other problems in geophysical fluid flows. The first results are published in Geophys. J. Int., 205, 1767–1779, 2016.
|
|
|
Numerical modelling of
lava flow with a ruptured crust
Although
volcanic lava flows do not significantly affect the life of people, its hazard
is not negligible as hot lava kills vegetation, destroys infrastructure, and
may trigger a flood due to melting of snow/ice. The lava flow hazard can be
reduced if the flow patterns are known, and the complexity of the flow with a
debris is analyzed to assist in disaster risk mitigation. In this paper we
develop three-dimensional numerical models of a gravitational flow of
multi-phase fluid with rafts (mimicking rigid lava-crust fragments) on a
horizontal and topographic surfaces to explore the dynamics and the interaction
of lava flows. We have obtained various flow patterns and spatial distribution
of rafts depending on conditions at the surface of fluid spreading, obstacles
on the way of a fluid flow, raft landing scenarios, and the size of rafts.
Furthermore, we analyze two numerical models related to specific lava flows: (i) a model of fluid flow with rafts inside an inclined
channel, and (ii) a model of fluid flow from a single vent on an artificial
topography, when the fluid density, its viscosity, and the effusion rate vary
with time. Although the studied models do not account for lava solidification,
crust formation, and its rupture, the results of the modeling may be used for
understanding of flows with breccias before a
significant lava cooling. The first results are published in J. Geodyn., 97, 31-41, 2016.
|
|
|
On the use of multiple-site estimations in probabilistic seismic hazard assessment
We
analyze specific features of multiple-site (MS) probabilistic seismic hazard
assessment (PSHA), i.e. the annual rate of ground-motion level exceedance in at
least one site of several sites of interest located within in an area or along
a linear extended object. The relation between MS hazard estimates and strong
ground-motion data obtained during large earthquakes is discussed in the cases
of the 1999 Chi-Chi MW 7.6 earthquake and the 2008 Wenchuan MW
7.9 earthquake. The strong-motion records obtained in the epicentral zones may
be considered as an example of the ground motion exceeding the design level
estimated using the conventional point-wise (PW) PSHA. We show that the
MS-PSHA, when being performed for the standard return period 475 years,
provides reasonable estimations of the intensity level that may occur during
the earthquakes, parameters of which are close to the parameters of events with
maximum possible magnitude accepted in PSHA for the regions. The MS-PSHA
efficiency is discussed with respect to (a) evaluation of the performance of
probabilistic seismic hazard maps and (b) application of the MS hazard
estimates as a basis for design loads. Based on the results of this work, we
propose a multi-level approach to PSHA considering fixed reference probability
of exceedance (e.g. 10% in 50 years): (i) a standard
PW-PSHA to be performed in a seismic-prone region (first level), and (ii) this
analysis should be supplemented by a MS-PSHA for urban and industrial areas, or
zones of a particular economic and social importance (second level). The
results are published in Bull. Seismol. Soc. Am., 106(5),
2233–2243, 2016.
Relation between the multiple-site PGAMLT
and point-wise PGAPNT
hazard estimations for the 475-yr (a) and 2475-yr (b) return periods. The
dashed line shows ratio between the maximum PGA recorded during the Wenchuan
earthquake (about 810 cm s-2, geometric mean of two horizontal
components) and the estimated design ground motion PGA475 ~ 300 cm s-2
(a) and PGA2475 ~ 600 cm s-2 (b).
1: area 600 km2; 2: 400 km2; 3: 100 km2.
Seismic hazard from instrumentally recorded, historical and simulated earthquakes
We
study regional seismic hazard accounting for observed (instrumentally recorded
and historic) earthquakes, as well as for seismic events simulated for a
significantly longer period of time than that of observations. We use a
probabilistic seismic hazards analysis (PSHA) with elements of scenario-based
analysis for the Tibet-Himalayan region. The large magnitude synthetic events,
which are consistent with the geophysical and geodetic data, together with the
observed earthquakes are employed for the Monte-Carlo PSHA. Earthquake
scenarios for hazard assessment are generated stochastically to sample the
magnitude and spatial distribution of seismicity, as well as the distribution
of ground motion for each seismic event. The peak ground acceleration values,
which are estimated for the return period of 475 yr, show that the hazard level
associated with large events in the Tibet-Himalayan region significantly
increases if the long record of simulated seismicity is considered in the PSHA.
The magnitude and the source location of the 2008 Wenchuan M=7.9 earthquake are
among the range of those described by the seismic source model accepted in our
analysis. We analyze the relationship between the ground motion data obtained
in the earthquake’s epicentral area and the obtained PSHA estimations using a
deaggregation technique. The proposed approach provides a better understanding
of ground shaking due to possible large-magnitude events and could be useful
for risk assessment, earthquake engineering purposes, and emergency planning.
The results were published in Tectonophysics,
657, 187-204, 2015.
Comparison of instrumental intensity (USGS shake map for the 2008
Wenchuan MW 7.9
earthquake) with the PGAs (in cm s-2; for rock site conditions and
the return period of 475 yrs) taken from (a) the national Seismic Ground Motion
Zonation Map GB 18306-2001, values in parentheses are given for the
medium-stiff soil; (b) the
GSHAP map; and (c) from the results
of this study.
Mantle/lithosphere
dynamics beneath the Japanese islands
Recent seismic tomography studies image a low velocity zone
(interpreted as a high temperature anomaly) in the mantle beneath the subducting Pacific plate near the Japanese islands at the
depth of about 400 km. This thermal feature is rather peculiar in terms of the conventional
view of mantle convection and subduction zones. Here we present a
dynamic restoration of the thermal state of the mantle beneath this region assimilating geophysical, geodetic,
and geological data up to 40 million years. We hypothesise
that the hot mantle upwelling beneath the Pacific plate partly penetrated
through the subducting plate into the mantle wedge
and generated two smaller hot upwellings, which
contributed to the rapid subsidence in the basins of the Japan Sea and to back-arc spreading. Another part of the hot
mantle migrated upward beneath the Pacific lithosphere, and the presently
observed hot anomaly is a remnant part of this mantle upwelling. The results
were published in Nat. Sci. Rep., 3, 1137,
2013.
Crustal
block-and-fault dynamics, fault slip rates and earthquake flow in the
Tibet-Himalayan region
The
Tibetan plateau and Himalayans have resulted from the continuous Indian and
Eurasian plate convergence following their initial collision at about 55
million years ago. Earthquakes in the region occur mainly in response to the
crustal motion and stress localization associated with this convergence. To
understand the basic features of the motion and seismicity in the
Tibet-Himalayan region, I develop a model of its block-and-fault dynamics. The
model structure is composed of twelve interacting upper crustal blocks
separated by fault planes. I develop several sets of numerical experiments
constrained by the regional seismic observations, geodetic and geological data.
Synthetic large events in the numerical experiments are clustered mainly on the
fault segments associated with the Himalayan Frontal Thrust as well as at some
internal faults of the Tibetan plateau. The clustering of earthquakes at some
fault is a consequence of dynamics of the regional fault system rather than
that of the fault only. I show that variations in the relationship of magnitude
to frequency of the events depend on changes in the motion of the upper crustal
blocks and on the rheological properties of the lower crust and fault zones. I
demonstrate that the predicted crustal motion in the region is characterized by
the north-northeastern movement of India toward Eurasia. Fluctuations in
rheological properties of fault plane zones and/or the lower crust influence
displacement rates of the crustal blocks and hence slip rates at the faults
separating the blocks. This can explain the discrepancies in estimates of slip
rates at major faults in the region (e.g., Altyn Tagh, Karakorum) over short and long time scales.
Observed and modeled large earthquakes in the Tibet-Himalayan region. (a) Seismicity from 1902 to 2000. White bold lines (model faults)
delineate the structural geological elements, and the white arrow indicates the
motion of India relative to Eurasia. Earthquake epicenters are marked by colored
(depending on the depth of the hypocenters) circles. (b)-(d) Simulated
earthquakes for three BAFD catalogs (modified after Ismail-Zadeh et al., 2007).
Quantitative
approach to thermo-mechanical restoration of salt-bearing sedimentary basins
(SBSBs).
The
research activity will address most important processes occurring in SBSBs:
thermal and structural evolutions of the basins. In applications they are
referred to as a maturity of hydrocarbons, their migration paths, deformations
of sedimentary layers, and stresses regime. I intend to study thermo-mechanical
restoration models for Pricaspian salt basin (East
European platform), which plays an important role in hydrocarbon explorations
and exploitations. The goal of the intended research is to understand the interplay
between geodynamic, geothermal, and tectonic processes involved in the
evolution of the SBSB. The specific objectives of the research are (i) to develop a 3D numerical approach to structural
restoration of SBSBs; (ii) to develop a 3D numerical approach to
thermo-mechanical restoration of SBSBs, and hence to determine temperatures
within the basins in the geological past based on their present-day
estimations; and (iii) to analyze the evolution of the Pricaspian
basin using thermo-mechanical restoration models and geological and geophysical
observations.
Inverse Problem of
the Thermal Convection
To
restore quantitatively the seismically detected mantle structures and
temperature field, we need a numerical tool for solving an inverse problem of
thermal convection at infinite Prandtl number. I develop a variational
approach to three-dimensional numerical restoration of thermoconvective
mantle flow. The approach is based on a search for the mantle temperature and
flow in the geological past by minimizing differences between mantle
temperature derived from seismic velocities (or their anomalies) and
temperature predicted by forward models of mantle flow. The mantle temperatures
and flow fields in the past so obtained could be employed as constraints on forward
models of mantle dynamics.
Salt Tectonics and
Basin Evolution
Numerical
studies of ductile deformations induced by salt movements have, until now, been
restricted to two-dimensional modeling (2D) of diapirism.
I study several 3D models of salt diapirism and
deformation of overlying sediments. We analyze a model of salt diapirs that develop from an initial random perturbation of
the interface between salt and its overburden and then restore the evolved salt
diapirs to their initial stages. An evolutionary
model of a 2D salt wall loaded by a 2D clinoform of
sediments predicts a decomposition of the salt wall into 3D diapiric structures
when the overburden of salt is supplied by 3D syn-kinematic
wedge of sediments. Also we study how lateral flow effects the evolution of
salt diapirs. The shape of a salt diapir
can be very different, if the rate of horizontal flow is much greater than the
initial rate of diapiric growth solely due to gravity.
Geodynamics, Stress
and Seismicity
Tectonic stress
in the lithospheric slab descending beneath the Vrancea (SE-Carpathians)
To
understand processes of stress generation in the descending slab, I analyze
tectonic stress in the slab by means of analytical and numerical modeling. I
have found that the maximum shear stress migrates from the upper plane of the
Benioff-Wadati seismic zone to its lower plane in
course of changes in slab dynamics from its active subduction through roll-back
movements to passive sinking solely due to gravity. It can explain a location
of hypocenters of Vrancea events at the side of slab adjacent to the Eastern
European platform. To understand processes of stress generation, I analyze the
tectonic stress induced by the sinking Vrancea slab employing 3D Eulerian FEM
and using (as input model data) seismic tomographic high-resolution images of
the Vrancea lithosphere, results of refraction seismic study, and other
geophysical and geological data. On the basis of experimental data on elastic
parameters and anelasticity I obtain initially a
model of mantle temperature from the P-wave velocity anomalies; crustal
temperature is derived from heat flow data. Based on temperature- and
pressure-dependent viscosity and density, modeled tectonic stress predicts
horizontal compression beneath the Vrancea region, which coincides with the
stress regime defined from fault-plane solutions for intermediate-depth
earthquakes. The stress reaches maximum at depths of 70 to 110 km and 130 to
180 km and decays below being in a good agreement with the observed seismicity.
Stress and
buoyancy-driven deformations of the lithosphere beneath central Apennines
The
juxtaposed contraction and extension, a long-standing geological enigma
recognized in different geodynamic frameworks world-wide, observed at the surface
and active nowadays, hence better studied, only in the Italian Peninsula, has
for a long time attracted the attention of geoscientists. Several models,
invoking mainly external forces, have been put forward to explain the close
association of these two end-member deformation mechanisms clearly observed by
geophysical, geodetic and geological investigations. These models appeal to
interactions along plate margins or at the base of the lithosphere such as
back-arc extension or shear tractions from mantle flow or to subduction
processes such as slab roll back, retreat or pull and detachment. On the basis
of seismic images of the lithosphere beneath the central Apennines and
numerical modeling of lithosphere dynamics and in-situ stress, I find that buoyancy
forces can explain solely the present-day stress distribution in the Apennine
lithosphere, complex lithospheric deformations, and unusual earthquake
distribution.
Dynamic
quantitative restoration of profiles across diapiric structures
The
backstripping method that is widely used in basin
analysis sometimes fails for salt-bearing basins, because the highly mobile and
buoyant salt deforms its sedimentary overburden. I developed a new numerical
approach for two-dimensional dynamic restoration of cross-sections through
successive earlier depositional stages. Our models interpreted basin profiles
as multiple ductile layers with various physical parameters (e.g., density,
viscosity, Young’s modulus). The evolution of salt structures was modeled
backwards in time by removing successively younger layers and restoring older
layers and any diapirs to the stage they were likely
to have been. The applicability of the technique was demonstrated by
reconstructions of upbuilt and downbuilt diapirs. I used the technique to restore a depth-converted
seismic cross-section through the south-eastern part of the Pricaspian
salt basin [5]. Mature salt diapirs in the section
are shown to have been downbuilt from a salt layer
with an initially uniform thickness, as a result of differential sedimentary
loading until the end of the Triassic before one of the diapirs
was buried and actively upbuilt. Our approach is well suited for restoration of
cross-sections with ductile overburdens and now is being developed to
three-dimensional restorations and other rheologies.
Gravitational and
buckling instabilities of rheologically stratified geostructures
I
investigated the gravitational and buckling instabilities of a structure
consisting of a buoyant layer of viscous fluid overlain by a dense perfectly
plastic layer. The structure is subject to either horizontal extension or
shortening and models rocksalt under a brittle
overburden. Considering the viscosity of the buoyant layer to be much less than
the effective viscosity of the overlying layer, we obtained the following
results. The instability pattern of the structure is similar to that of a
perfectly plastic structure. The characteristic wavelength, corresponding to
the most unstable mode, increases initially with the thickness ratio between
the lower and upper layers, but then decreases by a series of abrupt jumps.
This is in contrast to the result that the characteristic wavelength approaches
to a constant value with the increasing thickness ratio in the case of viscous
layers. The buckling instability induced by rapid horizontal straining
overwhelms the gravitational instability, and the growth rate of the
instability depends linearly on the effective viscosity ratio. We studied
models of diapirism in the Great Kavir,
Iran. We show that a small interdiapir spacing in the
salt canopy province and a wide range of the spacings
in the salt pillow province of the region can be explained by the perfectly
plastic sedimentary overburden and horizontal shortening.
Stress in
descending lithospheric slabs
I
examined the effects of viscous flow, phase transition, and dehydration on
stress in a relic slab to explain the intermediate-depth seismic activity in
the Vrancea region. I developed a two-dimensional finite-element model of a
slab gravitationally sinking in the mantle which predicts (i)
downward extension in the slab as inferred from the stress axes of earthquakes,
(ii) the maximum stress occurring in the depth range of 70 km to 160 km, and
(iii) a very narrow area of the maximum stress. The depth distribution of the
annual average seismic energy released in earthquakes has a shape similar to
that of the depth distribution of the stress in the slab. Estimations of the
cumulative annual seismic moment observed and associated with the volume change
due to the basalt-eclogite phase changes in the
oceanic slab indicate that a pure phase-transition model cannot explain solely
the intermediate-depth earthquakes in the region. We consider that one of
realistic mechanisms for triggering these events in the Vrancea slab can be the
dehydration of rocks which makes fluid-assisted faulting possible. The approach
can be applied to other regions of intermediate-depth seismicity.
Block-and-fault
dynamics of the lithosphere and earthquakes in the SE-Carpathians
I
studied a numerical model of block-and-fault dynamics of the lithosphere
beneath the earthquake-prone Vrancea region. The model presents a system of
absolutely rigid blocks separated by infinitely thin plane faults. The
interaction of the blocks along the fault planes and with the surrounding
medium is assumed to be a viscoelastic. Motions of boundary blocks cause the
displacements of the block system. The velocities of the motions are found from
a model of mantle flows induced by a sinking slab beneath the Vrancea region.
When a ratio of stress to pressure for some portion of a fault plane exceeds a
certain strength level, a stress-drop ('earthquake') occurs. As a result of the
numerical simulations catalogues of synthetic earthquakes are produced. Several
numerical experiments for various model parameters illustrate that the spatial
distribution of synthetic events is significantly sensitive to the directions
of the block movements. Small variations in a slab rotation control the pattern
of the synthetic seismicity. The results of the analysis indicate that the
catalogues obtained by the simulation of the block structure dynamics have
certain features similar to those of the real earthquake catalogue of the
Vrancea region. The results are in a good agreement with recent seismic
tomography data.
Computational
Methodology to Study Geodynamic Problems
I
developed a new two-dimensional Eulerian finite element methodology to study
geodynamical problems where chemical composition changes abruptly across
interfaces. The method combines a Galerkin-spline
technique with a method of integration over regions bounded by advected interfaces to represent discontinuous variations
of material parameters. It allows to directly approximate a natural free
surface position, instead of a posteriori calculation of topography from the
normal stress at the upper free-slip boundary. The suggested approach avoids
the difficulties due to material discontinuities at intermediate boundaries,
like Moho or the earth's surface, and is also free from the deficiencies of the
Lagrangian approach always resulting in mesh
distortion. Also I developed new three-dimensional methods for study of
lithospheric deformation and mantle convection. Special computer codes were
written for parallel supercomputers with distributed memory.