Solving the conundrum of extensional folding in metamorphic core complexes

Folds with axial traces parallel to the extension direction are a common feature in continental detachment systems and metamorphic core complexes. Yet, how they form has been puzzling for many decades. Here, we show that the key to solving the conundrum lies in revising the long-held single-scale view toward natural deformation and application of kinematic models. We demonstrate that extensional folding can result naturally from the partitioned stress field in competent layers in plate-scale extension and transtension deformations. Competent layers that develop extension folds should be regarded as rheological inclusions in the lithosphere rather than infinitely extending plates clamped at system boundaries and subjected to system boundary conditions.

so is the western Grenville Province (Fig.1). The strains reflected by extension folds cannot be explained by the transtension model.
Here, we demonstrate that the difficulty in understanding extensional folding is due to our single-scale approach toward lithosphere deformation. We think it is unrealistic to regard layers that develop extension folds as infinitely-extending elastic or viscous plates that are subjected to the system boundary conditions. It is more realistic to regard them as rheologically heterogeneous inclusions embedded in the lithospheric mass undergoing macroscale transtension. It is the partitioned stress and strain field in these layers, rather than the homogenized stress and strain field of the macroscale transtension, that is relevant for extensional folding. We show that partitioned stress field can produce extension folds in competent layers in extension and transtension settings.

Partitioned stress field in a competent layer-like inclusion
Eshelby's inclusion solution [14][15][16] relates the stress and strain field inside an ellipsoid inclusion to the uniform macroscale field by: Where both the inclusion and matrix are elastic (an all-elastic system), ε and E are elastic strains.
Where both the inclusion and the matrix are Newtonian viscous (an all-viscous system), ε and E are strain rates. Where the matrix is viscous and the inclusion elastic, ε is elastic strain tensor of the inclusion and E the viscous strain rate tensor of the matrix. C is the elastic or viscous stiffness of the matrix depending on its rheology and S the Eshelby tensor for the inclusion.
We use the fundamental interaction equation (1) and the equivalent-inclusion approach 14,17 to derive analytical solutions for the partitioned stress in competent layer-like elements under macroscale transtension. The following assumptions are made. First, we regard the competent layer as a flat oblate inclusion with principal semi-axes 1 2 3 a a a =  . By varying the aspect ratio 3 1 a a  = , the inclusion may represent different bodies including a uniformly-thick plate ( 0  → ).
We believe that, on the macroscale suitable for continental transtension, no geological unit is rheologically continuous throughout the deformation zone and subjected to the system boundary conditions. Second, for simplicity we assume that the rheology of the lithosphere on the macroscale and that of the inclusion are linear, isotropic, and incompressible. Third, on the macroscale appropriate for transtension, the lithosphere is approximated by a Homogeneous-Equivalent Medium (HEM) whose rheology is represented by the homogenized rheology of all the rheological elements in the macroscale representative volume 15,16,18,19 . We use the convention of tensile stress being positive so that tensile normal stress corresponds to positive extension strain.
To model extension folds in North America Cordillera metamorphic core complexes, we consider stress partitioning in an all-elastic system. Equation (1) leads to a strain partitioning 15 ) which, upon using Hooke's law, is rewritten in the following stress partitioning form: where r is the ratio of the elastic shear modulus of the element  to that of the HEM s  and To model extension folds in the high metamorphic grade rocks in the Grenville Province, we regard the HEM as a Newtonian viscous material. The competent layers that develop extension folds may also be viscous, with a higher viscosity than the HEM. In such a case, the above results (equations (3)) can be applied if r in the equations is taken as the viscosity ratio between the inclusion and the HEM. We consider an additional situation where the element is instantaneously elastic; infinitesimal elastic strains of the layer are converted to permanent strains continually to give rise to the observed folds. In such a scenario, we apply the equivalent-inclusion approach to equation (1) by substituting the elastic inclusion with a fictitious inclusion of the same rheology as the HEM but with the same internal stress as the elastic inclusion. The interaction equation (1) becomes ( ) , where, as in equation (2)

Folding of a horizontal competent layer in transtension
The macroscale stress tensor for transtension, expressed in the coordinate system xyz ( Fig.2a), is of the following form (see Supplementary Information): Taking the eigenvalues, the corresponding principal stresses for an all-elastic system are: and those for an elastic layer in a viscous HEM are: 6 Note equations (8)  . The case of  =0 means that the competent layer is an infinitely-extending sheet, and the solution converges to the single-scale solution of a horizontal sheet under macroscale transtension. In reality, a competent layer is an inclusion which means that  is small but not zero.
It can be readily confirmed that the principal partitioned stress axes are parallel to the macroscale principal stresses for both equations (7) and (8). Therefore, whether treated as an allelastic system, an all-viscous system, or as elastic layers in a viscous HEM, the partitioned principal stresses in a horizontal layer are parallel to the macroscale principal stresses but are increased in magnitude by approximately a factor of r or 1  − . Significantly, the horizontal 2  is always compressive and does not vanish (equations (7b) and (8b)) even if the system is in pure As there are no available analytical solutions for the buckling of a flat oblate element embedded in an elastic or viscous matrix, we use the theory of cylindrical buckling of a rectangular elastic plate under in-plate loading 21,22 to provide an approximate analysis of a horizontal elastic layer element in an elastic or viscous HEM under transtension.
The tendency for a competent elastic layer to buckle is determined by its flexural rigidity 21 , which for an incompressible elastic plate is ( h the thickness of the plate). The analysis in refs. 21,22 shows that the critical load P for the layer to buckle is . This yields an estimate for the critical 2  for buckling instability to develop in the layer: According to equation (8b), 2  increases in magnitude asymptotically to ( ) The critical aspect ratio for a horizontal layer is ( ) Similarly, in an all-elastic (or all viscous) system, we can use equation (7b) to find the critical aspect ratio. In the pure extension situation, the critical aspect ratio for a horizontal layer to develop buckling is . Thus, all horizontal layer elements with aspect ratios c   or ' c   , depending on the rheologies considered, will eventually reach the critical compression for buckling and develop extension folds in transtension (Fig.3). Folding of a layer increases its effective aspect ratio. The folding of a layer will continue until the effective aspect ratio is below the critical value.

