S. N. Korobeynikov

Computer Modeling of Granite Gneiss Diapirism in the Earth's Crust: Controlling Factors, Duration, and Temperature Regime

This paper is devoted to the modeling the granite gneiss formation by means of diapiric upwelling. The natural examples of granitic diapirism in the Precambrian granite-greenstone belts and complexes of metamorphic cores are described. A new approach is proposed to describe the partial melting and development of gravity instability in the crustal granitic layer, which experienced heat impact and melting during intrusion of basaltic melt. Rheology of partially melted material and surrounding medium is regarded to be temperature-dependent, following either plasticity or creep (non-Newtonian viscosity) law. Modeling results show that crustal rheology plays a significant role in the character of diapirism (shape of upwelling bodies, duration of the process, and width of thermal aureole). The rates of upwelling within the crust behaving as elastoplastic body are orders of magnitude higher (meters to tens meters per year) than those obtained for creep (viscous) liquid model (0.8 cm/yr). Modeling results revealed that the limiting depth of upwelling of partially crystallized melt, with allowance for temperature dependence of creep, corresponds to the isotherm of 400°C.

Formation and upwelling of mantle diapirs through the cratonic lithosphere: Numerical thermomechanical modeling

This paper reports the results of the numerical modeling of gravitationally instable processes in the lithospheric mantle of ancient cratons. The gravitational instability is considered as a result of melting at the lithosphere base owing to its local heating by anomalous mantle. Modeling was based on a finite element method in 2D formulation and took into account the geological structure and thermomechanical parameters of the lithosphere of the Siberian platform. Numerical results revealed the main tendencies in the mantle diapirisim of the mafic and ultramafic magma ascending through the “cold” high-viscosity lithosphere. It was shown that the shape of diapiric magmatic bodies is controlled by realistic visco-elastic-plastic rheology of lithosphere. The ascent of diapir in lithosphere was modeled for diverse regimes differing in duration, temperature field, and upwelling depth. It was concluded that the ascent of melt through lithosphere to the crust-mantle boundary is mainly controlled by rheology, and conditions of oscillatory diapirism with recurrent magma replenishments were modeled. Modeling results may shed light on some features related to the trap magmatism of the Siberian igneous province. The duration and rate of magma upwelling as well as the parameters of periodical magma upwelling were estimated and attempt was made to explain the high-velocity seismic anomalies that were recorded in the subcrustal regions of the Siberian platform.

Computer Modeling of Granite Magma Diapirism in the Earth's Crust

A new point of view describing processes of partial melting and development of gravitational instability in a thickening crust with increased thickness of the granite layer is suggested. Numeral experiments support the following main conclusions. The critical volume of partially melted material should be formed for the beginning of flotation in a gravitational field. Due to model estimations, the height of the melting area in the granite crust should be not less than 6–7 km. A mushroom-shaped form of the floating body was observed in all models regardless of the thermal source size (fixed or variable width): the high temperature channel (magma leader) and head body of the diapir are formed. The height of diapir floating depends on rheological features of the surrounding crust: 10 times increase in the yield strength (from 1 to 10 MPa) while temperature decrease confines the possible level of rising to a depth of 15–16 km. An elevation of about 750 m is formed in the day surface relief above the axis part of the diapir.

Computer Modeling of Underthrusting and Subduction under Conditions of Gabbro-Eclogite Transition in the Mantle

The rheology of Newton’s viscous or nonlinear liq-uid with viscosity dependent on the temperature anddeformation rate is often used in modeling of subduc-tion (see, for example, [1]). The temperature in the sub-ducted plate and accretionary and mantle wedge, aswell as the velocity of motion in the suprasubductionmantle wedge, is commonly calculated in this case [2].The motion of the plate itself is set with a constantvelocity [3] or is deduced from the convection of thelithospheric mantle [4, 5]. The subduction is oftenreproduced as plunging flow of viscous and cold mate-rial (see, for example, [6]). The interpretation of numer-ical modeling is based on the viscosity or temperaturefield. The region with elevated viscosity (lowered tem-perature) relative to the adjacent medium is regarded asa subducted material (see, for example, [7]). In suchmodels, the subducted plate has a typical shape of thesinking plume or drops [8].We used a different approach to modeling of geody-namic and tectonic processes elaborated previously in[9, 10]. A model with elastoplastic rheology of the sub-duction plate and lithospheric mantle with delineationof the slab boundary and the adjacent mantle is pro-posed in terms of the solid mechanics as an alternativeto models with viscous liquid. The conditions of fittingthe physicomechanical laws are strictly fulfilled at theplate–mantle interface: we set up the frictionless condi-tions in the case of serpentinization of the mantlewedge or friction conditions according to the Cou-lomb–Mohr law. The setting of the problem assumesthe existence of interacting geological bodies: subduc-tion slab, underlying mantle, and rigid continentalcrust. In contrast to the aforementioned models, inwhich the subduction plate may be determined at adepth only rather conditionally, the proposed approachallows exact location of all interfaces between geologi-cal bodies and provides tracing of their progressivedeformation. Thus, the problem of determining theshape of the plunging slab, its boundaries with the adja-cent mantle, and the stress state of the subducted mate-rial may be solved within the framework of the newapproach.We consider two (oceanic and continental) platesthat were in equilibrium with the mantle at the time ofcollision and have an inclined interface. The plates areequilibrated by contact forces owing to the pressure ofthe right plate on the left one (Fig. 1). It is suggestedthat the oceanic plate collides with the continental platewith a rate of

Numerical modeling of mantle diapirism as a cause of intracontinental rifting

