Tectonophysics, 78 (1981) 651-660 Elsevier Scientific Publishing Company, Amsterdam - Printed in The Netherlands
651
STRESS TENSOR DETERMINATION FROM THE STUDY OF e TWINS IN CALCITE: A LINEAR PROGRAMMING METHOD
Ph. LAURENT
l, Ph. BERNARD
2, G. VASSEUR
2 and A. ETCHECOPAR
’
’ Laboratoire de Ge’ologie Structurale, LA 266 U.S.T.L. 34060 Montpellier Cddex (France) ’ Centre de Ge’ologie et Gt!ophysique, LP CNRS, U.S.T.L. 34060 Montpellier Ce’dex (France) (Received January 26, 1981)
ABSTRACT Laurent, Ph., Bernard, Ph., Vasseur, G. and Etchecopar, A., 1981. Stress tensor determination from the study of e twins in calcite: a linear programming method. In: G.S. Lister, H.-J. Behr, K. Weber and H.J. Zwart (Editors), The Effect of Deformation on Rocks. Tectonophysics, 78: 651-660. A new method for the determination of the stress tensor from calcite e twin lamellae has been discovered. All e planes (twinned and untwinned) are taken into account to enable derivation of tensor characteristics. The basic assumption is that the resolved shear stress, a linear function of the unknown stress tensor components, is greater (respectively lower) than the critical resolved shear stress for the twinned (respectively untwinned) planes; the relevant inequality constraints can be solved using linear programming algorithms. The application to experimental (Carrara marble) and natural samples (a goniatite from Morocco) has come up against one major difficulty - polyphase deformation. Thus, before applying this method, it is necessary to operate some selection of the data. One procedure based upon linear programming is proposed. The results are consistent with those expected both for the experimental and natural cases. For the Carrara marble, an estimate of the critical resolved shear stress (tc) for twinning gave a value around 300-400 bar. The method can give, in principle, the stress tensor orientation and the relative magnitude of principal stresses simultaneously.
INTRODUCTION Dynamic analysis of deformed rocks has been carried out from the study of numerous different crystals (Carter and Raleigh, 1969). Amongst them, the calcite is one of the most useful. Common in crustal rocks and belonging to a high symmetry group (R~c), this mineral is experimentally well documented (Turner et al.,.1954; Friedman, 196’7). Moreover, the gliding systems 0040-1951/81/0000-0000/$02.50
@ 1981 Elsevier Scientific
Publishing Company
652
and twins are directly measurable with an optical microscope equipped with a four-axes U-stage. Under low temperature and anhydrostatic pressure conditions, it is classically established that strain in calcite is easily achieved by twinning on e(Oli2). A critical resolved shear stress value (t,) of 100 bar (Jamison and Spang, 1976), probably less (M.S. Paterson, pers. commun. 1979) in the positive sense, is sufficient for twinning to occur. This value seems to be independent of temperature in the range 0-400°C and at hydrostatic pressure (Turner et al., 1954). PREVIOUS
METHODS
Previous methods for estimating strain (Groshong, 1974) and stress (Turner, 1953) analysis have led to quite consistent results (e.g. Friedman et al., 1976). A brief review follows: (a) Stress orientation is determined statistically: the most favourable orientation of the principal stress axes for twinning to have occurred in each twinned crystal is defined, considering separately or together one, two or three e systems (e.g. Nissen, 1964; Friedman and Conger, 1964). (b) Differential stress value (Jamison and Spang, 1976) can be obtained by measuring the percentage of grains which are twinned on one, two or three different e planes. The main basic assumptions are an initial isotropic fabric of c-axes and an “irrotational uniaxial stress field”; (it corresponds to what we call an axially symmetric stress field). (c) Principal stress ratio (Friedman and Heard, 1974) has been determined in naturally deformed samples by the measurement of the e twin lamellae abundance and by comparison with experimental data. These methods do not allow simultaneously the definition of the principal stress directions and ratio. The stress tensor directions are determined for each crystal and the final result is gained by the statistical analysis of the whole. This assumption is not mechanically justified, particularly if the intermediate principal stress is different from the two other ones. Similar problems arise for the interpretation of slickensides data (Carey, 1979). Moreover, these methods do not make a total use of the information that some e planes are twinned and others not. STATEMENT
OF THE METHOD
For a given medium composed of numerous crystals, twinning is assumed to result from the application of one homogeneous stress tensor ts . No further assumption about u is required and all the twinned or untwinned e planes will be taken into account. Thus for a given e plane with unit normal vector n, twinning occurs (respectively does not occur) if the resolved shear stress t,, i.e. the component of the pressure force exerted up on this plane along the twinning
653
direction with unit vector 111,is greater (respectively lower) than some positive critical yield value, t,, given experimentally. The mathematical expression of this resolved shear stress is:
(1)
t,=m.a-n This assumption
is consistent with the fact that normal stress has no effect on twinning (Friedman, 1964). Therefore, only the deviatoric part (rn of the stress tensor c has to be considered and can be characterized by 6 unknowns, since it is symmetric :
6,
=
x1
x*
x3
xz
xl
x5
x5
x5 I
; x3
(2)
with xl+xq+x6=0
(3)
Each of the data consists of: ni, the unit vector normal to the ith e plane, mi the corresponding unit vector in the twinning direction and a Boolean variable Vi which is respectively “true” or “false” if twinning has occurred or not along mi. With respect to the previous assumption, each of the N independent data yields: mi * CD ’ TZi< t,
if
Vi = “false”
.K?li’ by . ?Zi> t,
if
Vi = “true”
(4)
with i = 1, . . . . N. For a given t, , the left-hand sides of eqs. 3 and 4 are linear in X, . . . X6. The inverse problem of deriving en from the data is reduced to that of finding X = X1 . . . X6 a vector of R6, such that eqs. 3 and 4 are satisfied. In general there are two possibilities: (a) The problem has a solution, but this solution is not unique: the set of solutions X is a linear manifold of R6 which can be shown to be limited by a convex polyhedron (Goldman, 1956). Various features of this convex set ban be described using the algorithms of linear programming (Da&zig, 1963); for example, max. (min.) value of Xi or of any linear function of X are obtained immediately using the simplexe algorithm: (b) The problem has no (feasible) solution: in such a case it is possible to obtain an approximate “best” solution by the use of positive slack variables t,, . . . . &,, such that: mi * CD * ni < t, + ti
if
Vi = “false”
mi.an
if
Vi = “true”
’ l2i > t, -
ti
(5)
654
The “best” solution is the one which minimizes the positive sum: X;j:T Ei (i.e. which gives the smallest residuals with respect to the data). This solution can also be obtained using the algorithms of linear programming. The inferred tensors are proportionally dependent upon the value of t,, for which an arbitrary value of 1 is chosen. The actual value of the deviatoric stress can then be obtained from t, values deduced from experimental deformation. PROBLEM
OF SUPERPOSITION
In actual cases, even from laboratory experiments, rocks have suffered different phases corresponding to different stress tensors. It is important to notice that, in accordance with the previous assumptions: (1) for each untwinned e plane, none of the stress tensors has given rise to a resolved shear stress t, greater than the critical one t, . (2) for each twinned e plane, at least one of the stress tensors has given a t, greater than t,. In actual deformation, one has to keep in mind that eventually a few twinned or untwinned e planes may not satisfy these relationships. For a polyphased deformation, the use of the previous method which yields to the definition of one tensor for all the data would lead to a meaningless solution. Therefore it is necessary to perform some selection of the input data before applying this method. One can propose several procedures for this selection. One of the most convenient ways can be described as follows; for each twinned plane, find the tensor (Jo such that: (1) the resolved shear stress t, for all the untwinned planes is smaller than the critical one. (2) the resolved shear stress t, for this twinned plane is as large as possible. Stated as above, this is also a standard linear programming problem which can be solved using classical algorithms. As will be shown later, the principal directions of the various tensors obtained from the various twinned planes generally cluster together into a few separated domains which leads to a recognition of the various phases. APPLICATIONS
A computer data. 1. Experimentally
program
was written
deformed
and tested
using monophased
synthetic
sample
A sample of carrara marble, kindly supplied by Dr. Chaye d’Albissin (University of Paris), has been studied. A core cylinder of this rock was subjected to an axially symmetric stress field (a, = 2460bar, u2 = u3 = 500bar) under room temperature for 45 min. The longitudinal shortening was 1%. The studied thin section is parallel to u,. 50 c-axes were determined (150 e planes) using a U-stage. The density diagram (Fig. 1) of all the e planes is subisotropic.
655
B Fig. 1. Density diagram of all the measured e planes (150 measurements) in the experimentally deformed sample of Carrara marble. T = top; B = bottom; N = north; S = south.
MAX.
28%
Fig. 2. Density diagram of ~1 directions (which correspond to the maximum of tm for each twinned plane) in the experimentally deformed sample of Carrara marble. I, ZZand ZZZcorrespond to clusters of twinned planes selected (cf. text). 2’ = top;Z? = bottom; N = north; S = south.
656
In a first step, the separation method described above was applied to each twinned plane: on Fig. 2 is shown the density diagram of u1 directions corresponding to each twinned plane. Obviously three clusters can be distinguished: in the first one, labelled I (73% of the twin planes), the u1 direction corresponds to the axis of the core (i.e. the direction of the experimental 0,) whereas the third one, labelled III (19% of the twin planes) shows a u1 direction subperpendicular to the previous one, in the plane of the thin section. This cluster number III is interpreted as a result of a previous natural deformation which will be evidenced in a check sample. The cluster labelled II (4% of the twin planes) corresponds to an intermediate u1 direction. A similar study of a check sample (150 e planes with 61 twinned planes) was performed. More than 70% of the u1 directions corresponding to each twinned plane obviously belong to the same class (Fig. 3) which emphasizes the existence of a previous deformation phase. Unfortunately, the relative orientations of the check sample and of the laboratory tested one are not known. However, it seems reasonable to consider that this class is related to the cluster number III of the experimentally deformed sample, i.e. to the natural phase. It is interesting to compare the obtained results with those gained using Turner’s (1953) method on the same deformed sample (parallel to the core axis). The ui directions deduced from this method are displayed on Fig. 4. The resulting pattern is not very different from the one of Fig. 2, but much more diffuse; hence the proposed method gives a selected cluster of uI direction with a better accuracy.
Fig. 3. Density diagram of ~1 directions (which correspond to the maximum of tm for each twinned plane) in the experimentally undeformed sample of Carrara marble. Density contours as for Fig. 2.
657
B Fig. 4. Density diagram of ~1 directions determined using Turner’s (1953) method, experimentally deformed sample of Carrara marble. Density contours as for Fig. top; B = bottom; N = north; S = south.
in the 2. T=
Next, the twinned planes corresponding to selected clusters of Fig. 2 are studied. These selected data together with all the untwinned planes are then re-introduced in the optimizing algorithm in order to obtain one tensor explaining all these data. Two solutions corresponding to two different selecTABLE
I
A comparison directions.
between
the
computed
(52
01
Expected azimuth Plunge Actual magnitude (bar) of deviatoric part Solution Computed Plunge Relative
and expected
9o” O0 1307
experimental
t* Cl
(53
O0 O-90” -653
principal
t* c2
stresses
and
t* c3
O0 o-9o” -653
(1) azimuth
magnitude
111° 0.5 O 4.14
21° 37O -1.40
202O 53O -2.14
315
466
238
6’ loo 4.24
197 80 -3.16
383
2720
206
(G = 1) Solution Computed Plunge Relative (k
(2) azimuth
magnitude
96’ 2O 3.41
= 1)
t=T is the value of the critical resolved shear stress t, computed [actual magnitude of Ui (bat)/relative magnitude of ~i].
using
the
values
of
Ui
658
Cons of twinned data have been computed: for solution (l), only the twin planes (58% of the twins), corresponding to those data of cluster number I (Fig. 2) which gives a solution in the vicinity of the peak (density higher than 10%) were selected. For solution (2), twin planes corresponding to clusters number I and III (Fig. 2) were selected (77% of the twins). In Table I, the results obtained with these two solutions (principal stresses and directions) are compared to the values expected on an experimental basis. In solution (l), for ul, a difference of 20” in azimuth between the expected direction and the computed one is obtained. In solution (2) the difference is only 6”. The difference between the two solutions (1) and (2) emphasizes the major importance of the selection of the data. A comparison between the values of the experimental stress (deviatoric part) with the computed one (scaled by t,) leads to an evaluation of actual magnitude of t, . Using only a, the obtained value of t, is 315 bar (= 1307 bar/4.14) for solution (1) and 383 bar (= 1307 bar/3.41) for solution (2). The estimate of t, from u2 and u3 which shows a much larger range is not really meaningful: this large range is due to the difference to axial symmetry in the computed solutions (es = a3). Such a high value of the critical resolved shear stress (300-400 bar) is astonishing by comparison with recent laboratory data. This difference may be related to the polycrystalline nature of the sample. 2. Naturally deformed
sample
A sample of fossiliferous Devonian limestone from Morocco has been chosen because : (1) It was expected to have been deformed by only one phase of deformation. This phase is responsible for large-scale upright folds, faults and stylolitization. (2) The natural penetrative deformation is slight (undistorted fossils imbedded into a stylolitized micritic matrix): the stratification plane (So) is horizontal at the point of sampling. (3) Stylolitization may give a good approximation of ui direction. (4) Numerous fossils are filled with large crystals of slightly deformed clear calcite. The c-axes are more or less randomly distributed. In a thin section parallel to So, 91 e planes (31 twinned planes) have been measured in a goniatite in order to determine the stress tensor. It was soon realized that all the twinned planes could not be explained by one unique tensor (no feasible solution for eq. 4). Therefore, a selection procedure was performed. Figure 5 presents the density diagram of u1 directions deduced from every twinned plane. More than 50% of the directions are E-W which are consistent with the CJ, directions deduced from stylolite peaks. However, other u1 directions are also observed which emphasizes the existence of other phases.
659
Fig. 5. Density diagram of (~1 directions corresponding to the maximum of t, for each twinned plane in the naturally deformed sample. Density contours as for Fig. 2. N = north;S = south; E = east; W = west.
The selection of the twinned planes corresponding to the E-W cluster for u1 yields to the solution given below: Ul =
1.91 t,
A, = 258”
plunge 4”
C72 = -0.14
t,
A, = 163”
plunge 45”
03 = -1.77
t,
A, = 351”
plunge 45”
The direction of u1 is found to agree with the direction of stylolite peaks (270-255”). The dips of u2 and u3 directions do not seem to be geologically meaningful. In this case, due to the low number of twinned data retained, the problem is not constrained enough. The solution is not unique and other solutions are likely to exist. More data are necessary to give further constraints. DISCUSSION
The problem of determination of the stress state from twin measurements is stated as an inverse problem using a few reasonable assumptions. The preliminary results presented here are consistent with the expected ones which, seems to justify the fundamentals of this method. However, the superposition of phases corresponding to different tensors has increased the difficulties of applying this method directly. Thus it was necessary to perform a selection of the data in order to deduce information about the
660
phases. Further studies concerning the selection method are on the way. At the present time, this method is valuable because it leads to a new quantitative estimation of the various stress components. REFERENCES Carey, E., 1979. Recherche des directions principales de contrainte associees au jeu d’une population de failles. Rev. Geogr. Phys. Geol. Dyn., 21 (1): 57-66. Carter, N.L. and Raleigh, C.B., 1969. Principal stress directions from plastic flow in crystals. Geol. Sot. Am. Bull. 80: 1231-1264. Dantzig, G.B., 1963. Linear Programming and Extension. Princeton Univ. Press, Princeton N.J., 341 pp. Friedman, M., 1964. Petrofabric techniques for the determination of principal stress directions in the rocks. In: W.R. Judd (Editor), State of Stress in the Earth’s Crust. Elsevier, New York, N.Y., pp. 450-550. Friedman, M., 1967. Description of rocks and rock masses with a view to their physical and mechanical behaviour. In: Proceedings of the First Congress of the International Society of Rock Mechanics. Lisboa, III: 182-197. Friedman, M. and Conger, F.B., 1964. Dynamic interpretation of calcite twin lamellae in a naturally deformed fossil. J. Geol., 72: 361-368. Friedman, M. and Heard, H.C., 1974. Principal stress ratios in cretaceous limestones from Texas Gulf Coast. Bull. Am. Assoc. Pet. Geol., 58 (1): 71-78. Friedman, M., Teufel, L.W. and Morse, J.D., 1976. Strain and stress analyses from calcite twin lamellae in experimental buckles and faulted drape-folds. Philos. Trans. R. Sot. London, Ser. A, 283: 87-107. Goldman, A.J., 1956. Resolution and separation theorems for polyhedral convex sets. In: H.W. Kuhn and A.W. Tucker (Editors), Linear Inequalities and Related Systems. Princton Univ. Press, Princeton, N.J., pp. 41-51. Groshong, R.H. Jr., 1974. Experimental test of least-square strain gauge calculation using twinned calcite. Geol. Sot. Am, Bull., 85: 1855-1864. Jamison, W.R. and Spang, J.H., 1976. Use of calcite twin lamellae to infer differential stress. Geol. Sot. Am. Bull., 87: 868-872. Nissen, H.U., 1964. Dynamic and kinematic analysis of deformed crinoid’s stems in a quartz graywacke. J. Geol., 72: 346-360. Turner, F.J., 1953. Nature and dynamic interpretation of deformation in calcite of three marbles. Am. J. Sci., 251: 276-298. Turner, F.J., Griggs, D.T. and Heard, H., 1954. Experimental deformation of calcite crystals. Bull. Geol. Sot. Am., 65: 883-934.