Folding of an inclined competent layer-like inclusion in transtension
Equations (3) and Equations (5) (3) and (5)) which remain at the same level or smaller than the macroscale components ( 33    ). Therefore, in a strongly competent layer, 3  may be negligible compared to other components.
This means that two of the three principal partitioned stress axes are always nearly parallel to the layer and the third principal stress is normal to the layer, despite of the layer orientation. This stress state is consistent with the field observation that competent layers are more prone to buckling and boudinage instabilities to develop folds and boudinage structures. It also implies that the classical elastic plate theory is a reasonable approximation for the deformation of thin competent layers in transtension.
We denote the strike and dip of an inclined competent layer by  and  respectively. It turns out that the orientation of the partitioned principal stresses in a competent layer is 8 independent of whether the system is all-elastic, all-viscous, or an elastic layer in a viscous HEM. and  is measured with respect to the -y axis (Fig.2a).
If buckling folds initiate in the layer parallel to 1  -axis, we can use equation (10) (2)). Equation (2) suggests that ij  in a mylonite shear zone are likely distinct from and smaller in magnitude than ij  of the bulk lithosphere making it unjustified to assume that ij ij   or that the invariants of the two are approximately equal.

Acknowledgments:
We are grateful to Lucy X. Lu, Ankit Bhandari, and Rui Yang for many discussions. We tank Lucy Lu for help with plotting Fig.4  Author contributions: DJ developed the theory, derived all the equations, and wrote the paper. CL conducted field structural analysis in the Grenville Front region and prepared Fig.1. Both were involved with the interpretation of field data and discussion of the ideas of the paper.

Competing interests: Authors declare no competing interests.
Data and materials availability: All data is available in the main text or the supplementary materials.

Additional information
Supplementary Text Figure S1 Supplementary References  a, Transtension geometry and coordinate system xyz used to define the deformation. The horizontal divergence vector v is at an angle  relative to the system normal (thex axis). An initial square is deformed progressively into a parallelogram with its lower-right corner (p) moving along v (dashed line to the final position p'). The horizontal principal stretching ( 1 E ) axis is at 2  relative to thex axis and the horizontal principal shortening 2 E -axis is normal to 1 E -axis. They also are parallel to the two horizontal principal stresses 1  and 2  (not shown in Figure but referred to in text). The two horizontal principal finite strains, measured by stretch, are 12 SS  and 1 S orientation is between 1 E -axis and v depending on the magnitude of finite strain. b, When  is small, 2 S remains close to 1 regardless of the finite strain rendering minor horizontal principal shortening (1-2 S ).

Transtension
The velocity gradient tensor for macroscale transtension considered in coordinate system xyz (Fig.2 in the text) can be expressed as: cos 0 0 where v is the divergence velocity, d the width of transtension, and  the angle of divergence velocity measured relative to the transtension normal (the x-axis). The The three principal stretches are obtained (see ref. 2 for more detail) by taking the eigenvalues of the tensor The three principal stretches are:  Fig. 2b in the text is based on Eq. s5b.

Inclined Competent Layer in Transtension
For an inclined inclusion with strike  and dip angle  , we set the inclusion coordinate system x'y'z' such that the x'-axis is along the dipline of the layer and pointing down, y'-axis is parallel to the strike of the layer, and z'-axis is normal to the layer and pointing up (Fig.S1). The rotation matrix ij Q relating xyz and x'y'z' coordinates are 3 : cos cos sin cos sin sin cos 0 cos sin sin sin cos The macroscale stress tensor expressed in the inclusion coordinate system are obtained by tensor The partitioned stress tensor components are obtained by inserting ' ij  in place of ij  in equations (3) and (5). In either case, the orientation of the 1 - axis relative to x'-axis is: With  obtained from Eq.s8, the trend and plunge for 1 - axis (equation (10) in text) can be obtained readily from simple trigonometry.

Derivation of Equations (5)
For shear stresses ij  ( ij  ), we expand equation (6)  Fig. S1. Coordinate system x'y'z' for an inclined layer-like inclusion. Coordinate system xyz is the one in which macroscale transtension is defined. The yellow plane is the inclined layer,  and  are respectively strike and dip angles. The x'-axis is along the dipline of the layer and pointing down, y'-axis is parallel to the strike of the layer, and z'-axis (not shown) is normal to the layer and pointing up.  is the angle of the partitioned principal 1  axis relative to the x'axis.