The numerical model of mantle diapirism and active rifting is developed. The model describes the possibility of extension and thinning of the Earth’s crust under the action of a local 100-km long heat source in the sublithospheric mantle, which causes melting and rising of the magmatic diapir through the cratonic lithosphere. The model combines the mechanisms of the uplifting of the anomalously hot material due to its gravitational instability, underplating of magma beneath the continental crust, and its extension by the forces of the convective flows at the base of the plate. The obtained results shed light on some geological features of the joint formation of the large Vilyui igneous province and Vilyui sedimentary basin.

Influence of the choice of a rheological law on computer simulation results of slab subduction

The influence of the choice of the type of yield surface for elastoplastic materials and material constants for the slab and the mantle on the scenario of mathematical modeling of the slab collision is investigated. Computer simulation is performed by the FEM numerical solution of nonlinear equations for deformable solid mechanics, using the MSC.Marc 2005 code. The simulation results essentially depend on the choice of material constants for the slab and the mantle, as well as on the type of the yield surface for the elastoplastic material of the subduction slab. The numerical simulations demonstrate that the primary driving mechanism of subduction can be a geometrical inhomogeneity of the subduction slab near the zone of the slab collision with simultaneous consideration of consolidation of the slabmaterial as the slab descends into the mantle.

Computer Simulation of Underthrusting and Subduction Due to Collision of Slabs

We come up with a mathematical simulation of a collision between lithospheric slabs (plates) where one slab is forced into the mantle beneath another. Problems of the Earth’s crust and mantle deformation are solved numerically: for spatial discretization of equations of deformable solid mechanics, a finite-element method is used, and for evolution of the collision process, a stepwise integration of quasistatic deformation equations is applied. Problems of plate motion are solved within a geometrically nonlinear setting in a two-dimensional approximation (plane deformation) with due regard for large deformations of bodies and contact interactions of slabs with the mantle. A numerical solution is obtained via a MSC.Marc 2005 code, encompassing formulations of equations with required types of nonlinearities. A part of the Earth’s crust that has no tendency to delving into the mantle is simulated by a prescribed motion of a rigid body. A part of the Earth’s crust that should sink by virtue of properties of initial geometry is simulated as a deformable solid made up of elastoplastic strain-hardening material. The mantle is simulated by an ideal elastoplastic material with a low yield stress value. We are concerned with parts of the Earth’s crust that have different geometric parameters. Computer simulation of plate collision shows that under standard conditions, underthrusting of one slab beneath another occurs; at sites of initial thickening of a slab in a contact zone, subduction (deep sinking) of the slab into the mantle is expected. In the latter case account should be taken of a well-known experimental fact, that of material compaction of the sunken piece of a slab.

The influence of crustal rheology on plate subduction based on numerical modeling results

Computer simulation of subduction was performed using nonlinear equations of deformable solid mechanics encompassing all types of nonlinearity: geometric, physical, and contact. This study presents a numerical model of subduction with allowance for the gabbro-to-eclogite phase transition. The model rheology is a plastic compressible material (Mohr-Coulomb law for a deformed rock material). It was shown that deep subduction can be modeled well with the selection of appropriate parameters of rock plasticity providing the initial thickening in the subducting slab nose.

Mathematical modeling of magma fracturing and dike formation

Magmatic phenomena play an important role inmass and heat transport within the Earth’s shells.Transport of magma in the lithosphere is determinedto a great extent by the force of gravity. Intrusive magmatism and diapirism are one of the manifestations ofgravity instability and heat and mass transport, i.e.,lifting of lighter matter in a viscoplastic medium anddescending of heavy matter. Diapirism in the lithosphere has been investigated in detail within physical,experimental [1], and mathematical [2, 3] models.The cause of magmatic intrusions is the pressure,which is arised in the gravity field mainly by the density difference between the melt and the host rocks.Horizontal tectonic motions also play a specific role inthis process. Magma intrusions are controlled by the“magma fracturing” mechanism, which has much incommon with hydraulic fracturing [4]. Magma fracturing is the process of the destruction of rocks bymeans of a fissure filled with melt under pressure.Intrusive magmatic bodies are formed as a result ofmagma fracturing. Their formation depends on twofactors: the fissuring rate of the host rocks and the ability of the magma to fill the channel in the crack tip [5].Many publications are dedicated to the mechanism ofmagma fracturing, in which the authors use laboratorymodeling [6] and analytical and numerical methods ofmathematical modeling [5, 7–9]. The influence ofmagma in these mathematical models was described ina simplified form: an internal constant pressure wasspecified in the growing fissure, while the motion ofthe magma was not considered.We developed a model of magma fracturing, inwhich we emphasized the interaction between mediawith contrasting rheological and thermophysicalproperties: a lowviscosity hightemperature melt andbrittleelastic rock at a temperature below melting.The melt material is modeled by a Newtonian fluid,while the rock material is an elastic medium with thepossibility of brittle damage. The breaking occurswhen the maximum stress at the solid body point σ

Computer modeling of folding in rocks

The results of numerical mathematical modeling of rock deformations under compression are given. The numerical solutions are obtained using discretization of the equations of the solid mechanics with the finite element method. The model of an ideal elastic-plastic material with a Huber-Mises yield surface was used in the calculations. The layered medium structure is taken into account in modeling of the compression of layers of the lower/middle crust on a local scale. The natural folds in strongly deformed metamorphic sequences were reproduced by the mathematical deformation models. It is shown that folding in the lower part of the Earth’s crust is possible, when the yield stresses of the host rocks are approximately two orders of magnitude less than those of the hard layers. The effect of the boundary conditions and thickness of the compressed rocks on the folding is shown.