A NEW APPROACH TO MASS BALANCE MODELING:APPLICATIONS TO IGNEOUS PETROLOGYBYPEIWEN KEB.Sc., THE UNIVERSITY OF SCIENCE AND TECHNOLOGY OF CHINA, 1989A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTERS OF SCIENCEINTHE FACULTY OF GRADUATE STUDIESDEPARTMENT OF GEOLOGICAL SCIENCESWe accept this thesis as conformingto the required standardTHE UNIVERSITY OF BRITISH COLUMBIAAugust 1992.© Peiwen Ke, 1992.In presenting this thesis in partial fulfilment of the requirements for an advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be granted by the head of mydepartment or by his or her representatives. It is understood that copying orpublication of this thesis for financial gain shall not be allowed without my writtenpermission.(Signature)Department of Geological SciencesThe University of British ColumbiaVancouver, CanadaDate ^A v [ 2_6 , t ,dy 3 DE-6 (2/88)ABSTRACTAn overdetermined linear system of equations is proposed to investigate mass balancerelationships in igneous systems. A chi square (x 2) minimization technique is used to solve theset of linear equations. The x2 technique is superior to R 2 techniques because it considers theanalytical errors on measurement data. Singular Value Decomposition (SVD) is used to searchfor the x2 optimal solution.The attributes of the x2 technique are discussed. The x2 minimization fitting technique,combined with SVD methods allows the confidence regions to be calculated and expressed bya geometric shape centered on the X2 optimal solution. The geometric shape varies with thedimension of the linear system: one-dimensional cases yield a line fragment; two-dimensionalcases yield an ellipse and the three dimensions require an ellipsoid. The size of the confidenceregions is determined by the confidence limits and is related to the coefficient matrix of thelinear problem. Model with smaller confidence regions are considered to be better thanmodels associated with larger confidence regions where the same confidence limit is used.Therefore, the confidence region can be used to select the model that best explains chemicaldifferentiation in igneous system.Synthetic data sets have been analyzed using the new program MASSBAL.0 whichsolves the mass balance problems using the principles of the X2 technique and SVD. The X2technique presents a solution closer to the true one than does the traditional residualminimization fitting (R2 technique). It overcomes the deficiencies that the R2 technique hasand appears to be a more robust and sensitive fitting technique.IIThe x2 technique has been applied to the apparent clinopyroxene paradox formed in lavasfrom Diamond Craters, Oregon. Three models are proposed to explain the chemicaldifference: a three phase model (01+P1+Px), a two phase model (O1+Pl) and a one phasemodel (01). The X2 mass balance modeling results show that the three phase model(O1+Pl+Px) has the smallest minimum )(2 value and the most compact confidence region;thus it is preferred and clinopyroxene did participate in the fractionation process.The )(2 technique has also been applied to study the assimilation problem in Paricutinvolcano and the fractionation problem in PLA group of Mount Meager volcanic complex. Themass balance calculations show that assimilation has play a major role in Paricutin volcano andabout 15% xenolith has contaminated the magma. In PLA fractionation processes, plagioclaseand amphibole are shown to have caused the chemical difference between lavas, although otherphenocrysts ( biotite and pyroxene ) did not participate.IIITABLE OF CONTENTSABSTRACT^ IITABLE OF CONTENTS^ IVLIST OF TABLES VILIST OF FIGURES^ VIIIACKNOWLEDGMENTS XCHAPTER 1 INTRODUCTION AND PRINCIPLES OF MASS BALANCEMODELING^ 11.1 Objectives^ 21.2 Previous Mass Balance Modeling Approaches^ 31.3 Fitting Techniques Used in Literature 9CHAPTER 2 MASS BALANCE MODELING ALGORITHM^ 112.1 Mass Balance Calculations in Igneous Differentiation 112.2 Fractionation of Three Phases^ 122.3 Optimization of Mass Balance Equations^ 162.4 Mapped Confidence Region for the X2 Optimization^ 202.5 Goodness-of-Fit for Optimal Solutions^ 262.6 Computer Program and Calculation Procedures^ 26CHAPTER 3 SIMULATED DATA SETS AND ALGORITHM TESTING^ 313.1 Calculation of Simulated Data Set.^ 313.2 Simulated Data Sets with Associated Uncertainties. ^ 323.2.1 Solutions for Synthetic Data Sets with 1% Relative Error.^ 343.2.2 The Effect of a on X2 and R2 Optimization^ 393.2.3 x2 Confidence Regions: A Measure of Goodness-of-Fit and ATest of Hypothesis^ 45CHAPTER 4 NATURAL DATA SETS 50IV4.1 Paricutin Volcano: An Example of Assimilation^ 514.2 Diamond Craters Volcano: An Example of Fractionation ^ 544.3 Mount Meager: Another Example of Fractionation 574.2.1 Introduction^ 574.2.2 Mineralogy and Rock Chemistry Data of the Mount MeagerVolcanic Complex^ 614.2.3 Mass Balance Calculations with the Plinth Assemblage (PLA) ^ 694.2.3.1 Comagmatic Relationships^ 694.2.3.2 Mass Balance Calculations 72CHAPTER 5 SUMMARY AND DISCUSSION^ 78REFERENCES^ 83APPENDIX 88Al Linear Algebra Theorem of SVD^ 88A2 Detailed Petrographic Descriptions. 99A3 Whole-Rock Chemical Analysis^ 104A4 Electron Microprobe Analysis . 105A5 MASSBAL.0 source code and test data. ^ 120vLIST OF TABLESTable 1.1 Summary of mathematical developments and applications of massbalance modeling in literature. ^ 5Table 1.2 Symbol notation for review of literature (Table 1.1)^ 8Table 2.1 Symbol notation for development of mass balance algorithm.^ 13Table 2.2 Critical R2 values for acceptance of solutions used in literature. ^ 19Table 2.3 AX2 values as a function of confidence level and degrees of freedom^ 23Table 3. la Synthetic chemical data set without assigned errors ^ 33Table 3.1b Results of calculation with MASSBAL.0 on synthetic data sets withoutsimulated errors^ 33Table 3.2 Random errors and the synthetic data.^ 36Table 3.3 Results for synthetic data set No. 1. 37Table 3.4 X2 and R2 fitting results for synthetic data sets involving three phasefractionation.^ 40Table 3.5 x2 results of synthetic data sets with inaccurate errors. ^ 43Table 3.6 x2 results for the examples with incorrect models. 46Table 4.1 Chemical data for assimilation problem. ^ 52Table 4.2a Comparison of solutions from MASSBAL.0 and the literature. ^ 53Table 4.2b Confidence regions around the X2 solutions on the assimilationproblem.^ 53Table 4.3a Chemical Data for Diamond Crater's problem.^ 58Table 4.3b The X2 optimal solutions and the confidence regions on around the X2optimal solutions for Diamond Crater's problem.^ 58Table 4.4 Phenocryst, microphenocryst and groundmass mineralogy in MountMeager rocks^ 62VITable 4.5 Whole rock chemistry of the Mount Meager volcanic rocks. ^ 65Table 4.6 Estimates of XRF machine precision and sample preparation uncertaintyfor sample MM-32.^ 68Table 4.7 Chemistry of minerals from PLA^ 70Table 4.8 Chemical data for PLA fractionation problems^ 73Table 4.9 Mass balance calculation results for PLA. 75Table A3.1 XRF operating conditions for the analysis of major and minorelements. ^ 106Table A3.2 XRF operation conditions for the analysis of trace elements.^ 107Table A3.3 The effect of Fe2+/Fe3+ ratio on CIPW norms.^ 108Table A4.1 Electron microprobe analyses of PLA plagioclases. 109Table A4.2 Electron microprobe analyses of PLA amphiboles .^ 114Table A4.3 Electron microprobe analyses of PLA pyroxenes. 115Table A4.4 Electron microprobe analyses of PLA biotites^ 116VIILIST OF FIGURES Figure 2.1 X2 probability density function^ 22Figure 2.2 Two dimensional solution spaces defined for two levels of confidence^ 25Figure 2.3 Comparison of the goodness-of-fit for different models. ^ 27Figure 2.4 Flowchart of the calculation procedures of the program MASSBAL.C. ^ 28Figure 3.1 Procedures to generate error-assigned synthetic data.^ 35Figure 3.2 Confidence region around the best solution for error-assigned syntheticexample 1^ 38Figure 3.3 Influence of standard deviation on x2 and R2 optimization .^ 41Figure 3.4 The effect of error estimation on the x2 technique optimization.^ 44Figure 3.5 Confidence regions on solution for model with "extra" phase. ^ 48Figure 3.6 Comparison between (false) 2-D confidence regions and projections of(true) 3-D confidence regions. ^ 49Figure 4.1 Two dimensional confidence regions around optimal solutions forassimilation problem (late stage).^ 55Figure 4.2 Three projections of the 3-D confidence regions around optimalsolutions for assimilation problem (middle stage). ^ 56Figure 4.3 Comparison between confidence regions for two phase model and theprojections on the corresponding plane for the three phase model for the DiamondCrater's differentiation problem ^ 59Figure 4.4 The location of the MM volcanic complex with the Garibaldi belt ofsouthwestern British Columbia^ 60Figure 4.5 The element ratios P/K to test the comagmatic relation in PLA rocks^ 71Figure 4.6 The projections of three dimensional confidence regions around theoptimal solution for the PLA early stage of chemical differentiation.^ 76VIIIFigure 4.7 The projections of three dimensional confidence regions around theoptimal solution for the PLA late stage of chemical differentiation. ^ 77IXACKNOWLEDGMENTSThis thesis is completed under the supervision of Dr. J. K. Russell, to whom Ihave much appreciation, for he has been continuously helping and discussing with meduring the three years program. His carefully reviewing different drafts of this thesis hasimproved the manuscript dramatically. Especially, his enthusiasm to science hasencouraged me through all the hard times in the research.I also wish to thank Dr. Tom Brown and K. Fletcher for many good suggestionsand comments during both the research and writing periods. I am grateful to the helpfuldiscussions with Bruce James, Yao Cui, Tom Clemo, Xiaolin Cheng, Bo Liang and MinSun. In addition, I wish to express my appreciation to Stanya for her assistance in majorand trace element measurements.NSERC operating Grant A0820 provided funding for this project.xCHAPTER 1 INTRODUCTION AND PRINCIPLES OF MASS BALANCE MODELINGThis thesis includes the mathematical development of chi square minimization fittingtechnique (x2 technique) in mass balance modeling and test runs of program MASSBAL.0using synthetic data sets, in addition to natural data sets. This chapter introduces previousapproaches of mass balance calculations in geology and the fitting techniques used in literature.In Chapter 2, a group of linear equations is developed to express mass balance relationsin three phase fractionation process with the amounts of crystallizing mineral phases as theunknowns. Both the residual square minimization fitting technique (R 2 technique) and the x2technique are discussed, including the limitations of the former and the attributes of the latter.Program MASSBAL.0 was developed to perform all the mass balance calculations, based onthe Singular Value Decomposition (SVD), in addition to X2 and R2 minimization principles.This program was written in C using a UNIX operating system. It provides two results: i) theoptimal solution for the mass balance equations and ii) the confidence regions around theoptimal solution.Besides the theoretical derivations, Chapter 3 discusses test runs for programMASSBAL.0 that have been performed to illustrate that the x2 technique provides morereliable solutions than does the R2 technique. Synthetic examples are used to carry out the testruns, since the expected solutions are known for synthetic data sets. Three groups of syntheticdata are involved. First, a synthetic data set without assigned errors is used to test thereliability of the program MASSBAL.C. The second group of data sets is synthetic data withassigned errors. This analysis is used to study that the influence of the standard deviation onthe x2 and the R2 optimizations when the data are subject to measurement errors. The lastgroup of synthetic runs involves incorrect hypotheses (adding the unwarranted phases in thesolution and dropping the essential phases from the solution) in order to test the goodness-of-fitof the x2 technique.In addition to the synthetic examples, three natural data sets were also applied todemonstrate that the X2 technique provides the better solutions than does the R 2 technique.These natural data sets are: assimilation problem in Paricutin volcano, Diamond Crater'sfractionation problem and the fractionation processes in Plinth assemblage of Mount Meagervolcanic complex.1.1 ObjectivesThis thesis is about modeling of geochemical data. The process of modeling data usuallyincludes three steps. The first step involves observation of geological phenomenon,recognition and definition of the problem. Secondly, a model is proposed to explain thephenomenon. The model commonly involves a number of parameters that can be relatedthrough mathematical equations. Thirdly, the parameters are fitted by finding the bestagreement between the observed data and the values predicted by the model (Press et al.,1987). The array of fitted parameters is termed the "best solution" or "optimal solution".The central part of this thesis emphasizes the third step and utilizes the x2 technique inmass balance modeling of igneous differentiation. The X2 technique represents a significantimprovement over the traditional R2 technique.Mass, balance modeling has been applied to interpretations of geochemical data since theearly 20th century (Bowen, 1928). It is still of importance in quantifying a variety ofgeochemical processes, especially where other approaches (e.g., thermodynamic modeling) are2not possible. Thus, developing a robust and sensitive fitting technique for mass balancemodeling remains a significant goal.Three criteria are used to judge a fitting technique. First, the technique must be able toidentify a "best solution" for the fitted parameters. If the observed data are subject tomeasurement errors, the measurement errors should be considered during the search for thebest solution. Secondly, the technique must provide a way of evaluating the goodness-of-fitagainst some useful statistical standards (Press et al., 1987). Lastly, it should also provide away to estimate the ranges of solution that share similar confidence levels. The optimalsolution is usually not the only possible answer since we are working on complicatedgeological systems.The R2 technique commonly used to search for the optimal solution does not satisfy allthree criteria listed above. The specific limitations to the R 2 technique are introduced in thenext section and discussed in detail in Chapter 2. The x2 technique developed in this thesiscan overcome the deficiencies of the R 2 technique and satisfies the three essential attributes ofa mathematical fitting technique.1.2 Previous Mass Balance Modeling ApproachesMass balance has been used to solve different geological problems. Metamorphicpetrologists employ mass balance and mass transfer concepts to the study of element mobility.They are interested in quantifying the gains and losses of mobile elements duringmetasomatism and other alteration processes (Myers et al., 1986). Igneous petrologists alsouse this universal concept to interpret the major element chemical patterns resulting frommagmatic differentiation processes, such as fractionation (Nesbitt & Cramer, 1981; Gerlach &Luth, 1980), magma mixing (Wright & Doherty, 1970), partial melting and anatexis (Bea,31989), and assimilation (DePaolo, 1981). Typical problems involve well-defined compositionsof the magma end members that can be related by a single process. The mass balancemodeling is a reasonable method by which the phases participating in the differentiationprocess are identified and their amounts are calculated. In the last three decades, there havebeen approximately 17 papers that have presented either new applications or new mathematicaldevelopments of mass balance modeling. The summary of applications, mass balanceformulas and fitting techniques in these papers are presented in Table 1.1 with symbol notationin Table 1.2.Previous applications of mass balance modeling in petrology can be divided into twomain groups. The first group of applications includes studies of open geological systemprocesses, such as metasomatic alteration (e.g., Gresens, 1966). The other is the study ofclosed system processes such as igneous differentiation within a single magmatic system. Byclosed we mean that material is conserved within the considered object. Fractionation, forexample, can be treated as closed system differentiation as long as the compositions of twomagma end members and involved mineral phases are known (Stormer & Nicholls, 1977;Wright & Doherty, 1970).Whereas it is straightforward to apply mass balance to closed system processes, opensystems present more difficulty. Where the system is open, we expect mass transfer of at leastthe mobile elements into or out of the system. However, even in open systems some elementscan be recognized as immobile or conserved and separated from those elements which aresubject to mass transfer. We consider the systems to be closed for these immobile elements.Mass balance relations of these immobile elements can be used as constraints to calculate themass gains or losses of the active elements. Therefore, mass balance modeling can still beapplied to study open geological system processes. Of course, if no immobile element can befound, the trace of the open system process is very hard to track.4Table 1.1 Summary of mathematical developments and applications of mass balance modeling in literature.Authors^Applications^ Formulas^ Fitting TechniquesGresens(1966)Chakraborty(1977)Gerlach & Luth(1982)Estimate element gains and lossesin metasomatic alterations.Estimate element gains and losses incore-mantle-matrix type metasomaticdifferentiation.Evaluation of the clinopyroxeneparadox.Implications for the origin of Hi-alumina arc basalt.Estimate element gains and lossesin metasomatic alteration, anothersolution for Gresens's equation.PDeX i = 100 [ f,(—) d i - p iPpP•pi =Dc •dci + Dm •dn:P = D c + Dm(Mp i =100(D c •dci + Dm •dmi - P•p i)/P•pino formulasno formulaseX i = [i; di - p i]Pn- 1 ,p i = E o i .yi + yn • pii= 1Francis(1986)Brophy(1986)Grant(1986)Wright & Doherty Computation for petrologic mixing(1970)^problems.Interpretation of differentiation model no formulasfor the 1959 Kilauea Iki.Estimates are based on volume ratios calculated with oneimmobile element, no 4nalytical errors are considered.Estimates are based on mass ratio of core and mantlecalculated from a pair f immobile elements. Massratios from different pa rs of immobile elements areassumed to be consistent, no consideration of analyticalor hypothetical errors.Unweighted least squars fitting.Unweighted least squard fitting, the solution with theR2 less than 0.01 is accepted.Unweighted least squart fitting, the critical R 2 valueis 0.1.Graphical solution, cu es of components in originalrock versus altered roc are made. Element gains andlosses are estimated by omparing curves of mobileselement with curve of immobile element or best fitof a straight line throug a series of points of immobileelements.Unweighted least squarn fitting, the accepted solutionsmust be non-negative and with the minimum R2.Bea^Partial melting and anatexis mixing by B=a0A.- aiAi(1989) mass balance among source, i=1S2 = It' (Si )2i= 1segregate and residue.Table 1.1 continuedAuthors^Applications^ Formulas^ Fitting TechniquesStout & Nicholls(1977)Michael & Chase(1987)no formulasEvaluation of REE mobility inlow-grade metabasalts.Identification of two scales ofdifferential element mobility inductile fault zone.Interpretation of fractionalcrystallization model for petrogenesisof the Snake River Plain Basalts.Studying the influence of primarymagma composition, H20 andpressure on Mid-Ocean Ridge basaltdifferentiation.Average volume ratios were estimated by assuming allREEs are immobile, gains and losses of REE werecalculated, one whose gain or loss more thananalytical error was considered to be mobile elementand taken off, new average ratio was calculated againuntil each REE's gain or loss is with analytical error,they were regarded as immobile. Gains and losses ofother mobile REEs were estimated by the averagevolume ratio of mobile REEs.Immobile elements were used to constrain and calculatethe gains and losses of mobile elements.Unweighted least squares fitting. The solution forThwhich each element satisfies S.< • was considered asthe best fit.Unweighted least squares fitting, acceptable solutionsare those with R2 value less than 0.10 and whichmakes residuals of each element smaller than analyticaluncertainty.Bartley(1986)Dipple et al.(1990)eX i = P(fvfpd i - p i)r PDeX i =100[fv(—d i - p i)]Ppn b i.M.^d i Unweighted least squares fitting.E ^ – — (100 +^- p ij^100^100 j=1--'Albarede & Provost Determinations of mode of lunar basalt D.. 0 .L u •j = 1Chi-square minimization fitting.(1977) and stoichiometric coefficients forcoronitic mineral reactions.Table 1.1 continuedAuthors^Applications^ Formulas^ Fitting TechniquesStormer & Nicholls Interactive testing of magmatic(1978)^differentiation models.Mj (b ij -d i)/100^d i -p ii= 1Unweighted least squares fitting. Solutions providingthe R2 less than 2.0 and with residual on each oxideless than analytical error are considered better. Thismethod is regarded to work better with the number ofcomponents at least as twice big as the number ofphases.Stout & Nicholls(1983)Study of origin of the hawaiites. no formula. Unweighted least squared fitting. The solution whichhas residual squares less than 2.0 is not rejected; noconclusion can be reached for the solution withresidual between 2.0 and 3.0 ; The solution with theresidual squares more than 3.0 is rejected.Myers et al.(1986)Determination of quartz eclogite as asource mixing with carbonate andpelagic sediment as parental Aleutianmagmas.p i = X.dei^+(1-X-Y)•dsi Graphic solution, Cpm-C i curves with different X weremade. The X value which can satisfy the variability ofall major, trace and rare element was considered as thebest fit.Myers et al.(1986/87)Determination of mixing proportionsof two and three components mixingproblems with end membercompositional variability.d.d i = Same as above.2^Id I dTable 1.2 Symbol notation for the summary of literature in Table 1.1.Symbols^NotationsP^Mass of parental or original rockD Mass of daughter or altered rockDm^Mass of mantle rockDo^Mass of matrix rockM. Masses of mineral phasesPi^Weight fraction of the ith oxide component in parental or original rockd i^Weight fraction of the ith oxide component in daughter or altered rockbij^Weight fraction of the ith oxide in the jth phasey i^Weight fraction of the jth mineral in source rockeX i^Element gains or losses from 100 or P grams parental or original rockd Weight fraction of the ith oxide component hybrid magmadm^Weight fraction of the ith oxide component mantle rockWeight fraction of the ith oxide component matrix rockde^Weight fraction of the ith oxide component quartz eclogitedei Weight fraction of the ith oxide component carbonated:^Weight fraction of the ith oxide component sedimentd. Weight fraction of the ith oxide in mixing component 1d2^Weight fraction of the ith oxide in mixing component 2PD^ Density of daughter or altered rockPp Density of parental or orginal rockfP^Density ratio between altered and original rocksfv Volume ratio between altered and original rocksX^Mixing proportion of quartz eclogite or parental magmaY Mixing proportion of carbonate or component 1The ratio of the i th element in between the source rock and the liquid.D id^The partition coefficient of the i th element in the j the phase.(Di The i th phase mass fraction in the whole rock.Ao^Compositional vector of source rockB Compositional vector of segregateA. Compositional vector of source minerals (j=1,2,...,n)Si^Residual square of oxide component iS Sum of residual squaresMI^Magnitude of imbalance1.3 Fitting Techniques Used in LiteratureA variety of fitting techniques have been applied to determine the parameters in massbalance modeling. These techniques include conserved element fitting (Bartley, 1986), fittingby minimization of sums of residual squares (Stormer & Nicholls, 1977), and fitting byminimization of chi square (Albarede & Provost, 1977).In order to calculate the element gains and losses between the altered rock and originalrock, the volume or mass ratio must be known. Mass balance equations of immobile elementswere applied to determine the volume (Bartley, 1986) or mass (Grant, 1986) ratio. Thisprocess, which includes: finding immobile elements, estimating the volume or mass ratiobetween, before, and after the alteration, and calculating element gains and losses, is definedto as conserved element fitting.Although the conserved element fitting technique is simple to use, it can be biased.There is less problem if there is only one immobile element during the alteration process.When more than one immobile element exists, it has been assumed that there is a singlevolume or mass change ratio. Generally, one immobile element is used to calculate thisvolume or mass change ratio (Hanor & Duchac, 1990). Actually, the volume ratios calculatedfrom different immobile elements can be very different due to analytical and hypotheticalerrors. To avoid bias, every immobile element must be considered during the volume ratioestimation. Actually, least square fitting could probably be used to estimate this ratio using allthe immobile elements simultaneously.Igneous petrologists are more interested in closed system reactions. The parameters ofinterest are different for different magmatic processes, but commonly include the amounts of9mineral phases involved in the fractionation process (Stout & Nicholls, 1977) or theproportions of end members for assimilation (DePaolo, 1981) and magma mixing (Wright &Doherty, 1970) processes. Assimilation and magma mixing can be treated as closed systemprocess, because the end-member chemical compositions are known, as is the composition ofthe added material.Mass balance equations for each oxide component are used as constraints to calculatethese parameters. Usually there are at least nine such mass balance equations. In mostigneous applications the number of unknown parameters is less than the number of equations.Thus, mass balance equations for differentiation processes define an overdetermined linearsystem of equations. Commonly, the R 2 technique has been used to find the "best solution"(e.g., Stormer & Nicholls, 1977).As compared with arbitrarily choosing a certain number of equations to estimate theparameters (Olsen, 1984), the R2 technique is a major improvement, however, it still haslimitations. First, it does not consider analytical errors. Secondly, when a solution isachieved by minimization R 2 , there is no criterion for establishing the goodness-of-fit. Last, itdoes not provide a way to estimate the confidence region around the optimal solution. Thus, anew fitting technique needs to be developed.The X2 technique developed in this thesis can be dated back to Albarede & Provost's(1977) first study. However, in their work, they mostly emphasize on theory and did notconsider the confidence region around the optimal solution. In this thesis, the applications ofX2 technique have been emphasized. The confidence regions around the optimal solutions havebeen applied to test different models for three natural data sets.10CHAPTER 2 MASS BALANCE MODELING ALGORITHMThis chapter introduces two fitting techniques for mass balance modeling of chemical data.The R2 technique is reviewed and the x2 technique is developed. The X2 is defined as the R2weighted by the variance (Press et al., 1987). The optimizations of these two techniques aresimilar where the variances for each oxide are equal. The advantages of the X2 technique overthe R2 technique lie in its consideration of analytical errors. The chi-square probabilityfunction follows the incomplete gamma function Q(x2 1 v), which makes it possible to estimatethe goodness-of-fit of the x2 optimization.The optimal solution to the x2 minimization problem is solved by the Singular ValueDecomposition (SVD). SVD is a powerful technique for solving a linear problem withsingular or very close to singular matrices. The standard deviations on fit parameters aremutual independent and can be expressed by the confidence regions around the optimalsolution that share similar significance (Press et al., 1987). SVD provides a simple way tocalculate the confidence regions.2.1 Mass Balance Relations in Igneous DifferentiationSuppose that two igneous rocks are related through magmatic differentiation. The parentmagma is represented by one rock composition pi and the daughter system is represented bythe other rock composition d i . The compositions of separating phases are also assumed to beknown. In fractional crystallization, for example, compositions of the involved mineral phasesmay be provided by the analysis of phenocrysts of these two rocks. The following equationexpresses the mass balance relation between the parent and the daughter magmas. The mass ofthe parent magma (P) is equal to the mass of the daughter magma (D) plus the mass(es) of theparticipating phases (Mi ).11P = D + mi .j=1(1)P, D and M. in equation (1) represent different phases for different processes. Firstly,for fractionation processes, P and D are masses of two cogenetic rocks. Mi are masses ofmineral phases responsible for the change in the system's compositions (Stormer & Nicholls,1977). Secondly, for assimilation processes, P and D are similarly masses of two cogeneticrocks, whereas Mi are masses of wallrock components (DePaolo, 1981). Thirdly, for magmamixing processes, Mi and D represent masses of end-member magmas (Wright & Doherty1970). Finally, equation (1) can also be used to describe partial melting and anatexis mixingprocesses. In this case, P, D and Mi indicate source, residual and segregates, respectively(Bea, 1989).For all the geological problems being studied, geologists are concerned about the amountsof phases (M.) in order to explain the chemical difference between two magmas. Mass balancecalculations are a means to estimate these amounts. A three-phase fractionation case is used asan example to illustrate how mass balance calculations help to solve geologic problems.Symbols used in mathematical development are listed in Table 2.1.2.2 Fractionation of Three PhasesTwo rocks are assumed to be related through fractional crystallization; from petrographicexamination, mineral phases Mi (j =1,2,...,m) are considered to have been involved. Thequestion is whether all of the mineral phases or any of the combinations of them can explainthe compositional differences between the cogenetic rocks. If so, we also want to know howmuch of each mineral is needed. In order to answer these two questions, the compositions ofthe parent and the daughter magmas, as well as the compositions of minerals must be known.12Table 2.1 Symbol notations used in development of mass balance modelimz alLwrithm.Symbols^NotationsA unweighted compositional matrixA' [—aij ]9x3, weighted compositional matrixaiunweighted compositional arrayB' bi[—] 97, 1 , weighted compositional arrayaiX^[Mi]3xi, mass array of participating mineral phasesaii^mii - d i , chemical difference between the j th mineral and daughter maumabi^pi - d i , chemical difference between parental and daughter mamasM. Mass of fractional crystallized mineral phases8mj^ranges of weight of fractional crystallized mineral phasesP Mass of parental magmaD Mass of daughter magmaPi^Weight fractions of the i th oxide component in parental mamadi^Weight fractions of the i th oxide component in dauuhtermii Weight fractions of the i th oxide component in the j th mineral phasex2^Summed chi squared of nine oxide componentsR2^Summed residual squared of nine oxide componentsRR Relative error6.^Standard deviation of analytical error for the i th oxide componentQ(X21 V)^x- probability functionQ(01 V) Limiting value of Q(X 2 1 V) functionQ(00 I V)^Limiting value of Q(X 2 I V) function^ i - j, degree of freedomAX2^Difference between certain chi square value and the minimum chi squared valuep Probability of ex2 being less than a certain valueU & U'^ixj column-orthouonal matrixVT & V'T^Transposes of an jxj orthogonal matrix V and V'W & W'^jxj diagonal matrix with positive or zero elementNP (i)^Column vectors of V'Diagonal elements of W'13Usually, the compositions of rocks are determined by X-ray fluorescence analysis and arereported as weight percents of nine oxides Si0 2 , Ti02 , Al203 , FeO, MgO, CaO, Na20, K20and P205 . The analytical errors of the XRF results are also estimated during the analysis. Thecompositions of the minerals Mi are given by electron microprobe. Equation (1) can bewritten in terms of the individual oxide constituent:P•pi = D•di + il [ Mi•mii ],^ (2)j= 1for i=1,2,...,n, where- pi are weight fractions of oxide component i of the parent composition,- di are weight fractions of oxide component i of the daughter composition, and- mii are weight fractions of oxide component i of mineral phase j.If the mass of parent system P is defined as 1 mass unit, (1) implies:D = 1- ti mj .j=1Equation (3) together with equation (2) produces:rii ivy m i; _ d i) = pi - di.j= 1Equation (4) is the final mass balance equation for fractionation. It can be written in anabbreviated form as (e.g., Stout & Nicholls, 1977; Stormer & Nicholls, 1977):ril (Miaii) = bi.^ (5)j=1(3)(4)14where- aij are the difference of the weight fractions of oxide component i betweenthe mineral phase j, mij, and daughter's, d i ,- b i is the weight fraction difference between parent and daughter magmas,- n is the number of oxide components, and- m is the number of involved mineral phases.Commonly there are at least nine oxide components in the system that provide nine massbalance equations similar to equation (4). The following group of linear equations expressesthe three phase fractionation.all^a12 a13'.\a21 a22 a23a31 a32 a33a41 a42 a43a51 a52 a53a61 a62 a63a71 a72a73a81 a82 a8391 a92 a93/(M1\M2 =\M3)(6)In matrix notation:A * X = B,^ (7)where A =[aii ]9x3 , X = [MORI, and B=Ibiii9 x i.The matrix A and array B comprise the composition vectors of the two magmas (pi andd) and the involved mineral phases (m ij). B is the measured chemical difference pi - dibetween the two rock compositions representing the two magmas. The array X is the unknownvector representing the amounts of minerals M j participating in the fractionation. The left-15hand sides of (6) represent the fractionation model subject to adjustable parameters M. Theproblem is to fit parameters Mi so that the best agreement between the observed data and themodel can be obtained.There are a larger number of equations than unknowns in the linear set of equations (6),thus, it becomes an overdetermined linear problem and can not be solved explicitly. The "bestsolution" is an approximation that can be obtained by some fitting techniques.2.3 Optimization of Mass Balance EquationsThe most common approach that has been used to find the optimal solution is theunweighted R2 minimization fitting technique (e.g., Stout & Nicholls, 1977). It has beenapplied to geology for at least two decades (e.g., Wright & Doherty, 1970).The principle of the R2 technique is to search for the optimal solution which produces theleast R2 value (R2min), where R2 (sum of squares of residual) is given by:9^3R2 = E (bi - E (aijMi))2i =1^j =1(8)The R2 technique represents a quantitative tool in the interpretation of geochemical data,however, it has several limitations. One of its limitations is that the R 2 function does notinvolve any parameters representing analytical uncertainty. Analytical errors associated withXRF analysis can be as high as 2% relative error for major elements and 8% for traceelements. These errors influence the compositional matrix A and array B in equation (6).Thus, they can affect the "optimal solution" found by R 2 optimization.16Recognizing the effects of analytical errors on R2 optimization, previous workers havetried to solve this problem. In Bea's application of mass balance to partial melting andanatexis processes (1989), he invented a parameter Thj which was related to analytical errors.The solutions which could make the residual for each oxide component less than Thj wereconsidered acceptable. In Stormer & Nicholls's (1977) early mass balance modeling, theysimply replaced Thj by analytical errors in their mass balance modeling. They considered thesolutions which make the residual for each oxide component less than analytical error werebetter.All these previous trials were good. However, they did not solve the fundamentalproblem; the measurement errors should be considered during the optimization procedure. x2minimization techniques are designed to handle variances associated with measurement errorand are therefore an improvement over the traditional R 2 minimization techniques.The x2 function can be expressed as (9). The principle of the )(2 technique is to searchthe optimal solution which can provide the least )(2 value (x2iiii.).7C2 =^[ (bi -^(aijMj))2 /05i2^(9)i = 1^j=1where- )(2 is the sum of the chi squares of each equations, and- ai is the standard deviation of measurement error of pi.The )(2 technique is related to the R2 technique. Its formula is R2 weighted by thevariance ail. The x2 technique actually deals with a new linear combination (10):A' • X ---- B',^ (10)17r aij^ r bijwhere , = 1.-01 19x3, X = [Mil3x1, and B ., =L— i9xi.One of the advantages the X2 technique has over R2 technique is that the solution foundby X2 is more sensitive to the constraints imposed by the precise data than to imprecise data.The x2 technique provides a way to treat good and bad analyses differently. The bigger theanalytical uncertainty associated with the oxide component measurement, the smaller itscontribution to the x2 function (9). Thus, the optimal solution from X2 technique appears to bemore reliable solution than the R2 optimal solution.In order to get a good solution from the X2 technique fitting, analytical error must bemeasured precisely. Both overestimation and underestimation for analytical errors influencethe X2 optimization adversely. The overestimation of error reduces the effectiveness of thedata and will make it very difficult to recognize bad hypotheses. Conversely, theunderestimation of error gives the data more weight than it deserves and will cause X2solutions to tend towards the unweighted R 2 solution. Both poor estimates of error willproduce poor solutions.Not only analytical uncertainty influences the optimal solution, but also errors associatedwith hypothesis. For example, the number of mineral phases used as the fractionationassemblage can cause a problem. If we keep adding participating mineral phases in theequation (8), we can always gain a smaller and smaller R 2 value. When should we stopadding mineral phases and accept the solution? Different geologist set their own criteria to dothis. Their philosophy is that only those solutions which produce a R 2 value less than somecritical R2 value can be accepted. Table 2.2 gives the critical R 2 values used by different18Table 2.2 Critical R2 values for acceptance of solutions.AuthorsR2Critical valuesFrancis (1986) 0.01Brophy (1986) 0.1Michael & chase (1987) 0.1Stormer &Nicholls (1978) 2.0Olsen (1984) 2.0Stout & Nicholls (1977) 2.0Bea (1989) Thi ** Thi is a parameter which is subject to analytical error (Bea 1989).19geologists. Compared with these arbitrarily set values, the X2 technique has another advantagediscussed below.The chi-square probability function P(x 2 1v) is defined as the probability that the observedchi-square for a correct model should be less than a value X2 (Press et al., 1987). Q(x 2 I v) isits complement and currently used in this thesis. This theorem makes it possible to evaluatethe validity of the X2 optimal solution, including the calculation of the confidence regionsaround the optimal solution and the goodness-of-fit (Press et al., 1987).2.4 Mapped Confidence Regions for x2 OptimizationIgneous systems are very complicated and errors inevitably exist in the analytical data weuse. Therefore, the solution selected by any optimization techniques should not be regarded asthe sole answer for the geological phenomenon we observe. There can be many othersolutions which are reasonable estimates of the true solution. We use the confidence regionaround the best solution to represent all the solutions with the same confidence level. Regionsdefined by the "critical" confidence level should be considered to contain equally validsolutions.Optimization provides the minimum value of X2 at a particular set of fit parameters.Obviously, the solutions around this point (x 2r in) has larger x2 values. Let Ax2 represent thedifference between x2min and the x2 value of the observed solution. The 0x 2 value willincrease with the observed solution perturbing further away from the optimal solution (Press etal., 1987). The confidence region is the zone around the optimal solution and with a constant0x2 value as the boundary (Press et al., 1987). The larger the 0x2 value is, the larger theconfidence region, and the larger the acceptable solution space (e.g. amounts of mineralphases) becomes.20Two factors comprise the confidence region. One is the geometric shape, for example,one dimensional problems involving a single fit parameter have confidence regions which areline fragments. Two dimensional problems involving two fit parameters have confidenceregions which are ellipses. Three dimensional problems involving three fit parameters haveconfidence regions which are ellipsoids. Higher dimensional problems involving highernumbers of fit parameters have confidence regions which are polygonal geometry. Theconfidence regions are centered on the optimal solution. The second characteristic of theconfidence region is the size which reflects the confidence level corresponding to a constantvalue of Ax2 . Three methods can be used to obtain the Ax2 value corresponding to a desiredconfidence level. These methods are:1) Estimating the delta-chi-square probability density function curves, which areexpressed in Figure 2.1;2) Looking up the Ax2 distribution tables, which is given as an example in Table 2.3;3) Computing by program GAMMQ.0 (Press et al., 1987). The source code of programGAMMQ.0 is included in the disk attached at the back of this thesis. This method isintroduced below.The probability distribution of the variable Ax2 follows the incomplete gammadistribution (Bury, 1975). Gammq function Q(Ax 2 1v) is defined as the probability that theobserved Ax2 will exceed the Ax 2 for a correct model (Press et al., 1987).gammq(v/2, Ax2/2) = 1-p,^ (12)where- p is the probability or the desired confidence level.2112 16 20 24 28 32 36A X 20^4Figure 2.1 The Chi square probability density functions, for several degrees of freedom.22Table 2.3 dx2 values as a function of confidence level and de2rees of freedom.Vp 1 2 3 4 5 668.3% 1.00 2.30 3.53 4.72 5.89 7.0490% 2.71 4.61 6.25 7.78 9.24 10.695.4% 4.00 6.17 8.02 9.70 11.3 12.899% 6.63 9.21 11.3 13.3 15.1 16.899.73% 9.00 11.8 14.2 16.3 18.2 20.199.99% 15.1 18.4 21.1 23.5 25.7 27.823This gamma function (14) has limiting values:Q(0 1 v) = 1,^Q(...1v) = 0,^(13)where- v is the degrees of freedom, which is equal to the number ofequations (one equation for each oxide) minus the number of fitparameters (mineral phases in fractionation modeling), v =n-m.The critical confidence level used in this thesis is 95 %. With the 0x2 estimate as theboundary, the size of the confidence region can be determined. SVD provides an easier wayto estimate the confidence region around the optimal solution. The linear algebra theorem,that SVD is based on, is presented in Appendix Al. For the three phase fractionation case, thesize of ellipsoidal confidence region is related to the parameters of A's decomposed matricesW and V. The boundary of the ellipsoid is given by:k2 = w112(v(i) • 8N4 1)2 + w222(v(2) • 6M2)2 + w332(v(3) • om3)2^(14)where w 11 , w22 and w33 are the diagonal elements of W, V(1) ,V(2) and V(3) are thecolumn vectors of V, 81\4 1 , 61■42 and 3M3 are the ranges of the weight of three mineral phases.While Ax2 is the boundary of the confidence region, the vectors V (i) are unit vectorsalong the principal axes of the confidence region ellipsoid. The semi-axes have lengths equalto I Ax I /w i • The confidence region for the two dimensional case is illustrated by Figure 2.2.24Figure 2.2 Sketches of two dimensional solution spaces defined by two different levels ofconfidence. The regions are centered on the optimal solution (double circle) and the largerregion is associated with higher confidence. Major and minor axes are parallel to unitvectors V(i) and V(2) respectively, and have magnitudes of a/w i and a/w2 for largerellipse, b/w i and b/w2 for smaller ellipse.252.5 Goodness-of-Fit for Optimal SolutionsIn the study of magmatic differentiation processes, different hypotheses can be used toexplain the geological phenomena. To study fractionation, for instance, we can form eithertwo-phase or three-phase fractional crystallization models to explain chemical differencesbetween parent and daughter magmas. The confidence region can be used to analyze thegoodness-of-fit of the x2 optimal solution, and further test the different models. Optimalsolutions with more compact 95 % confidence regions are regarded as better solutions. Figure2.3 shows that two different models provide different confidence regions. Model b's 95%confidence region is more compact, thus, the optimal solution of model b is considered betterthan that of model a (Press et al., 1987).The zero line of one mineral phase is defined as the line on which the amount of thismineral phase involved in is zero (Figure 2.3). If the 95% confidence region includes the zeroline of one mineral phase, the solution with zero amount of this mineral phase, is included inthe 95% confidence region, and it is possible that this mineral phase did not participate in themagmatic differentiation process.2.6 Program MASSBAL.0 and Calculation ProceduresProgram MASSBAL.0 was written to search for the fitted parameters of igneousdifferentiation models which purport to explain the chemical differences between two rocks.For example, it calculates the amounts of fractional mineral phases. It does two majorcalculations. First, the best solutions are calculated using both X2 and R2 techniques.Secondly, it calculates the confidence regions around the optimal solutions. The calculations26Figure 2.3 The different solution spaces used to compare the goodness-of-fit for differentmodels. The double circles represent the optimal solutions and are centered to theconfidence region ellipses. The dashed lines are zero lines (see text).27are based on the compositions of the two magma end members and the mineral phases and theprocedures can be described as the following steps:1. Calculate the compositional matrix A =[a ii]9x3 , and the compositional arrayB=[bi]9xi;2. Convert the relative errors into the absolute values, ai = RR x pi; if RR=O, 6i=1;3. Calculate the weighted compositional matrix A' =f aij--19x3 , and the weighteddicompositional array B''[ bi ]9x1;4. Decompose the weighted compositional matrix A' by Singular Value Decomposition(SVD), A' =1.P.WWIT;5. Calculate the fitted parameters by the X2 technique;6. Calculate the minimum X2 value (x2inin);7. Calculate the significance level of the X2 optimization through gammp function;8. Decompose the compositional matrix A by SVD, A =U-W•VT;9. Optimize the fitted parameters by the R2 technique;10. Calculate the minimum R2 value (R21nin);11. Calculate 0x2 value which is corresponding to certain probability p, such as 95% bygammp function;12. Calculate the orientations and lengths of the confidence region.A flowchart in Figure 2.4 is used to illustrate the calculation procedures of the programMASSBAL.C.28compositions ofmineral phasescompositions ofmagma end members1XRF analytical errorsmatrix A — — array BSVDweighted matrix A'U, V, Wgammq functionweighted array B'X 2II^optimalsolutionsiconfidenceregion aroundoptimal solutionFigure 2.4 Flowchart of the calculation procedures of the program MASSBAL.C.29This program was written in UNIX based C. The SVD and x2 routines are from Press etal. (1987). The source code, the format of input and output are given in the Appendix. Theexecutable code is included in the attached disk at the end of the thesis.30CHAPTER 3 SIMULATED DATA SETS AND ALGORITHM TESTINGThis chapter demonstrates how to apply the program MASSBAL.0 to the study ofigneous differentiation and compares solutions derived by x2 and R2 optimization usingsynthetic data sets.3.1 Calculation of Simulated Data Set.The synthetic data set comprises four synthetic basalts P, D1, D2 and D3, and threestandard minerals An70 , Fo80 and Di80 (Table 3. la). The composition of basalt P wassynthesized by modifying a natural basalt. The compositions of An 70 , Fo80 and Di80 werecalculated from their chemical formulas Ca0.7Na0.3A1 1.7Si2.308 , Mg 1.6Fe0.4SiO4 andCaMg0.8Fe0.2Si206 respectively. The derivative basalt compositions are created by simulatingcrystallization between a parent (P) and the daughter (D) compositions. The daughter magmaD 1 represents crystallization of 5g An70 out of 100g P. Crystallization of 5g of each An 70 andFo80 out of 100g P produces D2. Crystallization of 5g of each An70 , Fo80 and Di80 producesD3.The normalized compositions of magmas D1, D2 and D3 were calculated from thefollowing formulas:WtDi = ( 100*Wtp - 5*WtAn70 )/95.0 (14)WtD2 = ( 100*Wtp - 5*Wt An7c, - 5*WtFoso )/90.0 (15)WtD3 = ( 100*Wtp - 5*WtAn70 - 5*WtFoso - 5*WtDiso )/85 . 0 (16)where Wt represents weight fraction of oxide component.31Program MASSBAL.0 was used to calculate the relationships between the parent magmaP and each daughter magma and thereby to explore the effectiveness of the x2 optimizationtechnique. Three models were formed to explain the chemical difference between themagmas. These models were the same as the calculation procedure used to create the datasets: one phase fractionation from P to D 1 , two phase fractionation from P to D2 and threephase fractionation from P to D3.The comparison between the calculated and the expected solutions are presented in Table3.1b. The former are close to the true solutions, which are 5.0g for any of the participatingmineral phases, within computer roundoff errors. The X2 and R2 optimal solutions areidentical, and the minimum values x2„,in and R2min are both zero because there is no error.All the results are consistent with the expectations for error-free data.3.2 Simulated Data Sets with Associated Uncertainties.This section investigates the X2 solutions of synthetic data sets that have associateduncertainties. These synthetic examples are created by adding random errors to the syntheticdata set of Table 3.1a.Equation (18) is used to convert the error-free data into synthetic data with error.Wt' = Wt + fa^ (17)where- Wt' is the unnormalized new weight fraction of oxide component,- Wt is the error-free weight fraction of oxide component,- a is the simulated standard deviation of weight fraction of oxide32Table 3. la Synthetic chemical data set without assigned errors.Oxide Parent Minerals Daughters*An7O Fo80 Di80 Di D2 D3SiO2 50.6374* 50.5439 39.1914 53.1679 50.6423 51.2785 51.1673TiO2 2.6470 0.0000 0.0000 0.0000 2.7863 2.9411 3.1141Al203 13.6785 31.6982 0.0000 0.0000 12.7301 13.4373 14.2278FeO 11.2741 0.0000 18.7454 0.0000 11.8675 11.4854 12.1610MgO 7.3983 0.0000 42.0632 7.1330 7.7877 5.8835 5.8100CaO 11.2453 14.3576 0.0000 39.6992 11.0815 11.6971 10.0500Na2O 2.3189 3.4003 0.0000 0.0000 2.2620 2.3876 2.5281K2O 0.5206 0.0000 0.0000 0.0000 0.5480 0.5784 0.6125P2O5 0.2799 0.0000 0.0000 0.0000 0.2946 0.3110 0.3293* See text for detailed calculation procedures for daughter compositions.Table 3. lb Results of calculation with MASSBAL.0 on synthetic data sets without simulated errors .Processes* Expected resultsx2 techniqueCalculated resultsR2 techniqueW. of minerals xlmin W. of minerals R2minP to D 1 M1 = 5.00000g M 1 = 5.00004g 0.00000 M1 = 5.00004g 0.00000P to D2 Mi = 5.00000g M 1 = 5.00014g 0.00000 M1 = 5.00014g 0.00000M2 = 5.00000g M2 = 4.99998g M2 = 4.99998gM 1 = 5.00000g M I = 4.99981g MI = 4.99981gP to D3 M2 = 5.00000g M2 = 4.99994g 0.00000 M2 = 4.99994g 0.00000M3 = 5.00000g M3 = 4.99997g M3 = 4.99997g33component, and- f is a random number valued between -00 and +00.The standard deviation (a) is calculated from the relative error by (e):G = eWt^ (18)To simplify the modification, f was assigned one of three values: -1, +1 and O. Auniform random number generator was used to generate random numbers ranging from 0 to 1(Press et al., 1987). Random numbers from 0 to 1/3 caused addition of one standard deviationerror to the original weight fraction, 1/3 to 2/3 subtraction of one standard deviation, and 2/3to 0 no change.The procedures to generate the synthetic data are illustrated as a flowchart in Figure 3.1.Table 3.2 is an synthetic example with assigned 1% relative errors.3.2.1 Solutions for Synthetic Data Sets with 1% Relative Error.The calculated results of the synthetic example from Table 3.2 are presented in Table 3.3.It includes two parts, one part includes the optimal solutions for three fractionation processes,the other part describes the confidence regions around the optimal solutions. The bestsolutions obtained by the X2 and R2 techniques are different and they are no longer identical tothe true solutions. The X2 technique provides a better optimal solution closer than does the R 2technique. The confidence regions for the two phase fractionation process are illustrated inFigure 3.2. The true solution is included within the 95% confidence region.34Wt x 1 %perfect dataWt0- 1/3 1/3-2/3 2/3-1data with errorWt' = Wt + ff=0uniform randomnumber generator0 to 1Figure 3.1 Procedures to generate the synthetic data sets with assigned errors.35Table 3.2a Random numbers to generate error-assigned synthetic data sets.Oxides Parent DaughtersP D1 D 3SiO, -a -a 0 0TiO2 +a +a 0 +aA1,03 -a 0 +a 0FeO -a +a +a +aMgO -a -6 0 -aCaO 0 +a 0 0Na2 O -a -a 0 -aK2O 0 +a -a +aP2O5 -a +a 0 -aTable 3.2b Standard deviation calculated from 1% relative error on data in Table 3.1a.Oxide Parent DaughtersP D 1 D, D3SiO2 0.5064 0.5064 0.5128 0.5117TiO, 0.0265 0.0279 0.0294 0.0311A1,03 0.1368 0.1273 0.1344 0.1423FeO 0.1127 0.1187 0.1149 0.1216MgO 0.0740 0.0779 0.0588 0.0581CaO 0.1125 0.1108 0.1170 0.1005Na.,0 0.0232 0.0226 0.0239 0.0253K,0 0.0052 0.0055 0.0058 0.0061P,05 0.0028 0.0029 0.0031 0.0033Table 3.2c Assigned error synthetic data set No. 1.OxidesParent Minerals DaughtersP An70 Fo80 D i 80 D/ D, D3SiO, 50.1310 50.5439 39.1914 53.1679 50.1359 51.2785 51.1673TiO2 2.6735 0.0000 0.0000 0.0000 2.8142 2.9411 3.1452A1,03 13.5417 31.6982 0.0000 0.0000 12.7301 13.5717 14.2278FeO 11.1614 0.0000 18.7454 0.0000 11.9862 11.6003 12.2826MgO 7.3243 0.0000 42.0612 7.1330 7.7089 5.8835 5.7519CaO 11.2453 14.3576 0.0000 39.6992 11.1923 11.6971 10.0500Na20 2.2957 3.4003 0.0000 0.0000 2.2358 2.3876 2.5028K,0 0.5206 0.0000 0.0000 0.0000 0.5535 0.5726 0.6186P205 0.2771 0.0000 0.0000 0.0000 0.2973 0.3110 0.326036Table 3.3 Results of MASSBAL.0 calculations usint! synthetic data set No. I includin i optimal solutions andrelated confidence bounds.Processes Phases Solutions (g)Expected^x-^R2R2- min Confidence regionsSemi -axes90%^95%^99%Axis orientationsP to D i An70 5.0 5.459 4.960 9.026 0.081 1.30 1.40 1.59 (-1.00)P to D, An70 5.0 4.778 5.048 12.566 0.361 1.41 1.30 1.62 (-0.992,-0.121)Pogo 5.0 4.850 5.366 0.66 0.61 0.76 (0.121,-0.993)P to D3 An70 5.0 5.121 5.735 5.452 0.273 1.45 1.58 1.82 (0.866,0.156,-0.475)Pogo 5.0 4.966 5.565 0.80 0.87 1.01 (0.497,-0.160,0.853)Di80 5.0 5.270 5.326 0.56 0.61 0.71 (0.057,-0.975,-0.216)37Figure 3.2 Confidence regions around the optimal solution for the error-assigned syntheticdata set No.1. Bouble circle and cross hair represent the optimal solutions obtained by x2and R2 techniques, respectively. The open diamond is the true solution and is includedwithin the 95% confidence region.383.2.2 The Effect of a on z2 and R2 OptimizationA group of synthetic data sets with increasing assigned uncertainties were used tocompare solutions obtained by the X2 technique versus the R2 solutions. The method to createthe new synthetic data sets with larger assigned uncertainties is the same as the methoddescribed previously except the relative errors used are different.These simulated data sets are based on a three phase fractionation model which producesD3 from P. Table 3.4 presents the calculated amounts of An 70 , Fo80 and Di80 by both fittingtechniques for 11 data sets based on relative errors of 0.1 to 3%. The deviation from the truesolution is monitored by the relative deviation (Devi), defined by the following equation:Devi = (^- 5.0 )/5.0,^ (19)for j = An70 , Fog) and Di80 where- 5.0 represents the expected amount of mineral phases,- M. are the calculated amounts of mineral phases.Figure 3.3 illustrates the influence of standard deviation on the optimization of the x2 andthe R2 techniques. It shows that as analytical uncertainties increases, the optimal solutionsfrom both fitting techniques deviate from the true solutions. However, the X2 techniquecontinuously presents a much better and more stable solution than the R 2 technique. Within3% assigned relative error, the optimal solutions from the X2 technique deviate less than 8%from the expected solution, whereas the solutions obtained with the R 2 technique deviate by upto 42%.39Table 3.4 Results for both X 2 and R 2 fittinil of synthetic data sets for three phase fractionation.£* Phases Solution (g) Deviation (%) X2min R2minx2^R2 X2 R2An70 5.015 5.074 0.306 1.4820.10 Fo80 4.998 5.057 -0.048 1.136 5.3184 0.0027Di80 5.027 5.033 0.542 0.658An70 5.025 5.148 0.496 2.9540.20 Fo80 4.993 5.114 -0.134 2.270 5.3986 0.0109Di80 5.053 5.066 1.060 1.310An70 5.065 5.370 1.298 7.3920.50 Fo80 4.984 5.284 -0.304 5.670 5.4149 0.0682Di80 5.135 5.164 2.696 3.274An70 5.089 5.517 1.782 11.3320.70 Fo80 4.977 5.397 -0.454 7.932 5.4178 0.1339Di80 5.188 5.228 3.768 4.566An70 5.121 5.735 2.428 14.7081.00 Fo80 4.966 5.565 -0.688 11.306 5.4518 0.2732Di80 5.270 5.326 5.600 6.504An70 5.146 5.881 2.920 17.6201.20 Fo80 4.959 5.678 -0.828 13.550 5.4578 0.3932Di80 5.325 5.391 6.502 7.420An70 5.186 6.098 3.622 21.9621.50 Fo80 4.949 5.846 -1.014 16.912 5.4796 0.6145Di80 5.409 5.488 8.182 9.758An70 5.211 6.242 4.216 24.8441.70 Fo80 4.942 5.957 -1.160 19.142 5.4917 0.7891Di80 5.465 5.552 9.292 11.044An70 5.247 6.457 4.940 29.1382.00 Fo80 4.931 6.124 -1.366 22.476 5.5239 1.0922Di80 5.549 5.648 10.988 12.968An70 5.310 6.812 6.150 36.2442.50 Fo80 4.913 6.400 -1.736 28.008 5.5478 1.7064Di80 5.693 5.808 13.852 16.160An70 5.370 7.164 7.002 43.2803.00 Fo80 4.895 6.676 -2.140 33.314 5.5892 2.4567Di80 5.838 5.966 16.756 19.328* E is the relative error in percent.4050s.400300a>7F-E000.0^0.5^1.0^1.5^2.0^2.5^3.0^3.5Model Relative Error (%)Figure 3.3 Influence of standard deviation a on the optimization of the x2 and R2 techniques.41The x2 technique is a more robust fitting method than the R2 technique in light of theeffects of random error, however, good estimates of analytical error are essential to proper useof the x2 techniques. Having good analyses implies both precise analyses and good estimationof measurement error. Nevertheless, over- and underestimating the measurement errorsinevitably happens in practice. The next exercise investigate the consequences of poorestimation of analytical error on the x2 optimization. Again the effects are demonstratedthrough a set of synthetic data which are characterized by inaccurate estimates of error.The data in Table 3.3 were created with 1% assigned relative errors. When calculatingthe amounts of An70 , Fo80 and Di8 ,3 needed to produce the magmas D 1 , D2 and D3 from P,instead of using the true relative error of 1%, the calculations assumed other relative errors(e.g. 0.5%). Not all of the relative errors for the oxides were adjusted because the x2 and theR2 techniques are identical when the relative variances of all the oxide components are thesame. The results are in Table 3.5 for data sets with assigned relative errors for Si02 , Al203 ,FeO and CaO of between 0.1 and 5%.Figure 3.4 demonstrates the effects on the optimal solution and x2min values with varyingrelative error. It shows that underestimation of errors causes the solutions to deviate from thetrue solution much more than does overestimation. Different phases are influenced by thepoor estimation of errors differently. However, for An 70 , the solution deviates up to 12%with severe underestimation of error and only 2% for five times overestimation of error.Underestimating errors is the same as ignoring errors in the mass balance calculation, and thusis equivalent to applying the R2 technique to solve the mass balance problem. On the otherhand, overestimation of errors reduces the contribution of certain oxide components to theresiduals and thus affects the x2 optimization. Therefore, the effects of overestimating errorare less than the corresponding effects of underestimating error. The mean for a chi-squaredistribution function is approximately equal to the degrees of freedom of the problem (Bury,42E*(%) Phases Solution (g) Dev. (%) x2minAn70 6.464 29.2720.10 Fo80 6.772 35.442 205.58Di80 5.547 10.938An70 5.630 12.5980.20 Fo80 5.650 12.996 98.40Di80 5.383 7.662An70 5.359 7.1860.30 Fog° 5.286 5.720 49.19Di80 5.329 6.588An70 5.249 4.9760.40 Fog° 5.138 2.758 29.14Di80 5.307 6.114An70 5.194 3.8880.50 Pogo 5.065 1.304 19.26D i 8 0 5.295 5.894An70 5.164 3.2820.60 Fo80 5.025 0.494 13.72Digo 5.287 5.742An70 5.146 2.9160.70 Fo80 5.000 0.000 10.33Di go 5.287 5.742An70 5.134 2.6820.80 Pogo 4.984 -0.318 8.10Di go 5.277 5.546An70 5.126 2.5300.90 Fo80Din4.9735.274-0.5345.4706.56An70 5.121 2.4281.00 Fo80 4.966 -0.688 5.55Digo 5.270 5.400An70 5.118 2.3621.10 F080 4.960 -0.784 5.45Di80 5.267 5.336E( % ) Phases Solution (g) Dev. (%) x2m inAn70 5.116 2.3201.20 Fo80 4.956 -0.882 4.00Di so 5.264 5.274An70 5.115 2.2941.30 Fo- -80 4.954 -0.914 3.51Di x° 5.261 5.?17An 70 5.114 7.2801.50 Foso 4.949 -1.018 2.81Di so 5.255 5.096An 70 5.115 2.2941.70 F080 4.947 -1.096 2.34D i80 5.241 4.824An70 5.115 2.3081.80 F080 4.946 -1.080 2.16Di so 5.246 4.928An 7o 5.117 2.3462.00 Foso 4.945 -1.082 1.88Di so 5.241 4.824An 70 5.119 2.3882.20 Fo50 4.945 -1.108 1.67Di so 5.236 4.726An 70 5.120 2.4102.30 Fo50 4.944 -1.110 1.58Di so 5.234 4.678A n70 5.123 2/4562.50 Fo io 4.944 -1.132 1.44Di go 5.230 4.590An70 5.128 2.5663.00 F080 4.945 -1.104 1.21Di80 5.220 4.396An 70 5.143 ").86?5.00 Foso 4.947 -1.064 0.86Diso 5.197 3.936Table 3.5 x2 results of synthetic data sets with inaccurate errors replacirwg the true 1% relative error.* E are the relative errors which replace the true value I %.4330025225,•1751-420^P-r;co^125 I- 115075>.< 1072 5fr 5^•0'-250 1 2^3 4 5^6^1 2^3^4 5Model Relative Error (7.)^Model Relative Error (%)Figure 3.4 The effects of poor estimates of standard deviation a on X2 optimization.441975). The calculated minimum x2 value can be compared against the expected mean as a testof the quality of error estimation. In the problem represented in Figure 3.4 the degrees offreedom are 6. Note that the overestimation of a gives minimum x2 value of between 0.86 to4, whereas underestimation of a yields values of between 6 to 205, as the estimates of a getpoorer, the x2min rapidly deviates from the expected value of F_-- 6.3.2.3 x2 Confidence Regions: A Measure of Goodness -of-Fit and A Test of HypothesisThe error estimation is not the only factor to influence the X2 mass balance calculation,errors in the model or hypothesis also play a major role. In some instances, the number of theparticipating mineral phases can be uncertain. An extra and unwarranted mineral phase maybe added to the postulated fractionating mineral assemblage or, through oversight, an essentialphase may not be included in the hypothesis. The calculated confidence region around theoptimal solution can be used as a measure of goodness-of-fit and can test the model or thepetrologic hypothesis. The effect of error in the model is demonstrated by another set ofsynthetic test runs.Two types of synthetic examples based on the models in Table 3. la were made: i) byadding an unrequired phase to and ii) by dropping an essential phase from the mineralcrystallizing assemblage. The results are given in Table 3.6. It shows that models thatdropped one or more essential phases are associated with large x2 values (166 to 1422)compared to the x2 values associated with the true solution (5.45). Additionally, the optimalsolutions deviate significantly from the true solution. On the other hand, the models withadded, and unwarranted, extra phases have much smaller minimum x2 values (around 7) andthe optimal solutions are very small for the extra added phase.45rnTable 3.6 Results of the x2 technique modeling for the examples with incorrect models.True Solution Incorrect Hypothesis Confidence regionsModel x2 sol. xzmin Model^x2 sol. xzmin Orientations Semi-axes90%^95% 99%Three phase An70 = 5.121 Two phases An70 =8.378 (-0.995^-0.100) 1.25 1.35 1.55P to D3 Fo80 =4.966 5.45 P to D3 F080=5.935 313.73 (-0.100^-0.995) 0.61 0.66 0.75Digo = 5.270Three phase An70 = 5. 121 Two phases An70 =2.256 (-0.894^0.447) 1.45 1.57 1.79P to D3 Fogo =4.966 5.45 P to D3 Di80=7.634 713.60 (-0.447 -0.894) 0.84 0.91 1.05Disci = 5 - 27°Three phase An70 = 5.121 Two phases Fo80 =4.334 (0.973^0.230) 0.60 0.65 0.75P to D3 Fogo =4.966 5.45 P to D3 D180=7.020 166.48 (-0.230^0.973) 0.94 1.02 1.16Diso= 5 - 27°Three phase An70 = 5.121 One phase An70 =6.537 1422.48 (1.000) 1.29 1.39 1.59P to D3 Fo80 = 4.966 5.45 P to D3Digo =5.270One phase An70 =5.459 9.03 Two phases An70 =5.577 7.77 (0.981^0.195) 1.30 1.41 1.61P to D I P to D I Fo80=0.219 (0.195^-0.981) 0.64 0.69 0.79One phase An70 =5.459 9.03 Three phases An70 =5.509 (0.874^0.192^-0.447) 1.43 1.56 1.80P to D I P to D i Fo80=0.201 7.59 (0.450^0.032^0.893) 0.85 0.93 1.07Digo =0.130 (-0.185^0.981^0.059) 0.60 0.65 0.75The zero line is used to test the unwarranted mineral phase. Basalt D 1 is produced for Pby fractionation of An70 . If we add the unwarranted phase Fo80 to the mineral assemblage,the confidence regions in Figure 3.5 shows that Fo 80 is not needed, because the Fo80 zero lineis included within the 95% confidence region. If all solutions within the 95% confidenceregion are considered acceptable, then Fo 80 is not needed to produce D 1 from P.Compactness of confidence regions can be used to evaluate the models which havedropped an essential phase. For example. basalt D3 is produced from P by fractionation ofthree phases An70 , Fo80 and Di80 . Three two phases models including An70 and Fo 80 , An70and Di80 , and Fo80 and Di80 were used to reproduce D3 from P. Figure 3.6 is the comparisonof each these two phase models' confidence region with the corresponding two phaseprojection from the three phase model's confidence region. It shows that the projections aremore compact than the two phase confidence region. Consequently, the three phase model isbetter than two phase models.47Figure 3.5 Confidence regions around the X2 optimal solution (double circle) for a (false)model with an extra phase (Fo 80). Cross hair represents the optimal solutions obtainedwith R2 technique. The open diamond is the true solution. Note that the 95% confidenceregion includes the Fo 80 zero line.48Figure 3.6 The comparisons between the (false) two phase confidence regions (dashedellipses) and 2-D projections of the (true) three phase (An 70 + Fo80 + Di80) confidenceregions (solid ellipses) with 90%, 95% and 99% confidence limits, respectively forsynthetic fractionation. The double circles are the projections of the optimal solution forthe (true) three phase model. The stars are the optimal solutions for the (false) two phasemodel.49CHAPTER 4 NATURAL DATA SETSThree natural data sets were studied using the X2 technique. These data sets include oneassimilation and two fractionation problems. The Paricutin volcano assimilation problem andDiamond Crater fractionation are somewhat controversial topics and are taken from literature(McBirney et al., 1987 and Russell & Nicholls, 1987).4.1 Paricutin Volcano: An Example of AssimilationParicutin Volcano erupted in Michoacan, Mexico and its chemical differentiation hasbeen studied by different petrologists using different methods. Wilcox (1954), using agraphical method, originally interpreted its chemical variations as resulting from acombination of An70 and Fo80 fractional crystallization, and concurrent assimilation of up to25 weight percent continental crust. Miesch (1979) used a method of vector analysis toexamine variations within the same series of analyses and concluded that the compositionsevolved through three stages. The first could be explained by fractionation of olivine andplagioclase, the second by fractionation of Ca-poor pyroxene and plagioclase, and the third byassimilation of granite wall rocks. McBirney et al. (1987) used unweighted residual squarefitting technique to assess the effects of xenolith assimilation with fractionation of specificcombinations of minerals for each of the three major stages. Their calculations suggestedolivine fractionation for first stage; loss of orthopyroxene and plagioclase with addition ofxenolith for the middle stage; and the late stage is explained by a combination oforthopyroxene crystallization and assimilation of xenolith.The X2 technique developed in this thesis is used to study the differentiation processes ofParicutin volcano's middle and late stages chemical variation. It is applied to test the50contribution of xenolith to the evolution of Paricutin basaltic andesites. The chemical data andthe combined models of fractionation and assimilation are the same as McBirney's (1987) data.The difference between W-47-9 and W-48-5 represents the middle stage differentiationinvolving plagioclase, orthopyroxene and xenolith. The difference between W-48-5 and W-16-5 represents the late stage differentiation involving orthopyroxene and xenolith. Table 4.1is the chemical data of the assimilation problems in Paricutin Volcano. This paper did notreport measurement errors, therefore a 1% relative error and the average XRF error wereused. The comparison between McBirney et al.'s (1987) R 2 calculations and MASSBAL.0 x2calculations is given in Table 4.2a.Table 4.2a shows that the X2 optimal solutions are very different from McBirney's (1987)R2 optimal solutions. The assigned analytical errors used also influences the solutions. In thegroup using an average XRF error, the needed amount of xenolith to contaminate the originalmagma from the x2 fitting is more than double that from the R 2 fitting. Additionallyplagioclase is predicted to be accumulated rather than fractionated. On the other hand, for 1%relative error, the xenolith needed is a little smaller but still larger than the R 2 solution. Theplagioclase is fractionated but the amount is about four times smaller than that predicted byMcBirney et al. (1987). The amount of xenolith calculated by MASSBAL.0 is more similarto Wilcox (1954)'s conclusion that 25wt. % xenolith have been added, however, the amountsof olivine and plagioclase are quite different from Wilcox (1954)'s answers (3wt. % 01 andlOwt.% P1).Generally, the )(2 fitting suggests that we need a larger amount of xenolith to play a rolein the calc-alkaline magma's evolution, and the amounts of plagioclase and orthopyroxene are,relatively, not as important as xenolith.51Table 4.1^Chemical data of field assimilation from literature.Oxides Parent Daughter Plag Opx Xeno^Calc. daughterW-47-9 W-48-5SiO2 57.76 58.65 53.35 53.86 72.01 58.67TiO2 0.89 0.77 0.04 0.20 0.94Al2O3 17.18 17.42 29.51 2.46 14.44 17.29FeO 6.46 6.06 0.93 13.64 1.78 5.92MgO 5.61 4.00 0.10 28.53 0.89 4.03CaO 6.90 6.41 11.80 1.48 2.61 6.65Na2O 3.69 3.89 4.42 3.96 3.93K2O 1.22 1.49 0.23 3.00 1.52P205 0.29 0.30 0.13 0.31R2 = 0.122Oxides Parent Daughter Opx Xeno Calc. daughterW-48-5 W-16-5Si02 58.65 59.69 53.65 72.73 59.69TiO2 0.77 0.80 0.20 0.73Al203 17.42 17.17 1.75 14.58 17.17FeO 6.06 5.59 14.16 1.80 5.66MgO 4.00 3.71 27.87 0.90 3.66CaO 6.41 6.12 1.56 2.63 6.09Na20 3.89 3.97 4.00 3.90K20 1.49 1.66 3.03 1.62P205 0.30 0.28 0.13 0.28R2 = 0.02252Table 4.2a Comparison of solutions from MASSBAL.0 and the literature.W-47-9 to W-48-5Ave. XRF ErrorsMASSBAL.0 SolutionsNo ErrorsLiterature1% Relative ErrorsPlag^-2.78g 1.37g 6.37 6.7gOpx 3.17g 4.55g 7.00 6.6gXeno^-18.65g -13.02g -2.76 -8.1g2X min^53.808 196.324R2min^0.229 0.229 0.229 0.122W-48-5 to W-16-5MASSBAL.0 Solutions LiteratureAve. XRF Errors 1% Relative Errors No ErrorsOpx 0.86g 0.69g 0.21 0.55gXeno -6.23g -7.03g -8,.33 -10.67gx2min 91.108 118.598R2min 0.036 0.036 0.036 0.022Table 4.2b The confidence regions for the x2 solutions of the assimilation problem.W-47-9 to W-48-5PhasesSemi-axes1% relative errorSemi-axesaverage XRF errorAxis orientations Axis orientations90% 95% 99% 90% 95% 99%Plag 0.63 0.68 0.79 (0.102,-0.986,0.128) 0.89 0.96 1.11 (0.265,-0.964,-0.034)Opx 1.48 1.61 1.87 (-0.993,-0.094,-0.067) 1.46 1.59 1.87 (-0.779,-0.235,0.582)Xeno 1.88 2.04 2.36 (-0.054,-0.134,-0.989) 2.69 2.93 3.38 (0.568,0.127,0.813)W-48-5 to W-16-5Phases^1% relative error^ Average XRF errorSemi-axes90%^95% 99%Axis orientations90%Semi-axes95%^99%AxisorientationsOpxXeno0.522.200.562.380.642.72(0.994,-0.105)(0.105,0.994)0.782.240.842.430.972.78(1.000,-0.016)(-0.016,1.000)53The possible ranges of the amounts of mineral and xenolith are also calculated andpresented in Table 4.2b. The confidence regions calculated from 1% relative error are used inthe related figures. Figure 4.1 is the 2-D confidence region for the late stage differentiation(Xeno + Opx) in Paricutin volcano. Figure 4.2 is the three projections of the 3-D confidenceregions for the middle stage differentiation (Pl+Opx+Xeno) in Paricutin volcano. One of thephenomenon from Figure 4.1 and 4.2 is that the 95% confidence region around the X2 optimalsolutions for both the three and two phases assimilation problems do not encompass thesolutions given by McBirney et al., 1987 .4.2 Diamond Craters Volcano: An Example of FractionationThe Diamond Crater volcanic complex (DCV) is located within the high-alumina basaltpetrographic province in Southeastern Oregon and comprises Quaternary alkali olivine basalts(Russell & Nicholls, 1987). The rocks from the Diamond Crater are porphyritic withplagioclase and olivine phenocrysts; plagioclase dominates the phenocryst assemblage (16% -24%). Titanaugite is an important phase but is reported only within the groundmass (Russell& Nicholls, 1987). The chemical compositions of the DCV rocks are very similar, with SiO2around 47% weight fraction (Russell & Nicholls, 1987). The slight diversity in chemistrybetween the DCV rocks is thought to derive from the fractionation of both plagioclase andolivine by Russell and Nicholls (1987). The mass balance relationships are investigated hereto further test whether augite could play a role in the chemical differentiation of the DiamondCraters lavas.All the rock chemistry and mineral chemistry data come from Russell and Nicholls(1987), including the analytical errors for the whole rock wet chemistry. Based on the elementratios (Russell and Nicholls, 1987), the compositions of the two most divergent rocks DCVO4554Figure 4.1 Confidence regions around the x2 optimal solution (double circle) for the latestage differentiation of Paricutin volcano. The dashed line is the orthopyroxene zero line.McBirney et al. (1987)'s solution is not included within the solution spaces.55432oD1A.C o-22 3^4^5^5^7PI (g)e (1.37, 4.55)McBirney's sol. (6.7, 6.6)McBirney's sol. (6.7, -8.1)@ (4.55, -13.02)-10- 11-12o -13a)X -t4- 15-16-2^-1 0^1^2^3Opx (g)O (1.37, -13.02)McBirney's sol. (6.6, -ad)Figure 4.2 Three 2-D projections of the 3-D confidence regions (Pl+Xeno+Bi) with90%, 95% and 99% confidence for the middle stage differentiation of Paricutin volcano.Double circles are the projections of the X2 optimal solution. McBirney et al. (1987)'ssolution is not included within the solution spaces.56and DCVO41 were used to represent the compositions for the parent and derivative magmas.The chemical data are reproduced in Table 4.3a.Three models were formed to investigate the mass balance relation between DCVO45 andDCVO41. The first model is three phase fractionation which include Pl, 01 and Px. Thesecond model is two phase fractionation of P1 and 01. The other model is 01 one phasefractionation model. Table 4.3b gives the calculation results, including the optimal amounts ofmineral phases and confidence regions around the optimal solutions.Table 4.3b shows that the three phase (01+Pl+Px) model has the smallest x2 values.Figure 4.3 are the comparison between the two phase (01 + Pl) model confidence regions andthe two dimension projection (01 + Pl) of three phase model (01 + Pl+ Px) confidenceregions. It shows that the semi-axes of the projection from the 95% three phase confidencereigon are 0.29g and 0. 83g for 95% confidence region, respectively, while the semi-axes ofthe 95% confidence region for the two phase model are 0.28g and 1.27g, respectively, thus,the solution space of the two phase model is more compact. The pyroxene zero lines are notincluded in the projections on Pl+Px and 01+Px planes from the three phase (01+Pl+Px)confidence regions. Thus, the three phase fractionation model is better and augite is essentialto produce the chemical difference within Diamond Crater's lavas.4.3 Mount Meager: Another Example of Fractionation4.2.1 IntroductionThe Mount Meager (MM) volcanic complex, is located about 150 km north ofVancouver, British Columbia in the Coast Mountains (Figure 4.4), and seated between theLillooet River and Meager Creek. The Mount Meager volcanic field, together with five other57Table 4.3a Chemical data for alkali olivine basalts and minerals at Diamond Craters, Oregon.Oxides Parent Daughter 01 P1 Px Relative errorDCVO45 DCVO41SiO2 47.38 47.55 38.65 49.20 49.64 0.106TiO2 1.20 1.46 2.13 1.667Al203 17.30 16.86 30.99 2.94 0.173FeO 9.93 10.30 17.24 0.48 9.76 2.163MnO 0.17 0.17 0.27 0.22 11.765MgO 8.81 8.19 42.77 0.07 13.09 0.341CaO 11.02 10.67 0.26 15.23 20.78 0.272Na2O 2.65 2.94 2.94 0.45 1.132K2O 0.27 0.41 0.10 7.407P2O5 0.21 0.26 9.524Total 99.21 99.05 99.38 99.01 99.18Table 4.3b The confidence regions for the x2 solutions of three models for Diamond Craterfrationation. These models are three phase model (01 + PI + Py), two phase model (01 + Pl)and one phase model (01).Models Phases Optimal Sol. (g) x2min Semi-axes90%^95% 99%Axis orientations01 4.07 0.24 0.26 0.29 (0.914 -0.346 0.214)01+Pl+Px P1 10.98 66.04 0.62 0.67 0.77 (-0.313 -0.263 0.912)Py 2.93 1.37 1.49 1.70 (-0.259 -0.901 -0.349)01+Pl 01 3.69 253.09 0.26 0.28 0.31 (0.938 -0.345)P1 8.33 1.18 1.27 1.45 (0.345 0.938)01 01 0.78 1000.57 0.28 0.31 0.35 (1.000)58m (4.07, 10.98)• (3.69, 8.33)8• •7 131^2^3^4^5^601 (g)7(10.98, 2.93)65el 4oe3219^10^11^1208 13^4PI (g)7 1 I^..m (4.07, 2.93)654a 3211^2^3^4^5^601 (g)Figure 4.3 The comparison between the confidence regions for the two phase (01+Pl) model(dashed ellipses) and the 2-D projections of the three phase (01+Pl+Px) confidence regions(solid ellipses) for the Diamond Crater's fractionation. The confidence regions are with 90%,95% and 99% confidence limits, respectively. The open circles are the projections of the threephase optimal solutions. The star is the two phase optimal solution.591250Queen^ -1-BRITISH COLUMBIA♦ "OSound SILVERTHRONEN,%,^• ill -FRANKLIN GLACIERcetV.bALAL GLACIER&Cs^e.\ v.\^MEAGER CREEK• ELAHO VALLEYIILMOUNTotc, CAYLEYIPGARICIALDI LAKE1300 1 20 °52 °I *MOUNT \GARIBALDIEXPLORER WATTS POINTPLATEPACIFICPLATE48o13004, 'I'JUAN DE FUCAPLATE125 °•MOUNT BAKERWASHINGTONGLACIER PEAK12ID °Figure 4.4 The location of the Mount Meager volcanic complex with the Garibaldi belt of southwestern British Columbia (Greenet al., 1990).volcanic fields including Mount Garibaldi, Garibaldi Lake, Watts Points, Mount Cayley andSalal Glacier, forms the north-northwestern extension of Garibaldi volcanic belt. MountMeager occurs at the northern end of this belt (Green et al., 1988). The Garibaldi belt is thenorthern continuation of the Cascade volcanic belt and is believed to derive from thesubduction of the Juan de Fuca plate beneath the North American plate (Riddihough, 1978).The first study looking at the MM rocks involved the Bridge River pumice and ash(Nasmith et al., 1967). Later, the Meager Creek hot springs were investigated as a potentialsource of geothermal power, which resulted in detailed mapping (Anderson, 1975 and Read,1977) of the Mount Meager area.Read (1977) produced a 1:20,000 geologic map of the Meager Creek Geothermal Area.He subdivided the Mount Meager volcanic complex into seven major assemblages based onage, field relations, outcrop and geography (Read, 1977; 1991). In order of decreasing age,they are the Devastator assemblage, Pylon assemblage, Job assemblage, Capricornassemblage, Plinth assemblage, Mosaic assemblage and Bridge River assemblage. Twentyseven samples were collected by Stasiuk and Russell (1988) which cover these assemblages.Mass balance calculations are applied to investigate the chemical differentiation processeswithin one of these assemblages, namely, the Plinth assemblage.4.2.2 Mineralogy and Rock Chemistry Data of the Mount Meager Volcanic ComplexThe Mount Meager volcanic rocks include a wide variety of felsic through intermediate tobasaltic rocks. However, andesite dominates the volcanic edifice. The mineralogy andpetrography of each sample was studied in both hand sample and thin section. All samples arehighly porphyritic. The phenocrysts observed in the MM rocks are summarized in Table 4.4.Detailed mineralogy descriptions are given in Appendix A2.61Table 4.4 Phenocryst, microphenocryst and groundmass mineralogy in the Mount Meager volcanic rocks*.Assemblage Name^Phenocrysts^Xenocryst^Microphenocrysts^Groundmass crysts01 Opx Cpx Bi Am P1 Ap^Q^01 Opx Cpx Bi Am P1 Ap Q Z 01 Opx Cpx Bi Am P1 Ap Q OpPliocene MM-34^X X X X^x^x^X X X^X^X XMM-36 x^X X x^X X x^X XDevastator MM-27^X X^ X^ XMM-33 X^X X^X^X X^XPylon^MM-28^X^X X^X X^X X X^ X^XMM-35^X X X X^X X X^X^XMM-20^X^X X X X X X X X XMM-21^X X X^ X^X X^X^X XMM-25^X X^X X X X^X X X^XCapricorn MM-34^X^X X X^X^X^X X X^ X X XMM-38 X X X X X^X XPlinth^MM-26^X X^X^X^X X X^ X X X X XMM-29 X X X X X X X X X X XMM-30^X X^X X X^X XMM-40 X X X^X X^X^ X X^X XMM-39^X X^X X X^X X X^X XMosaic MM-37 X^X^X^X^X^X^X^X^X^XMM-17 X X X X X X X^XMM-18 X X X X^X X XMM-19 X^X^X^X X^X^XTable 4.4 continuedAssemblage Name^Phenocrysts^Xenocryst^Microphenocrysts^Groundmass crysts01 Opx Cpx Bi Am P1 Ap^Q^01 Opx Cpx Bi Am P1 Ap Q Z 01 Opx Cpx Bi Am P1 Ap Q Op^Bridge river MM-kl^X X X X^X X^MM-16 X^X^ X^X^XMM-2B^X X X X^X^X X^ X^X XMM-32 X X X X X^X X XMM-K3^X X X^X X^X X XMM-7 X X X X X^X XMM-2A^X X X X X^X X* 01 - Olivine, Opx - Orthopyroxene, Cpx - Clinopyroxene, Bi - Biotite, Am - Amphibole, PI - Plagioclase, Ap - Apatite, Q - Quartz, Z - Zircon, Op -Opaque.Major, minor and trace element concentrations of 27 whole rock samples from MMvolcanic complex were made with pressed crushed glass pellets and pressed powder, andmeasured by X-Ray fluorescence spectrometry. The analyses were made with a Philip PW1410 X-Ray wavelength-dispersion spectrometer in the Department of Geological Sciences.Major and minor elements include Si0 2 , Al203 , Ti02 , MgO, Fe203 , MnO, CaO, Na20, K20and P205. The analyzed trace elements are Nb, Y, Zr, Sr and Rb. Five duplicate sampleswere run to estimate the precision of sample preparation. XRF "machine precision" wasestimated by running the same pellet five times. The operating conditions for the analysis ofmajor, minor and trace elements are given in Appendix 3.Ferrous iron was measured volumetrically and ferric iron was calculated from the totaliron determined by XRF. Duplicate analyses of ferrous iron were performed to ensure areproducibility of ± 0.05%. Finally, the absorbed and structurally bound water weremeasured gravimetrically at 100°C and by the Penfield method, respectively. The operationalprocedures are given in Appendix 3.Chemical compositions of the volcanic rocks are listed in Table 4.5. The data includemajor, minor and trace element concentrations, H20+, H20- , loss on ignition (LOI) andCIPW normative mineralogy. Estimates of sample preparation and machine precision arereported in Table 4.6.Much attention has been paid during the XRF analysis in order to reduce the uncertaintyrelated to the machine instability. The counts on the sample used as a monitor were strictlyinvestigated. Only when the counting number of the monitor was within 1% of the meanvalue were the counts on the three unknowns recorded. Otherwise, this group of samplewould be reanalyzed. The results in Table 4.6 also show that uncertainties due to samplepreparation are about the same magnitude as the uncertainties due to machine instability,64Table 4.5 Whole rock chemistry of Mount Meager volcanic complex.Sample No. MM-2A MM-2B MM-7 MM-16 MM-20 MM-21 MM-24 MM-25 MM-26Assemblage* BRA BRA BRA BRA PYA PYA PCA PYA PLAUTM 684121 684121 682124 664124 628038 628038 604034 642038 649085SiO2 (wt. %) 64.35 67.69 68.33 68.64 68.68 67.32 62.40 62.76 67.10TiO2 0.57 0.46 0.47 0.44 0.42 0.52 0.62 0.71 0.54Al203 16.57 16.10 16.45 16.18 16.19 16.60 17.41 17.66 16.02Fe2O3 2.09 1.06 1.37 1.44 1.42 1.64 3.33 3.09 1.82FeO 2.32 1.86 1.68 1.43 1.36 1.62 1.06 1.54 1.40MnO 0.08 0.07 0.07 0.07 0.07 0.07 0.07 0.09 0.08MgO 1.92 1.36 1.46 1.29 1.25 1.71 3.21 2.67 1.74CaO 4.17 3.06 3.33 3.07 3.03 3.67 5.00 4.79 3.28Na2O 5.10 5.16 5.48 5.32 5.42 4.78 4.73 5.17 4.99K2O 2.23 2.41 2.44 2.46 2.50 2.43 1.65 2.07 2.83P2O5 0.19 0.15 0.15 0.15 0.14 0.18 0.26 0.24 0.23Total 99.59 99.39 101.22 100.49 100.48 100.55 99.76 100.81 100.03LOT 1.66 2.14 0.42 0.56 0.50 1.92 2.62 1.04 0.29H20+ 1.34 1.96 0.31 0.40 1.04 2.15 1.26 0.59 0.23H2O' 0.28 0.19 0.11 0.08 0.06 0.03 1.56 0.58 0.07Nb (ppm) 7 7 5 12 14 23 4 26 30Y 18 14 15 17 17 18 14 19 19Zr 180 153 158 157 152 180 197 224 202Sr 570 449 506 477 0.2 595 1120 769 712Rb 41 41 50 44 52 31 20 53 5Q (wt. %) 13.78 18.10 16.49 18.40 17.97 18.35 13.07 10.86 17.05Or 13.16 14.24 14.13 14.37 14.59 14.23 9.70 12.03 16.62P1 61.45 60.00 61.04 60.10 60.40 59.06 63.58 64.21 57.44Di 2.97 0.44 1.93 0.99 1.31 0.29 1.36 2.62 1.44Hy 5.28 5.15 4.07 3.73 3.39 5.28 8.14 5.94 4.33Mt 2.18 1.11 1.40 1.49 1.47 1.70 1.33 2.27 1.89Hm 0.00 0.00 0.00 0.00 0.00 0.00 1.43 0.61 0.000.79 0.64 0.64 0.61 0.58 0.72 0.86 0.97 0.75Ap 0.40 0.31 0.31 0.31 0.29 0.37 0.54 0.49 0.48BRA: Bridge River Assemblage; PYA: Pylon Assemblage; PCA: Pliocene Assemblage; CPA:Capricorn Assemblage; DGA: Devastator Assemblage; PLA: Plinth Assemblage; MOA: MosaicAssemblage.65Table 4.5 continuedSample No. MM-27 MM-28 MM-29 MM-30 MM-32 MM-kl MM-k3 MM-33 MM-34Assemblage* DGA PYA PLA PLA BRA BRA BRA DGA CPAUTM 655079 666082 644089 634109 621119 684121 684121 621080 610080Si02 (wt.%) 60.85 60.96 68.74 60.50 68.82 66.90 68.11 68.22 66.64TiO2 0.85 0.79 0.47 0.84 0.44 0.50 0.46 0.55 0.53Al203 17.76 17.79 15.72 17.74 16.09 16.63 16.23 16.29 16.21Fe203 2.73 2.02 1.65 3.40 1.14 0.95 0.56 1.30 3.72FeO 2.77 3.08 1.32 2.02 1.70 2.75 2.63 1.17 0.01MnO 0.09 0.10 0.07 0.09 0.07 0.08 0.07 0.04 0.08MgO 2.82 2.89 1.42 2.73 1.31 1.59 1.39 1.27 1.94CaO 5.08 5.07 2.70 5.06 3.07 3.49 3.17 1.90 3.53Na20 5.56 5.08 5.10 5.43 5.30 5.32 5.28 3.02 4.77K20 1.76 1.58 3.06 1.74 2.46 2.22 2.37 5.86 2.48P205 0.36 0.28 0.17 0.36 0.15 0.17 0.15 0.21 0.18Total 100.63 99.62 100.44 99.89 100.55 100.61 100.43 99.82 100.11LOI 0.26 0.09 0.19 0.94 0.45 1.73 0.74 3.49 0.02H2O+ 0.40 0.28 0.28 0.90 0.28 1.48 0.75 1.90 0.10H20- 0.16 0.14 0.06 0.10 0.08 0.24 0.10 0.67 0.06Nb (ppm) 14 4.1 12 11 9.0 6 4 10 7Y 20 15 19 18 16 17 15 26 19Zr 231 177 166 163 167 170 157 155 163Sr 773 748 533 500 485 526 508 4 535Rb 29 24 53 58 41 38 44 116 48Q (wt.%) 7.24 9.36 18.39 8.44 18.36 15.3 16.89 21.21 18.14Or 10.22 9.28 17.91 10.19 14.36 12.95 13.85 34.85 14.61P1 67.06 66.34 56.24 67.07 59.74 61.91 60.35 35.40 58.17Di 3.57 1.79 1.05 3.06 1.11 1.00 0.92 0.00 0.00Hy 7.19 9.45 3.71 5.94 4.34 6.82 6.46 3.73 5.34Mt 2.81 2.10 1.71 3.12 1.18 0.98 0.58 1.37 0.00Hm 0.00 0.00 0.00 0.27 0.00 0.00 0.00 0.00 2.59II 1.16 1.09 0.65 1.16 0.61 0.69 0.63 0.77 0.14Ap 0.74 0.58 0.35 0.75 0.31 0.35 0.31 0.44 0.3866Table 4.5 continuedSample No. MM-35 MM-36 MM-38 MM-39 MM-40 MM-18 MM-17 MM-19 MM-37Assemblage* PYA PCA CPA PLA PLA MOA MOA MOA MOAUTM 608106 592101 572095 638097 639101 604991 642161 603994 5821095i02 (wt.%) 63.15 65.21 59.06 68.10 66.94 50.68 51.33 49.60 50.91TiO2 0.66 0.63 0.75 0.48 0.52 1.44 1.34 1.49 1.51Al203 16.80 16.69 18.88 16.15 16.75 15.56 18.06 15.31 16.57Fe2O3 4.15 3.24 3.40 3.34 1.96 2.53 2.47 2.40 3.80FeO 0.52 0.88 2.60 0.04 1.60 7.70 6.38 8.10 5.02MnO 0.09 0.08 0.12 0.08 0.09 0.16 0.15 0.16 0.14MgO 2.37 2.12 3.27 1.68 1.89 8.00 6.06 8.52 6.24CaO 4.57 4.26 5.85 3.28 3.62 8.70 8.63 8.55 9.17Na2O 5.20 4.97 5.15 5.16 4.90 3.97 4.50 3.78 4.42K2O 1.87 2.41 1.02 2.55 2.44 0.78 0.74 0.77 1.32P2O5 0.21 0.27 0.28 0.14 0.17 0.27 0.35 0.28 0.51Total 99.59 100.76 100.37 101.01 100.86 99.80 100.01 98.95 99.61LOT 0.95 1.10 1.08 0.14 0.95 0.69 0.72 0.00 0.071120+ 0.83 1.10 1.01 0.08 0.98 0.92 1.58 0.53 0.34H2O- 0.11 0.12 0.44 0.02 0.14 0.40 0.15 0.32 0.08Nb (ppm) 7 5 3 8 6 17 10 19 12Y 17 16 14 15 18 16 17 17 18Zr 202 240 148 138 155 149 196 150 266Sr 692 882 779 438 483 581 943 537 1575Rb 36 42 15 46 47 17 13 14 6Q (wt. %) 12.90 14.88 7.81 17.92 17.18 0.00 0.00 0.00 0.00Ne 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 1.87Or 11.03 14.07 5.95 14.85 14.22 4.57 4.32 4.55 7.76P1 63.56 60.02 70.71 58.87 59.70 57.57 66.55 56.52 57.71Di 3.45 2.57 1.43 0.62 0.27 15.14 10.79 14.46 16.46Hy 4.81 4.50 8.99 4.27 5.54 2.02 0.04 1.08 0.00Mt 0.04 0.91 3.51 0.00 2.02 2.63 2.55 2.51 3.95Hm 2.86 1.63 0.00 2.30 0.00 0.00 0.00 0.00 0.00Il 0.92 0.87 1.03 0.18 0.71 1.99 1.85 2.08 2.09Ap 0.44 0.56 0.58 0.29 0.35 0.56 0.72 0.59 1.0601 0.00 0.00 0.00 0.00 0.00 15.52 13.18 18.22 9.1067Table 4.6a Estimates of XRF machine precision for sample MM-32.Repetition 1 2 3 4 5 mean Std* E*Si02 67.99 67.96 67.99 68.04 67.99 67.99 0.029 0.042TiO2 0.43 0.42 0.42 0.42 0.42 0.42 0.0045 1.06Al203 16.11 16.11 16.13 16.10 16.13 16.12 0.013 0.083FeO 2.66 2.67 2.66 2.67 2.67 2.67 0.0055 0.21Fe203 --- 1.0*MnO 0.07 0.07 0.07 0.07 0.07 0.07MgO 1.31 1.30 1.29 1.27 1.28 1.29 0.016 1.23CaO 3.05 3.04 3.05 3.05 3.04 3.05 0.0055 0.18Na20 5.30 5.21 5.23 5.14 5.23 5.22 0.057 1.10K20 2.41 2.40 2.40 2.40 2.40 2.40 0.0048 0.19P205 0.14 0.14 0.14 0.14 0.14 0.14Total 99.38 99.26 99.31 99.23 99.30Table 4.6b Estimates of sample preparation uncertainty for sample MM-32.Duplicate pellet 1 2 3 4 5 mean Std* E*Si02 69.07 69.25 67.99 69.10 68.82 68.85 0.50 o.73TiO2 0.44 0.44 0.43 0.44 0.44 0.44 0.0045 1.02Al203 16.26 16.19 16.11 16.37 16.09 16.20 0.11 0.71Fe203 2.99 2.84 2.66 2.78 2.87 2.83 0.12 4.24FeO 1.0*MnO 0.07 0.07 0.07 0.07 0.07 0.07MgO 1.32 1.27 1.31 1.32 1.31 1.31 0.021 1.60CaO 3.03 3.07 3.05 3.04 3.07 3.05 0.018 0.59Na20 5.14 5.15 5.30 5.24 5.30 5.23 0.078 1.49K20 2.49 2.45 2.41 2.44 2.46 2.45 0.029 1.19P205 0.14 0.14 0.14 0.15 0.15 1.44 0.0055 3.80Total 100.93 100.84 99.38 100.90 100.55* Std = Standard deviation; E = the relative error as a percentage.The relative error of FeO was estimated by the titration analysis.68except for SiO2 , Al203 and Fe,203 , for which the former is about one to 10 times larger thanthe latter. In mass balance calculations, only the errors brought in during the pelletpreparation are used to represent the total analytical errors.4.2.3 Mass Balance Calculations with the Plinth Assemblage (PLA)Among the seven assemblages being studied in the Mount Meager volcanic complex, thelavas which comprise the Plinth assemblage (PLA) show potentially comagmatic relationshipsbased on the field relations and mineralogy. Mass balance calculations are applied to furtherstudy the relations between the PLA rocks and investigate the differentiation processes.4.2.3.1 Comagmatic RelationshipsFive samples were collected from the Plinth assemblage which are numbered as MM-26,MM-29, MM-30, MM-39 and MM-40. All five rocks are highly porphyritic and with largeamounts of phenocrystal plagioclase and quartz. The plagioclase is unzoned. Quartz occurs asbroken crystals and may be xenocrystal. Medium-sized biotite exists in all samples and variesfrom a few grains to 10% in volume. Biotite is breaking down to amphibole, opaque andpyroxene. There are also small amounts of pyroxene microphenocrysts in MM-29, MM-39and MM-40, and a small amount of amphibole microphenocrysts in the groundmass of MM-26and MM-29 (Table 4.4).The rock chemistry of PLA rocks is presented in Table 4.7. All the samples have veryhigh SiO2 (-2:-. 68 wt. %), except MM-30 which contains about 60 wt. % SiO2. They areclassified as rhyodacites, except for MM-30 which is andesite. The element ratios P/K areclose together and within analytical error of each other, except MM-30 (Figure 4.5). Thissimilar chemical signature, together with the similar mineralogy suggests that except MM-30,69Table 4.7a The average electron microprobe analyses of PLA biotites and amphibolse.Sample MM-26 MM-29 MM-30 MM-39 MM-40 MM-26 MM-29Bi Bi Bi Bi Bi Am AmSi02 37.16 36.97 36.57 36.69 40.83 48.94 44.88Al203 13.43 13.36 13.21 13.58 12.22 5.93 9.35FeO 17.37 16.72 16.60 17.05 13.40 11.46 10.40MgO 13.91 14.19 13.75 13.96 14.03 16.54 16.03TiO2 0.37 4.30 4.21 4.31 3.04 0.97 1.98MnO 4.27 0.33 0.40 0.35 0.32 0.78 0.62Na20 8.84 0.55 0.65 0.60 1.73 1.43 2.14CaO 0.03 0.07 0.14 0.01 7.75 11.28 11.87K20 0.56 8.56 8.32 8.88 3.20 0.32 0.46F 0.17 0.28 0.28 0.51 0.16 0.21 0.17Cl 0.09 0.10 0.10 0.08 0.01 0.07 0.01Total 96.20 95.43 94.23 96.02 96.69 97.93 97.91N* 5 8 4 3 6 1 4Table 4.7b The average electron microprobe analyses of PLA feldspars and pyroxenes.Sample MM-26 MM-29 MM-30 MM-39 MM-40 MM-29 MM-39 MM-40PI P1 P1 P1 P1 Py Py PySi02 62.57 61.44 60.62 61.44 62.52 51.12 54.28 53.82Al203 23.46 23.92 24.69 23.92 23.50 2.59 1.71 2.04FeO 0.22 0.18 0.21 0.18 0.20 5.68 13.74 13.66MgO 0.00 0.01 0.01 0.01 0.01 16.11 28.71 28.66MnO 0.14 0.37 0.38Ti 02 0.53 0.19 0.16Na20 7.92 7.63 7.10 7.63 8.01 0.32 0.04 0.02CaO 5.10 5.56 6.52 5.56 4.89 22.17 1.31 1.24K20 0.71 0.67 0.60 0.67 0.77Total 99.98 99.41 99.75 99.41 99.90 98.66 100.35 99.98N 4 12 16 25 12 2 3 3*N is the number of measured points.700.16 T 7-MM- 3 00.12"*--- 0 08MM-26MM-400.04MM- 29 MM-390.00 ^ 1^1^I16^18^20^22^24^26 28^30Si/KFigure 4.5 The element ratios of PLA rocks in Mount Meager volcanic complex. Foursamples have the P/K ratios close enough within the analytical errors and are comagmatic.71the PLA samples are comagmatic and related by fractionation. Although the presence ofquartz xenocrysts might suggest that assimilation has played a role, the lack of other evidenceprevents further study of assimilation. The mass balance calculations using x2 fittingtechniques is applied to study the chemical relationships between these rocks.The compositions of the rocks are presented in Table 4.5. The compositions of thephenocrysts appearing in these five rocks were analyzed by electron microprobe. Theanalyzed minerals include plagioclase and biotite for all samples, pyroxene for MM-29, MM-39 and MM-40, and amphibole for MM-26 and MM-29. The analyses, calculated structuralformula and mole% end members are presented in Appendix A4. Table 4.7 contains theaverage compositions of the phenocrysts in these five samples.4.2.3.2 Mass Balance CalculationsIn terms of SiO2 wt%, the chemical difference between MM-40 and MM-26 isconsidered as the early stage in the chemical differentiation of PLA group and MM-29 andMM-39 as the late stage. The fractionation of both plagioclase and biotite was considered toplay a major role in the evolution of these rocks. The other two phases, amphibole andpyroxene exist in the groundmass. We want to find out whether pyroxene and amphibole arecrucial to explain the observed chemical differences. Mass balance calculations are applied tocarry out this job. Two models are formed. The first one is to test whether amphibolefractionated at the early stages to produce MM-26 from MM-40. The second task is to test thecontributions of both the pyroxene and amphibole to the late stage chemical evolution in PLAgroup producing MM-29 from MM-39. Table 4.8 gives the chemical compositions of the endmember rocks and mineral phases for these two processes.72Table 4-8 Chemical data of the fractionation problems in MM modeled by mass balance calculations.a. MM-40 to MM-26Oxides ParentMM-40DaughterMM-26Mineral phasesplagioclase^biotite amphiboleSi02 (wt.%) 66.94 67.10 62.57 37.16 48.94TiO2 0.52 0.54 0.37 0.97M203 16.75 16.02 23.46 13.43 5.93Fe203 1.96 1.82FeO 1.60 1.40 0.22 17.37 11.46MnO 0.09 0.08 4.27 0.78MgO 1.89 1.74 0.00 13.91 16.54CaO 3.62 3.28 5.10 0.03 11.28Na20 4.90 4.99 7.92 8.84 1.43K20 2.44 2.83 0.71 0.56 0.32P205 0.17 0.23b. MM-39 to MM-29Oxides ParentMM-39DaughterMM-29 plagioclaseMineral phasesbiotite^amphibole pyroxeneSi02 (wt.%) 68.10 68.74 61.44 36.97 44.88 51.12TiO2 0.48 0.47 4.30 1.98 0.53Al203 16.15 15.72 23.92 13.36 9.35 2.59FeO 3.05 2.80 0.18 16.72 10.40 5.68MnO 0.08 0.07 0.33 0.62 0.14MgO 1.68 1.42 0.01 14.19 16.03 16.11CaO 3.28 2.70 5.56 0.07 11.87 22.17Na20 5.16 5.10 7.63 0.55 2.14 0.32K20 2.55 3.06 0.67 8.56 0.46P205 0.14 0.1773Plagioclase, biotite and amphibole were considered to participate in the early stage offractionation. Table 4.9 is the calculated results. It shows that about 6% plagioclase, 2%amphibole and 1% biotite crystallized from MM-40 to produce MM-26. Figure 4.6 containsthe three 2-D projections for the 3-D confidence regions. It shows that the zero line of biotiteis included in the 95 % confidence region around the optimal solution on the projections onPl+Bi and Bi+Am planes. Thus, biotite did not contribute to the early stage chemicaldifferentiation in PLA group. It was replaced by amphibole, which is consistent with thepetrographic observation that biotite is breaking down to amphibole.Plagioclase continues to dominate the fractionation assemblage during the late stagechemical differentiation in PLA group. Pyroxene microphenocrysts start to appear andscattered throughout the groundmass. The fractionation models formed by mineral assemblageof plagioclase, biotite, amphibole and pyroxene are used to explain the chemical differencebetween MM-39 and MM-29. The calculated results are also presented in Table 4.9. It showsthat the four phase model (Pl+Bi+Am+Px) has to be rejected since it produces a largeminimum X2 value (6199). Figure 4.7 shows the projections on Am+Px and Pl+Px planesfrom the three phase (Pl+Am+Px) confidence regions. The zero lines of pyroxene isincluded in the projections of 95% confidence regions. Figure 4.7 also shows the comparisonof two phase (P1+Am) confidence regions and the projection on Pl+Am plane from the threephase (Pl+Am+Px) confidence region. The semi-axes of the former are about 0.50g and1.60g, respectively, while the semi-axes of the latter are 0.48 and 1.81g, respectively,therefore, the confidence regions for the two phase model are a little more compact. Thus,pyroxene is not considered as a major mineral phase in the late stage of chemicaldifferentiation in PLA group. Only two phases, which are about 10% plagioclase and 3%amphibole caused the difference between MM-39 and MM-29.74Table 4-9 The mass balance results of PLA fractionation.a. Results of three fractionation models for MM-40 to MM-26 (one three phase model, two two phasemodel).Models Phases Optimal Sol. (g) X2min90%Semi-axes95% 99%Axis rientationsP1 5.05 1.86 2.00 2.28 (-0.984 -0.156 0.080)Pl+Bi +Am Bi -1.10 325.54 0.97 1.05 1.19 (-0,175 0.834 -0.523)Am 2.47 0.49 0.53 0.61 (0.015 -0.529 -0.849)Pl+Bi PI 6.06 504.46 1.91 2.05 2.32 (-0.996 -0.090)Bi 0.77 0.76 0.82 0.92 (0.090 -0.996)Pl+Am P1 5.58 345.59 1.87 2.01 2.27 (1.000 -0.007)Am 2.00 0.58 0.63 0.71 (-0.007 -1.000)b. Results of two models for MM-39 to MM-29 (two two phase model, two three phase model and onefour phase model).Models Phases Optimal Sol. (g) X2 m in90%Semi-axes95% 99%Axis orientationsPI -18.59 0.03 0.04 0.04 (0.062 -0.833 -0.494 -0.242)P1+Bi+Am+Px Bi -9.97 6199.72 0.27 0.30 0.34 (-0.139 0.414 -0.301 -0.847)Am -1.74 2.44 2.66 3.07 (0.870 -0.069 0.377 -0.311)Px 4.98 3.60 3.92 4.53 (-0.460 -0.360 0.723 -0.356)P1 12.153 0.25 0.27 0.31 (0.058 0.543 0.838)Pl+Am+Px Am 4.078 51.70 2.92 3.16 3.62 (0.884 0.362 -0.297)Px -0.741 0.92 1.00 1.15 (0.465 -0.758 0.459)PI -19.82 0.04 0.04 0.05 (0.064 -0.859 -0.508)Pl+ Am +Bi Am -13.83 6315.52 2.97 3.21 3.68 (0.961 0.191 -0.201)Bi 7.06 0.59 0.63 0.73 (0.270 -0.475 0.838)P1 -0.62 2.04 2.20 2.51 (-0.997 -0.075)P1+ Bi Bi -8.24 7313.89 0.05 0.05 0.06 (0.075 -0.997)P1+ Am P1 10.595 58.47 1.68 1.81 2.06 (1.000 -0.024)Am 3.178 0.46 0.48 0.56 (-0.024 -1.000)75Figure 4.6 Three 2-D projections of the 3-D confidence regions (Pl+Am+Bi) with 90%,95 % and 99% confidence for the PLA early stage differentiation of Paricutin volcano. Thedashed ellipses are the confidence regions with 90%, 95% and 99 % confidence limits forthe two phase (Pl+Am) model. Double circles are the projections of the x2 optimalsolution.763be I(11.. —1- 3• (4.08, -0 .74)-5^......... , -1^1^3^5^7Am (g)253-3- 5• ( 1 2.15. -0.74)7• (10.80, 3.18). (12.15, 4.08)8^10^12^14^lgPI (g)Figure 4.7 The 2-D projections (solid ellipses) of 3-D (Pl+Am+Px) confidence regions with90%, 95% and 99% confidence limits, respectively for PLA late stage differentiation. Thedashed ellipses are the two phase (Pl+Am) confidence regions with 90%, 95% and 99%confidence limits, respectively. The horizontal dashed lines are the pyroxene zero line.77CHAPTER 5 SUMMARY AND DISCUSSIONA new x2 minimization fitting technique has been developed in this thesis to study massbalance relations during igneous differentiation processes. The X2 fitting differs from theconventional R2 minimization by weighting the residual squared with the variance of themeasurement error. When analytical data are error free or the variances for all the data pointsare the same, these two fitting techniques provide identical results. However, when thechemical data (individual oxides) vary in analytical quality, the X2 is a better fitting techniquethan the R2 technique.The advantages the X2 fitting has over the R2 fit can be summarized on three points.First, in the X2 technique, the fitted parameters depend more on the better quality data than onthe poorer quality data , while the R 2 techniques treat each data point equally. Secondly, theX2 technique can estimate a confidence region around the optimal solution. And third, it ispossible to estimate statistically the goodness-of-fit in the X2 fitting, whereas in the R2technique, a critical R2 value is arbitrarily set to quantitatively assess the acceptability of theoptimal solution can be accepted or not.It was demonstrated using synthetic data that the X2 fitting is a more reliable techniquethan the R2 technique. The optimal solution found by the X2 technique is much closer to thetrue solution than that by the R 2 technique. With the relative error for all the data pointsincreasing up to 3% in the synthetic data, the X2 technique continuously produces a muchbetter and more stable solution than the R 2 technique. Within the range of 3% relative error,the optimal solutions from the X2 technique deviate from the true solution less than 8 %,whereas the optimal solutions from the R2 technique deviate more than 42%. Therefore, asmeasurement errors get greater, X2 fitting provides a better solution compared to that chosenby R2 minimization.78The importance of the X2 technique is also reflected in the calculated confidence regionsaround the optimal solutions. The geological systems we study are highly complicated and theanalytical data used to describe the system are subject to measurement errors. There is arandom component included in the measurement error, thus the optimal solution is not theunique realization of the true solution. Instead, there are many other solutions which can alsorepresent acceptable estimates of the true solution. The X2 technique brings in the confidenceregion to represent the possible realizations to the true solution, rather than using a uniqueoptimal solution as is done in R2 fitting.The confidence region is the zone of the possible solutions which share the sameconfidence limit. It includes two factors. One is the shape of the region, and the other is thesize reflecting the confidence level. The first one depends on the dimension of the solvedproblem, which is the number of fitted parameters. The shape of the confidence region is aline segment for one dimension centered on the optimal solution, an ellipse for two dimensionsand an ellipsoid for three dimensions. The confidence levels used in this thesis are 90%, 95%and 99%, which represent the probabilities of the true solution falling within the calculatedregions (90, 95 and 99%, respectively).The critical confidence level we apply is 95%. The solutions within the 95% confidenceregion around the x2 optimal solution are considered to be good estimates of the true solution.For all the synthetic data sets, the true solution is included in the 95% confidence region of theoptimal solution. The goodness-of-fit is also investigated through the 95% confidence region.The model with the optimal solution whose 95 % confidence region is more compact isregarded as better model to explain the differentiation processes. For example, if two phasefractionation model has an optimal solution whose 95% confidence region is more compact79than that of the three phase fractionation model, then the two phase model is better to explainthe chemical difference of two magma end members.The x2 minimization fitting appears to be a good technique for quantitative interpretationof chemical data. It is powerful to test hypotheses explaining chemical difference of twomagma end members. However, the estimation of analytical uncertainties is very important toget a good solution. The under- or overestimation of errors have a large effect on the optimalsolution and the confidence regions. The underestimation of error has larger effect on thesolution than the overestimation of error. The overestimation of error would bring us asmaller and unreliable confidence region.Program MASSBAL.0 was written to do all the x2 mass balance calculations. It usesSingular Value Decomposition (SVD) to search the optimal solution and estimate theconfidence region around the optimal solution.The newly developed x2 minimization fitting technique has been applied to study threefield problems. One is the assimilation in Paricutin volcano of Mexico, the other two arefractionation in Mount Diamond Crater of Oregon and Mount Meager Plinth Assemblage(PLA) group of British Columbia. The magmatic differentiation models of the first two wereinvestigated to clarify some debate on them. The chemical difference of PLA group was usedto provide some insight the magmatic evolution of the Mount Meager volcano.In the field assimilation problem, the x2 technique was applied to test the differentiationmodels explaining the chemical variation in the volcanic rocks of Paricutin volcano inMichoacan, Mexico. The evolution of this basaltic andesite to andesite volcano was dividedinto three stages. The early stage involved the fractionating of 01 and Pl, the middle stageinvolved the fractionating of Opx and Pl., and the last stage involved the fractionating of Opx.80However, the better explanation for the chemical variation can be obtained if a certain amountof xenolith is added in the magmatic system concurrent fractionation of mineral assemblage.The x2 mass balance calculations are applied to test the influence of xenolith assimilation.The X2 optimal solution calculated from MASSBAL.0 were different from the McBirneyet al.'s (1987) results. The X2 calculations suggest that for the middle stage a larger amount ofxenolith (about 20 wt. %) was needed to contaminate the magma on its way to the surface withabout 3% plagioclase added instead of taken away from the magma and 3% orthopyroxenefractionation. In the late stage, about 6% of xenolith was added in the magma while about onepercent of orthopyroxene crystallized from the magma. The difference between the x2calculations and the R2 calculations is so large that the solution space of the 95% confidenceregion around the X2 optimal solution does not include the R 2 optimal solution obtained byMcBirney et al. (1987). The x2 optimal solution for the amount of xenolith participating themiddle stage is consistent to Wilcox's (1954) explanation, which was about 20% xenolithcontaminating the magma.From the Diamond Craters x 2 calculation results, about 3% weight fraction of augite isessential in addition to 11 % plagioclase and 4% olivine, to produce the slight chemicaldifference within the DCV basalts. This result is different from Russell and Nicholls (1987)who concluded that only P1 and 01 were needed. The existing of augite in the groundmassactually infers that the augite has already participated in the fractionation process.The chemical difference of PLA rocks were divided into early and late two stages basedon the petrography and element ratios. Three phase plagioclase, amphibole and biotitefractionation model was used to explain the early stage chemical evolution. The zero line ofbiotite is included within the 95% confidence region around the x2 optimal solution, thus, verypossibly biotite has stopped crystallizing to produce MM-26 from MM-40. In the late stage,81only plagioclase and amphibole participated the fractionation process to produce MM-29 fromMM-39, while pyroxene is not needed.Although we can not prove that the X2 solution is better than the R2 solution, the x2technique provides the confidence region around the best solution. We can estimate thegoodness-of-fit by investigating the compatibility of confidence region, which has not beentried before. The further improvement in the x2 fitting can be obtained by more precisemeasurement errors and the representative compositions of mineral phases and xenolith.This thesis is the first try to use X2 as a maximum likelihood estimate in modelinggeological data. It reveals the importance of quality in analytical measurements where thosechemical compositions are to be used in modeling. Also it treats differently the good and badanalyses, thus accentuating the need for precise analyses and the reporting of all analyticalerrors.82REFERENCESAlbarede, F., and Provost, A., Petrological and geochemical mass-balance equations; analgorithm for least-square fitting and general error analysis, Comput. Geosci., 3, 309-326,1977.Anderson, R.G., The geology of the volcanics in the Meager Creek map area, southeasternBritish Columbia, B. Sc. thesis, 130pp., Univ. of British Columbia, 1975.Bartley, J.M., Evaluation of REE mobility in low-grade metabasalts using mass-balancecalculations, Norsk Geologisk Tidsskrift, 66, 145-152, 1986.Bea, F., A method for modelling mass balance in partial melting and anatectic leucosomesegregation, J. Metam. Geol., 7, 619-628, 1989.Bowen, N.L., The evolution of the igneous rocks, Princeton University Press, 332pp., 1928.Bury, K.V., Statistical models in applied science, New York, Wiley, 625pp, 1975.Chakraborty, K.R., A new method of material balance evaluation in metamorphicdifferentiated systems, Contrib. Mineral. Petrol., 65, 101-110, 1977.DePaolo, D.J., Trace element and isotopic effects of combined wallrock assimilation andfractional crystallization, Earth and Planetary Science Letters, 53, 189-202, 1981.Dipple, G.M., Wintsch, R.P., Andrews, M.S., Identification of the scales of differentialelement mobility in a ductile fault zone, J. Metam. Geol. , 8, 645-661, 1990.83Francis, D., The pyroxene paradox in MORB glasses -- a signature of picritic parentalmagmas, Nature, 319, 586-589, 1986.Gerlach, T.M., Luth, W.C., Mass balance differentiation models for the 1959 Kilauea lid lavalake, Eos Trans. AGU, 61, 1143, 1980.Grant, J.A., The isocon diagram; a simple solution to Gresens' equation for metasomaticalteration, Econ. Geol., 81, 1976-1982, 1986.Green, N.L., R.L. Armstrong, J.E. Harakal, J.G. Souther, and P.B. Read, Eruptive historyand K-Ar geochronology of the late Cenozoic garibaldi volcanic belt, southwestern BritishColumbia, Geol. Soc. Am. Bull., 100, 563-579, 1988.Gresens, R.L., Composition - volume relationships of metasomatism, Chem. Geol., 2, 47-65,1966.Hanor, J.S., and K.C. Duchac, Isovolumetric silicification of early Archean komatiites;geochemical mass balances and constraints on origin, J. Geol., 98, 863-877, 1990.McBirney, A.R., H.P. Taylor, and R.L. Armstrong, Paricutin re-examined: a classic exampleof crustal assimilation in calc-alkaline magma, Contrib. Mineral. Petrol., 95, 4-20, 1987.Michael, P.J., and R.L. Chase, The influence of primary magma composition, H 2O andpressure on Mid-Ocean Ridge basalt differentiation, Contrib. Mineral. Petrol., 96, 245-263,1987.84Myers, J.D., and C.L. Angevine, Mass balance calculations with end member compositionalvariability: application to petrologic problems, Eos Trans. AGU, 67, 404, 1986.Myers, J.D., C.L. Angevine, and C.D. Frost, Mass balance calculations with end membercompositional variability: applications to petrologic problems, Earth Planet. Sci. Lett., 81,212-220, 1987.Myers, J.D., C.D. Frost, and C.L. Angevine, A test of a quartz eclogite source for parentalAleutian magmas: a mass balance approach, J. Geol., 94, 811-828, 1986.Nasmith, H., and G.E.R. Mathews, Bridge River ash and some other recent ash beds inBritish Columbia, Can. J. Earth Sci., 4, 163-170, 1967.Nesbitt, H.W., and J.J. Cramer, Graphical representation of mineral equilibria and materialsbalances in igneous rocks, Contrib. Mineral. Petrol., 78, 136-144, 1981.Nicholls, J., The statistics of Pearce element diagrams and the Chayes closure problem,Contrib. Mineral. Petrol., 99, 11-24, 1988.Olsen, S.N., Mass-balance and mass-transfer in migmatites from the Colorado Front Range,Contrib. Mineral. Petrol., 85, 30-44, 1984.Press, W.H., B.P. Flannery, S.A. Teukolsky, and W.T. Vetterling, Numercal recipes in C,Cambridge University Press, 818pp., 1987.Read, P.B., Geology of Meager Creek geothermal area, British Columbia; Geological SurveyCanada, Open File 603, 1977.85Read, P.B., Mount Meager complex, Garibaldi belt, southwestern British Columbia,Geosci. Can., 17, 167-170, 1991.Riddihough, R.P., The Juan de Fuca plate, Eos Trans. AGU, 59, 836-842, 1978.Russell, J.K., and J. Nicholls, Early crystallization history of alkali basalts, Diamond Craters,Oregon, Geochim. Cosmochim. acta, 51, 143-154, 1987.Stasiuk, M.V., and J.K. Russell, Petrography and chemistry of the Meager Mountainvolcanic complex, southwestern British Columbia, Current Research, Part E, Geol. SurveyCan., Paper 89-1E: 189-196, 1989.Stasiuk, M.V., and J.K. Russell, Quaternary volcanic rocks of the Iskut River region,northwestern British Columbia, Current Research, Part E, Geol. Survey Can., Paper90-1E: 227-233, 1990.Stasiuk, M.V., and J.K. Russell, The Bridge River assemblage in the Meager Mountainvolcanic complex, southwestern British Columbia, Current Research, Part E, Geol. SurveyCan., Paper 90-1E: 153-157, 1990.Stormer, J.C. , and J. Nicholls, XLFRAC: a program for the interactive testing of magmaticdifferentiation models, Comput. Geosci., 4, 143-159, 1978.Stout, M.Z., and J. Nicholls, Mineralogy and petrology of Quaternary lavas from the SnakeRiver Plain, Idaho, Can. J. Earth Sci., 14, 2140-2156, 1977.86Wright, T.L., and P.C. Doherty, A linear programming and least square computer methodfor solving petrologic mixing problems, Geol. Soc. America. Bull., v. 81, 1995-2008, 1970.87APPENDICESAl Linear Algebra Theorem of SVDSVD is based on the following theorem of linear algebra: Any MxN matrix (MxN) Acan be written as the product of an MxN column-orthogonal matrix U, an NxN diagonalmatrix W' with positive or zero element, and the transpose of an NxN orthogonal matrix V (Press et al., 1987). This theorem can simply be expressed by:A =-- U•W•VT.88A2 Detailed Petrographic DescriptionsThe petrography of each MM sample is studied and the detail information is described below.A2.1 Pliocene assemblage (PLO) PLO rock members are: MM-24 and MM36.The Pliocene Assemblage was formed during the late Tertiary volcanic episode. As theoldest unit in this volcanic area, contacting with the basement metamorphic unit directly, it isalmost covered by the younger units. Only four outcrop traces of flows, intrusions andbreccias are exposed, which is mostly appears at the southern and southwestern edge of themap area, except one plagioclase andesite flow is exposed in the western Affliction Creek. Atleast three flow are recognized within this assemblage. The stratigraphy can be expressed asfollowing. The bottom layer of medium grey-green aphanitic flows and breccias is overlaid bythe volcanic breccia with dominant pluton and some aphanitic volcanic clasts which is as thickas up to 300m. The top two layers are plagioclase andesite flows and breccias, plagioclase,quartz and hornblende dacite flows and breccias in decreasing age order (Read, 1977).Samples yield the age about 1.9 to 2.2 Ma with an error of 0.2Ma by K-Ar dating (Green,Armstrong et al., 1988). Two samples were collected from the top two layers of this unit.One is from the plagioclase andesite flows (2.0+-0.1Ma), another from the quartz, plagioclaseand hornblende dacite flow (1.9 +- 0.2Ma).Although only two samples were collected from this area, they are quite different from themineralogy and petrology. MM-36 came from the andesite unit. The detailed informationabout the thin sections can be inferred from the appendix A4.1. The most distinguish feature89is: MM-24 has the quartz and biotite phenocrysts and MM-36 doesn't. They are common inthat both have the plagioclase and hornblende phenocrysts, clinopyroxene microphenocryst.Sample number: MM-24Thin section:Plagioclase phenocryst, quartz xenocryst, clinopyroxene, hornblende and biotite phenocrysts,minor apatite and zircon microphenocrysts, quartz glomerocryst are carried in the groundmassof plagioclase, hornblende, opaques and glass. Anhedral to subhedral plagioclase is zonedwith the size from 2mm to 0.5mm. Most of the pl grains have 0.01mm wide remelted fringeand very thin overgrowth rim, except one or two small pl phenocrysts without the overgrowthrim. All the quartz grains are embayed with the fractured surface, and the melt corroded thequartz from the cracks. The disequilibrium relation between the quartz and the melt and theabsence of the quartz in the groundmass infer that the quartz grains are the xenocrysts. Thebiotites were broken down to plagioclases and Fe oxides. The hornblende both as the twinnedphenocrystals and in the groundmass were altered with the opaque rim. Some of thehornblende microphenocrystals were totally replaced by the opaques. The number of theclinopyroxene is a small amount, only two were found in one thin section, but the size is about0.15mm, not so small. There is not any reaction rim and alteration phenomenon. One quartzglomerocryst about 0.07mm in diameter with minor hornblende is carried in the groundmass,the hornblende was almost altered to the opaques. One zircon microphencryst was noticed onthe groundmass with the unusual interference color.Mineral classification:Mineral^Mode^NamePlagioclase 5% AndesiteQuartz^2%Hornblende 1%Biotite^1%90Clinopyroxene^1%Groundmass 90%Sample numbler: MM-36Thin section:Plagioclase, hornblende, clinopyroxene and opaques phenocrysts, plagioclase, clinopyroxeneand opaques glomerocryst are carried in the groundmass of plagioclase, hornblende,orthopyroxene, clinopyroxene, opaques and glass. Plagioclase varies from 0.15mm to0.03mm in size. Two types of plagioclase exist in this thin section. One type with theremelted fringe and the overgrowth rim, while the other doesn't have, but both types arestrongly zoned. The hornblende appears the light yellow to deep reddish brown pleochroism.All the euhedral, subheral and anhedral hornblende grains are fresh without any reaction rimand alteration surface. While the orthopyroxene has the hornblende coronas. The irregularhabit of the opaques and their coexisting with the ferromagnesian mineral infer that they camefrom the alteration. The clot comprises the orthopyroxene, opaque and plagioclase. Theorthopyroxene has the hornblende reaction rim, the opaques keep the pyroxene shape wirh thereddish brown hornblende rim. These facts suggests that the clots are produced byorthopyroxene breakdown.Mineral classification:Mineral^Mode^NamePlagioclase 10% DaciteOrthopyroxene^2%Hornblende 3%Opaques^1%Clot 0.5%Groundmass^87.5%91A2.2 Pylon Assemblage (PYA) PYA rock members are: MM-20, MM-21, MM-25, MM-28 and MM-35.Sample number: MM-20Thin section:Plagioclase, hornblende, orthopyroxene, biotite and opaques phenocrysts, plagioclase andhornblende glomerocryst are carried in the groundmass of plagioclase, hornblende and glass.Subhedral to anhedral plagioclase varies from 0.2mm to 0.05mm in size, slightly zoned withthe interference color high to first-order yellow. Glass inclusion is rich in the plgioclasegrains. Euhedral to subhedral hornblende has the light greenish yellow to deep brownpleochroism with the apatite inclusion. Minor anhedral biotite was slightly remelted.Microphenocrysts of orthopyroxene were found in the groundmass. Most of the opaques havethe irregular shape, a few grains occur as the nice square, which are possible magnetite. Thegroundmass is glassy and rich in plagioclase. The plagioclase and hornblende in theglomerocryst are the same as the phenocrysts. And there isn't any reaction or alterationrelation is noticed. Therefore this crystal clot is formed by the minerals growing together.Mineral classification:Mineral^Mode^NamePlagioclase 10% AndesiteHornblende^3%orthopyroxene 1%opaque^1%glomerocryst^0.5%biotite 0.5%Groundmass^84%92Sample number: MM-21Thin section:Plagioclase, hornblende, orthopyroxene and opaque phenocrysts, plagioclase, orthopyroxene,hornblende, biotite and opaque glomerocrysts are carried in the groundmass of plagioclase,orthopyroxene, opaque and glass. The subhedral to anhedral plagioclases are common zoned,with the glass liquid inclusion. The rim is round because of the remelted. The hornblende hasthe light yellow green to deep green pleochroism. Most grains are andedral with the remeltedrim. A few have the orthopyroxene coronas. All the orthopyroxene are anhedral and small insize about 0.02mm. Two glomerocrysts are noticed, one with biotite. Both are common withhornblende, plagioclase, opaque and orhtopyroxene. The hronblende is the same in clots andas phenocrysts. Therefore these two clots are the production of the hornblende breakdown.Mineral classification:Mineral^Mode^NamePlagioclase 10% A ndesiteHornblende^3%Orthopyroxene 1%Opaque^0.5%Glomerocryst^1.5%Groundmass 84%Sample number: MM-25Thin section:Plagioclase, hornblende, orthopyroxene, clinopyroxene and opaques phenocrysts, plagioclase,orthopyroxene, clinopyroxene and opaques glomerocryst are carried in a glassy groundmass ofplagioclase and opaques. Subhedral to anhedral plagioclase zonedly grew with the needle likeapatite, second-order interferene color clnopyroxene, first-order interference orthopyroxene,and higher-order white interference calcite. The melt embayed the plagioclase from the (010)93direction. The hornblende does't appear as plgioclase inclusion, but euhedral to subhedralphenocryst with size about 0.05mm. With the light yellow to deep brown pleochroism incommon, all the hornblende grains were oxided with the Fe-oxide rim. Some of the smallgrains were totally replaced by the opaques. The euhedral twinned clinopyroxene andsubhedral elongate clinopyroxene almost have the equal amount as the orthopyroxene. Theminerals in the glomerocryst are the same as phenocrysts, but disequilibrium is obviousexisting within this crystal clot. This infer that the glomerocryst is probably the production ofsettlement of plagioclase and pyroxene at the bottom of the magma chamber.Mineral classification:Mineral^Mode^NamePlagioclase 15 % AndesiteClinopyroxene^3%Orthopyroxene^3%Hornblende 2%Opaques^1%Groundmass 76%Sample number: MM-28Thin section:Plagioclase, hornblende, orthopyroxene and opaques phenocrysts, clinopyroxene and apatitemicrophenocrysts, plagioclase and orthopyroxene glomerocryst, plagioclase and opaquesglomerocryst are carried in the groundmass of plagioclase, orthopyroxene, clinopyroxene,opaques and glass. The plagioclase has the minor apatite and pyroxene inclusion, showingsieved structure. The sieved plagioclase varies with the remelted fringe wide from 0.1mm to0.01mm. Hornblende is altered to opaques. Euhedral to subhedral orthopyroxene is melted toplagioclase.Mineral classification:94Mineral^Mode^NamePlagioclase 10% AndesiteHornblende^1%Orthopyroxene^2%Opaques 1%Groundmass^86%Sample number: MM-35Thin section:Plagioclase, hornblende and orthoroxene phenocrysts, plagioclase, hornblende andorthopyroxene glomerocryst are carried in the groundmass of plagioclase, hornblende,orthopyroxene and glass. This sample is very vesicular. Plagioclase grew in the centre of avesicle. Hornblende has the light yellow to deep red pleochroism. All the ferromagnesianmineral, like hornblende and orthopyroxene were altered to opaques.Mineral classification:Mineral^Mode^NamePlagioclase 15 % AndesiteHornblende^3%Orthopyroxene^2%Opaques 1%Groundmass^79%A2.3 Capricorn Assemblage (CPA) CPA rock members are: MM-34 and MM-38.Sample number: MM-3495Thin section:Plagioclase, quartz, biotite and hornblende pheonocrysts, apatite microphenocryst are carriedin the groundmass of plagioclase, hornblende, orthopyroxene, opaques and glass. Theplagioclase shows two different groups. One group has very thin multiple twin lathe. Anothergroup does't have the obvious mutiple twin, but nice zoned. Two types of plagioclase havethe sieved texture, shown by the remelted fringe and overgrowth rim. Nice glass inclusion isfound in a few plagioclase grains. One plagioclase ia found including zircon microphenocryst.The quartz is embayed and corroded by the fluid and magmatic liquid. Fratured surface hintthe strong strain. Small quartz piece is also found in the groundmass. This infer that quartz isone of the early crystallizing mineral phase. The hornblende and biotite both are altered toopaques. The microphenocryst and groundmass orthopyroxene is also oxidized at the rim.Mineral classification:Mineral^Mode^NamePlagioclase 7% RhyodaciteQuartz^2%Hornblende 1%Biotite^1%groundmass 89%Sample number: MM-38Thin section:Plagioclase and hornblende phenocrysts are carried in the highly vesicular groundmass ofplagioclase. Subhedral plagioclase has an average size of 0.03mm with the largest up to0.1mm. This sample is hightly vesicular. Plagioclase is slightly zoned, embayed at the rim.Sieved texture is obviously in this sample, early disequilibrium shown as the irregular remeltedinside zone, later stable condition brought an overgrowth rim. The yellow green hornblendewas resorbed, with an opaque rim. The lackage of hornblende in the groundmass suggests that96the hornblende is possibly xenocryst. The crystallizing sequence is probably as plagioclase,then amplibole followed by magnetite.Mineral classification:Mineral^Mode^NamePlagioclase 3% RhyoliteHornblende^0.5%Vesicle 10%Groundmass^86%A2.4 Devastation Assemblage (DVA) DVA rock members are: MM-27 and MM-33.Sample number: MM-27Thin Section:Plagioclase, hornblende and orthopyroxene phenocrysts, hornblende, plagioclase and opaquesglomerocryst are carried in the groundmass of plagioclase, orthopyroxene, clinopyroxene andglass. The subhedral plagioclase has the size from 0.2mm to 005mm, slightly zoned, sieved,embayed. The hornblende is reacting to plagioclase and opaque with liquid, While it keeps thehornblende's shape. Orthopyroxene has slight pleochroism, from pale to pale pink, very lowinterference color.Mineral classification:Mineral^Mode^NamePlagioclase 15%Hornblende^3%Orthopyroxene^2%Groundmass 80%97Sample number: MM-33Thin section:Plgioclase, quartz and biotite phenocrysts are carried in the groundmass of plagioclase, biotiteand glass. Two groups of plagioclase appear. One show the sieved texture, with the remeltedfringe and overgrowth rim, slightly zoned, slightly twinned. Another shows the thin multipletwinning lathe, and zoned, but without the remelted fringe and overgrowth rim. Biotite isbeing broken down by the melt. Subhedral to anhedral quartz is remelted, embayed, by thedisequilibrium liquid.Mineral classification:Mineral^Mode^NamePlagioclase 12%Quartz^2%Biotite 1%Groundmass^85 %A2.5 Plinth Assemblage (PLA) PLA rock members are: MM-26, MM-29, MM-30, MM-39 and MM-40.Sample number: MM-26Thin section:Plagioclase, quartz, hornblende and biotite phenocrysts, apatite microphenocryst are carried inthe groundmass of plagioclase, hornblende, clinopyroxene, opaques and glass. Plagioclase iszoned and the grains of plagioclase tend to attach together by the sysnesis process. Thefractured quartz is embayed, corroded and even with the clinopyroxene coronas. Apatiteinclusions are visible in some grains of quartz. Anhedral biotite has the reaction rim of98hornblende and minor plagioclase. Hornblende is replacing the biotite. Microphenocrysts ofclinopyroxene and apatite are rich in the groundmass.Mineral classification:Mineral^Mode^NamePlagioclase 15% RhyodaciteQuartz^5%Hornblende 1%Biotite^3%Clinopyroxene^0.5%Apatite^0.1 %Groundmass 75.4%Sample number: MM-29Thin section:Plagioclase, quartz, biotite and hornblende phenocrysts, hornblende, plagioclase andclinopyroxene microphenocrysts are carried in the groundmass of plagioclase, hornblende andglass. There are two types of plagioclase. One is zoned with the thin lathe of multipletwinning. Another is zoned but with the coarse core. The quartz has the fractured surface,embayed and corroded by the melt. Euhedral to subhedral biotite has the size about 10.5m,some grains have the hornblende reaction rim. Hornblende has a much smaller size relativelyto biotite. Twinned clinopyroxene microphenocryst scattered in the groundmass.Mineral classification:Mineral^Mode^NamePlagioclase 20% RhyodaciteQuartz^3%Biotite 5%Hornblende^2%99Clinopyroxene^1%Groundmass 69%Sample number: MM-30Thin section:Plagioclase, quartz, biotite and hornblende phenocrysts, clinopyroxene and orthopyroxenemicrophenocrysts are carried in the groundmass of plagioclase, hornblende and glass. Theplagioclase has the sieved texture, shown as the remelted fringe and overgrowth rim. Thequartz has the clinopyroxene reaction rim. Hornblende is replacing the biotite, shown as thebiotite has the hornblende reaction rim.Mineral classification:Mineral^Mode^NamePlagioclase 12% AndesiteQuartz^7%Biotite 5%Hornblende^2.5%Groundmass 73.5%Sample number: MM-39Thin section:Plagioclase, quartz and biotite phenocrysts, clinopyroxene, opaque and zirconmicrophenocrysts are carried in the groundmass of plagioclase, opaques and glass. Apatiteinclusion in plagioclase is visible. The groundmass is in red representing the oxidized status ofthis sample. Biotite is breaking down to opaques.Mineral classification:Mineral^Mode^NamePlagioclase 12% Rhyodacite100Quartz^2%Biotite 3%Clinopyroxene^1%Opaques 1%Groundmass^81 %Sample number: MM-40Thin section:Plagioclase, quartz and biotite phenocrysts, hornblende microphenocryst are carried in thegroundmass of plagioclase, hornblende and glass. There are very small amount of pyroxene,mostly clinopyroxene on the groundmass. Biotite is breaking down to opaque. Bothplagioclase and quartz are remelted and with the overgrowth rim. The quartz is regarded asxenocrystal.Mineral classification:Mineral^Mode^NamePlagioclase 17% RhyodaciteBiotite^4%Quartz 4%Hornblende^1%Groundmass 74%A2.6 The Mosaic Assemblage (MOA) MOA rock members are: MM-17, MM-18, MM-19 and MM-37.Sample number: MM-17Thin section: (located at north exposure)101Plagioclase and olivine phenocrysts are carried in the groundmass of plagioclase, olivine andtitanaugite. Plagioclase is slightly zoned. Titanaugite fill between the plagioclase and olivinein the intergranular texture. Plagioclase is slightly zoned.Mineral classification:Mineral^Mode^NamePlagioclase 15 % Alkali olivine basaltOlivine^15 %Titanaugite 10%Groundmass^60%Sample number: MM-18Thin section: (located at south exposure)Big olivine and plagioclase phenocrysts are carried in the groundmass of plagioclase, olivineand clinopyroxene. Olivine is slightly zoned and with alteration rim.Mineral classification:Mineral^Mode^NamePlagioclase 15 % Alkali olivine basaltOlivine^15%Groundmass 70%Sample number: MM-19Thin section: (located at south exposure)Plagioclase and olivine phenocrysts are carried in the groundmass of plagioclase,clinopyroxene and olivine. Olivine has the clean surface and clear rim. Magnetite appears onthe surface of olivine.Mineral classification:Mineral^Mode^Name102Olivine^15%^Alkali olivine basaltPlagioclase 10%Ground mass^75%Sample number: MM-37Thin section: (located at west exposure)Plagioclase, olivine and clinopyroxene phenocrysts are carried in the groundmass ofplagioclase, olivine and clinopyroxene. Minor anorthoclase or sanidine with melting rim occurin the groundmass. Plagioclase has the melting overgrowth surface.Mineral classification:Mineral^Mode^NamePlagioclase 10% HawaiiteOlivine^5%Augite 10%groundmass^75 %103A3 Whole-Rock Chemical AnalysisThe operation conditions for the XRF analysis of major, minor and trace elements are inTable A3.1 and A3.2. Table 3.3 is the CIPW calculations for the four rocks in Mosaicassemblage in order to present the influence of iron status on CIPW results.104A4 Electron Microprobe AnalysisMineral chemistry was studied on the scanning electrical microprobe. The chemicalcomponents of phenocrystals provide the essential data for mass balance calculation. Thecompositions of plagioclase, biotite, amphibole and pyroxene were measured and the analyticalresults are presented in Table 4.1, 4.2, 4.3 and 4.4.105Table A3-1. XRF operating conditions for the analysis of major and minor elements.Fe pk Fe bkg Ti pk Ti bkg Ca pk Ca bkg K pk K bkg Si pk Si bkgLine^Ica^--^ka^--^ka^--^ka^--^ka^--20 57.49 56.00^86.12 91.00^113.12 110.40^136.70 132.15 109.24 113.30Target^Cr^Cr^Cr^Cr^Cr Cr^Cr^Cr^Cr^CrCrystal LIF200 LIF200 LIF200 LIF200 LIF20 LIF200 LIF200 LIF200 PET PET0Icv/mA^60/40 60/40^60/40 60/40 60/32 60/32^60/40 60/40^60/40 60/40Collimator^C^C^C^C^C^C^C^C^C^CCounter F^F F^F^F^F F^F F^FVacuum^On On^On On^On On^On^On^On OnGain 128^128^128^128^128^128^128^128^128^128Counterkv^830x2 830x2^830x2 830x2 830x2 830x2^830x2 830x2 830x2 830x2Lower level^200^200^200^200^200 200^200^200^200 200Window^800^800^800^800^800 800^800^800^800 800Count time^lOsec lOsec^lOsec lOsec^lOsec lOsec^lOsec lOsec^100sec lOsecElement^Al pk Al bkg P pk P bkg Mg pk Mg bkg Na pk Na bkg Mn pk Mnbkg Line^ka^--^ka^--^ka^--^ka^--^ka^--20 145.02 139.00^89.53 92.60 45.15 44.00^55.20 53.30^62.97 65.00Target^Cr^Cr^Cr^Cr^Cr Cr^Cr^Cr^Mo MoCrystal PET PET^PET PET TIAP TIAP^TIAP TIAP LIF20 LIF2000kv/mA^60/40 60/40^60/45 60/45 60/45 60/45^60/45 60/45^60/45 60/45Collimator^C^C^C^C^C^C^C^C^C^CCounter F^F F^F^F^F F^F^F^FVacuum^On^On^On On^On On^On^On^On OnGain 128^128^128^128^128^128^128^128^128^128Counter kv^830x2 830x2^830x2 830x2 830x2 830x2^830x2 830x2 830x2 830x2Lower level^200^200^200^200^200 200^200^200^200 200Window^800^800^800 800^800 800^800^800^800 800Count time^100sec lOsec^100sec lOsec^100sec lOsec^100sec lOsec^100sec lOsec* The XRF analyses were done on XRF spectrometer PW 1410 in the Department of GeologicalSciences of UBC.106Table A3-2 XRF operating conditions for the analysis of trace elements.Element Nb pk Nb -bkg Nb+bkg Y pk Y -bkg^Y+bkgLine La La --20 21.40 21.10 21.70 23.80 23.50^24.10Target Cr Cr Cr Cr Cr^CrCrystal LIF200 LIF200 LIF200 LIF200 LIF200^LIF200kv/mA 60/40 60/40 60/40 60/40 60/40^60/40Collimator F F F F F^FCounter S.I. S.I. S.I. S.I. S.I.^S.I.Vacuum On On On On On^OnGain 128 128 128 128 128^128Counter kv 1044 1044 1044 1044 1044^1044Lower level 280 280 280 280 280^280Window 420 420 420 420 420^420Count time 40sec 40sec 40sec 40sec 40sec^40secElement Zr pk Sr pk Rb pk Zr,Sr,Rb -bkg Zr,Sr,Rb+bkgLine La La La20 22.51 25.15 26.60 19.50 24.55Target Cr Cr Cr Cr CrCrystal LIF200 LIF200 LIF200 LIF200 LIF200kv/mA 60/40 60/40 60/40 60/40 60/40Collimator F F F F FCounter S.I. S.I. S.I. S.I. S.I.Vacuum On On On On OnGain 128 128 128 128 128Counter kv 1044 1044 1044 1044 1044Lower level 280 280 280 280 280Window 420 420 420 420 420Count time 40sec 40sec 40sec 40sec 40sec* The XRF analyses were carried on the XRF spectrometer in the Department ofGeological Sciences of UBC.Table A3-3 Effect of Fe2 +/Fe3 + ratio on the computed CIPW norm for the rocks of Mosaic assemblageFe state Fe total as Fe70/ Fe total as FeO Fe total as Fe.,0/ + FeOSample 17 18 19 37 17 18 19 37 17 18 19 37Ne 0.00 0.00 0.00 0.00 1.90 1.21 1.48 4.83 0.00 0.00 0.00 1.87Or 4.32 4.57 4.55 7.76 4.32 4.57 4.55 7.76 4.32 4.57 4.55 7.76P1 66.55 57.57 56.52 60.83 63.39 55.54 54.06 52.77 66.55 57.57 56.52 57.71Di 7.56 11.66 10.80 12.71 10.79 15.14 14.46 16.46 10.79 15.14 14.46 16.46Hy 8.11 12.11 11.97 0.11 0.00 0.00 0.00 0.00 0.04 2.02 1.08 0.0001 3.49 2.99 4.63 8.00 17.03 20.97 22.79 15.03 13.18 15.52 18.22 9.10II 0.23 0.25 0.25 0.22 1.85 1.99 2.08 2.09 1.85 1.99 2.08 2.09Hm 6.59 7.67 7.95 6.50 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00Tn 2.42 2.61 2.74 2.81 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00Ap 0.72 0.56 0.59 1.21 0.72 0.56 0.59 1.06 0.72 0.56 0.56 1.06Mt 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 2.55 2.63 2.51 3.95108Table A4-1 Electron microprobe analyses of PLA plagioclases.Sample MM-26A MM-26A MM-26B MM-26B MM-26B MM-26Oxides core rim core core rim averageSiO2 61.05 61.17 63.68 63.57 63.36 62.57Al203 24.40 24.31 22.78 22.86 22.94 23.46FeO 0.15 0.27 0.20 0.25 0.22 0.22MgO 0.01 0.01 0.00 0.00 0.00 0.00Na2O 7.55 7.48 8.32 8.14 8.11 7.92K2O 0.51 0.59 0.85 0.79 0.79 0.71CaO 6.28 5.98 4.13 4.40 4.69 5.10Total 99.95 99.81 99.96 100.00 100.11 99.97Structural Formula (Based on 8 Oxygens)Si4 + 2.7174 2.7249 2.8173 2.8119 2.8028 2.7750A13 + 1.2798 1.2763 1.1879 1.1917 1.1956 1.22623.9972 4.0012 4.0052 4.0036 3.9987 4.0012K 4 0.0291 0.0337 0.0478 0.0444 0.0445 0.0399Na+ 0.6514 0.6463 0.7132 0.6976 0.6958 0.6810Ca2 + 0.2995 0.2856 0.1958 0.2085 0.2225 0.2423mg2+ 0.0008 0.0004 0.0000 0.0003 0.0000 0.0003Fe2 + 0.0054 0.0099 0.0074 0.0091 0.0082 0.00800.9862 0.9759 0.9642 0.9600 0.9710 0.9715End MembersK 0.0291 0.0337 0.0478 0.0444 0.0445 0.0399Na 0.6541 0.6463 0.7132 0.6976 0.6958 0.6810Ca 0.2995 0.2856 0.1958 0.2085 0.2225 0.24230.9801 0.9656 0.9568 0.9506 0.9628 0.9632Others 0.0062 0.0104 0.0074 0.0094 0.0082 0.0083109Table A4.1 continued.Sample MM-29A MM-29A MM-29A MM-29A MM-29A MM-29C MM-29C MM-29B MM-29B MM-29B MM-29B MM-29B MM-29Oxides 1 2 3 4 5 core rim 1 2 3 4 5 averageSi02 57.45 57.36 60.26 61.37 62.66 61.99 63.10 62.34 62.56 61.34 63.43 63.49 61.44Al203 27.17 26.79 24.76 23.94 23.41 23.25 22.85 23.44 23.28 23.54 22.39 22.26 23.92FeO 0.22 0.18 0.31 0.16 0.15 0.22 0.15 0.23 0.14 0.09 0.17 0.17 0.18MgO 0.01 0.01 0.01 0.00 0.00 0.00 0.00 0.00 0.01 0.01 0.01 0.00 0.01Na2O 5.97 6.05 7.14 7.46 8.10 8.00 8.22 7.92 8.05 7.70 8.47 8.53 7.63K20 0.35 0.28 0.45 0.64 0.73 0.72 0.71 0.82 0.71 0.74 0.95 0.94 0.67CaO 9.12 8.83 6.58 5.79 4.94 4.84 4.39 4.93 4.87 4.97 3.76 3.67 5.56Total 100.29 99.49 99.51 99.38 99.99 99.02 99.43 99.68 99.62 98.38 99.19 99.06 99.42Structural Formula (Based on 8 Oxygens)Si4 + 2.5690 2.5818 2.6957 2.7420 2.7780 2.7756 2.8070 2.7739 2.7822 2.7625 2.8279 2.8332A13 + 1.4319 1.4211 1.3055 1.2608 1.2235 1.2269 1.1978 1.2291 1.2205 1.2394 1.1764 1.17084.0009 4.0029 4.0012 4.0027 4.0015 4.0025 4.0048 4.0030 4.0027 4.0119 4.0043 4.0039K+ 0.0198 0.0164 0.0259 0.0368 0.0414 0.0409 0.0402 0.0465 0.0404 0.0428 0.0542 0.0538Na+ 0.5179 0.5279 0.6196 0.6467 0.6963 0.6946 0.7090 0.6835 0.6938 0.6723 0.7319 0.7382Ca2 + 0.4369 0.4256 0.3156 0.2773 0.2346 0.2323 0.2091 0.2350 0.2320 0.2397 0.1798 0.1754mg2+ 0.0004 0.0008 0.0008 0.0000 0.0002 0.0003 0.0000 0.0000 0.0008 0.0004 0.0009 0.0000Fe2 + 0.0084 0.0066 0.0116 0.0061 0.0054 0.0083 0.0059 0.0087 0.0053 0.0034 0.0062 0.00620.9834 0.9770 0.9734 0.9669 0.9780 0.9765 0.9642 0.9737 0.9723 0.9587 0.9730 0.9737End MembersK 0.0198 0.0164 0.0259 0.0368 0.0414 0.0409 0.0402 0.0465 0.0404 0.0428 0.0542 0.0538Na 0.5179 0.5276 0.6196 0.6467 0.6963 0.6946 0.7090 0.6835 0.6938 0.6723 0.7319 0.7382Ca 0.4369 0.4256 0.3156 0.2773 0.2346 0.2323 0.2091 0.2350 0.2320 0.2397 0.1798 0.17540.9746 0.9696 0.9611 0.9608 0.9723 0.9678 0.9583 0.9650 0.9662 0.9548 0.9659 0.9674Others 0.0089 0.0074 0.0123 0.0061 0.0057 0.0086 0.0059 0.0087 0.0060 0.0039 0.0071 0.0062110Table A4.1 continued.SampleOxidesMM-30B MM-30B MM-30B MM-30B MM-30B MM-30B MM-30B MM-30B MM-30B MM-30B MM-30B MM-30A1^2^3^4^5^6^7^8^9^10^11^coreMM-30AMM-30AMM-30C MM-30C MM-30middle^rim^core^rim^averageSiO2 63.77 54.88 56.42 56.48 60.86 62.87 56.83 62.48 59.35 61.86 63.74 62.67 62.67 59.52 62.97 62.60 60.62Al203 22.64 28.39 27.62 27.43 24.71 23.21 26.96 23.51 25.73 24.07 22.82 23.14 23.43 24.99 22.93 23.48 24.69FeO 0.24 0.18 0.18 0.20 0.22 0.19 0.15 0.22 0.22 0.12 0.22 0.17 0.23 0.45 0.13 0.17 0.21MgO 0.01 0.00 0.00 0.02 0.00 0.02 0.01 0.01 0.02 0.00 0.02 0.01 0.02 0.05 0.00 0.00 0.01Na2O 8.18 5.03 5.62 5.74 7.34 7.92 5.97 7.77 6.72 7.48 8.15 8.17 8.15 5.61 7.98 7.74 7.10K2O 0.82 0.28 0.251 0.29 0.53 0.73 0.30 0.75 0.45 0.58 0.96 0.82 0.76 0.51 0.89 0.69 0.60CaO 4.42 10.99 9.68 9.76 6.11 4.96 8.63 5.14 7.77 5.60 4.06 4.74 4.56 8.35 4.62 4.92 6.52Total 100.08 99.69 99.75 99.92 99.78 99.91 98.85 99.88 100.24 99.71 99.96 99.71 99.83 99.48 99.52 99.60 99.74Si4 + 2.8192 2.4822 2.5389 2.5405 2.7112 2.7878 2.5732 2.7733 2.6450 2.7494 2.8191 2.7868 2.7814 2.6696 2.8013 2.7812 2.7045A13 + 1.1798 1.5135 1.4646 1.4541 1.2973 1.2128 1.4386 1.2300 1.3515 1.2610 1.1893 1.2125 1.2255 1.3209 1.2023 1.2295 1.29823.9991 3.9957 4.0035 3.9946 4.0085 4.0006 4.0118 4.0032 3.9965 4.0104 4.0083 3.9994 4.0070 3.9904 4.0036 4.0107 4.0027Structural Formula (Based on 8 Oxygens)K+ 0.0462 0.0159 0.0144 0.0167 0.0302 0.0413 0.0176 0.0422 0.0254 0.0329 0.0544 0.0465 0.0428 0.0292 0.0503 0.0393 0.0342Na+ 0.7013 0.4412 0.4901 0.5006 0.6342 0.6807 0.5241 0.6691 0.5802 0.6445 0.6985 0.7042 0.7014 0.4881 0.6881 0.6669 0.6139Ca2 + 0.2093 0.5301 0.4665 0.4703 0.2915 0.2358 0.4186 0.2445 0.3710 0.2667 0.1926 0.2258 0.2169 0.4013 0.2203 0.2344 0.3115Mg2 + 0.0003 0.0002 0.0000 0.0016 0.0000 0.0015 0.0007 0.0003 0.0010 0.0000 0.0010 0.0004 0.0016 0.0032 0.0000 0.0000 0.0007Fe2 + 0.0089 0.0066 0.0066 0.0077 0.0081 0.0071 0.0058 0.0084 0.0081 0.0046 0.0080 0.0062 0.0085 0.0167 0.0048 0.0062 0.00770.9661 0.9941 0.9775 0.9968 0.9641 0.9664 0.9668 0.9645 0.9857 0.9486 0.9546 0.9831 0.9712 0.9385 0.9635 0.9467 0.9680End MembersK 0.0462 0.0159 0.1444 0.0167 0.0302 0.0413 0.0176 0.0422 0.0254 0.0329 0.0544 0.7042 0.7014 0.4881 0.6881 0.6669 0.6139Na 0.7013 0.4412 0.4901 0.5006 0.6342 0.6807 0.5241 0.6691 0.5802 0.6445 0.6985 0.2258 0.2169 0.4013 0.2203 0.2344 0.3115Ca 0.2093 0.5301 0.4663 0.4703 0.2915 0.2358 0.4186 0.2445 0.3710 0.2667 0.1926 0.9765 0.9611 0.9186 0.9586 0.9405 0.95960.9568 0.9873 0.9710 0.9876 0.9560 0.9578 0.9603 0.9558 0.9766 0.9440 0.94550.0067 0.0101 0.0200 0.0048 0.0062 0.0084Others 0.0093 0.0068 0.0066 0.0093 0.0081 0.0086 0.0065 0.0087 0.0091 0.0046 0.0090111Table A4.1 continued.SampleOxidesMM-39B1MM-39B2MM-39B3MM-39B4MM-39B5MM-39B6MM-39B7MM-39B8MM-39C1MM-39C2MM-39C3MM-39C MM-39C4^5Si02 62.37 62.08 61.30 61.35 61.03 60.54 59.46 62.41 62.91 62.43 62.80 63.23 61.63Al203 23.51 24.16 24.04 24.73 24.32 24.48 25.47 23.64 23.33 23.27 23.27 22.84 24.18FeO 0.20 0.12 0.17 0.267 0.19 0.28 0.23 0.15 0.19 0.14 0.15 0.18 0.15MgO 0.00 0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00Na2O 7.98 7.84 7.65 7.89 7.81 7.52 7.55 7.64 7.94 8.05 8.25 8.40 8.08K20 0.65 0.67 0.76 0.84 0.81 0.72 0.80 1.01 0.79 0.79 0.77 0.82 0.81CaO 5.18 5.43 5.10 4.62 4.92 5.44 5.32 4.80 4.74 4.67 4.51 4.08 4.31Total 99.90 100.32 99.01 99.70 99.08 98.98 98.84 99.66 99.89 99.34 99.76 99.56 99.17Structural Formula (Based on 8 Oxygens)Si4 + 2.7693 2.7465 2.7458 2.7301 2.7343 2.7185 2.6775 2.7747 2.7885 2.7837 2.7882 2.8097 2.7529A13 + 1.2304 1.2596 1.2681 1.2968 1.2840 1.2955 1.3515 1.2388 1.2187 1.2228 1.2179 1.1963 1.27303.9996 4.0061 4.0149 4.0269 4.0183 4.0140 4.0290 4.0135 4.0071 4.0065 4.0060 4.0060 4.0259K+ 0.0366 0.0379 0.0433 0.0479 0.0465 0.0410 0.0462 0.0576 0.0448 0.0452 0.0436 0.0466 0.0461Na+ 0.6874 0.6727 0.6644 0.6811 0.6783 0.6547 0.6591 0.6585 0.6825 0.6960 0.7101 0.7241 0.6998Ca2 + 0.2465 0.2574 0.2448 0.2204 0.2361 0.2619 0.2565 0.2289 0.2249 0.2230 0.2147 0.1943 0.2064mg2 + 0.0000 0.0004 0.0000 0.0000 0.0000 0.0000 0.0003 0.0001 0.0000 0.0000 0.0000 0.0000 0.0003pe2+ 0.0076 0.0046 0.0063 0.0100 0.0070 0.0103 0.0086 0.0057 0.0070 0.0053 0.0057 0.0067 0.00550.9781 0.9730 0.9587 0.9593 0.9680 0.9679 0.9708 0.9505 0.9591 0.9694 0.9742 0.9717 0.9580End MembersK 0.0366 0.0379 0.0433 0.0479 0.0465 0.0410 0.0462 0.0573 0.0448 0.0452 0.0436 0.0466 0.0461Na 0.6874 0.6727 0.6644 0.6811 0.6783 0.6547 0.6591 0.6585 0.6825 0.6960 0.7101 0.7241 0.6998Ca 0.2465 0.2574 0.2448 0.2204 0.2361 0.2619 0.2565 0.2289 0.2249 0.2230 0.2147 0.1943 0.20640.9705 0.9679 0.9525 0.9493 0.9610 0.9575 0.9618 0.9447 0.9521 0.9641 0.9684 0.9650 0.9522Others 0.0076 0.0051 0.0063 0.0100 0.0070 0.0103 0.0090 0.0059 0.0070 0.0053 0.0057 0.0067 0.0058112Table A4.1 continued.SampleOxidesMM-39C6MM-39C7MM-39C MM-39D MM-39D MM-39A MM-39A8^core^rim^rim^rimMM-39A MM-39A MM-39A MM-39A MM-39A MM-391^2^3^4^5^averageSiO2 61.33 62.31 61.94 62.74 62.44 62.48^62.30 62.77 60.54 62.68 61.77 61.69 61.94Al203 23.73 22.87 23.30 23.730 23.06 23.32^23.82 23.15 25.01 23.54 24.23 24.08 23.80FeO 0.16 0.24 0.20 0.22 0.18 0.11^0.15 0.24 0.12 0.23 0.18 0.20 0.19MgO 0.00 0.00 0.00 0.00 0.00 0.01^0.00 0.02 0.00 0.00 0.00 0.00 0.00Na2O 7.99 8.34 7.95 7.81 7.81 7.96^7.89 8.01 7.09 8.15 7.66 7.76 0.81K2O 0.69 0.81 1.12 0.68 1.20 0.86^0.86 1.00 0.55 0.78 0.64 0.80 0.81CaO 5.29 4.38 4.580 5.19 4.88 4.97^5.42 4.85 6.85 4.61 5.51 5.06 4.99Total 99.14 98.93 99.11 100.37 99.56 99.72^100.45 100.04 100.16 100.00 100.00 99.59 99.61Structural Formula (Based on 8 Oxygens)Si4 + 2.7480 2.7922 2.7750 2.7704 2.7852 2.7789^2.7559 2.7854 2.6913 2.7779 2.7413 2.7485Al3 + 1.2530 1.2079 1.2301 1.2348 1.2122 1.2221^1.2421 1.2106 1.3106 1.2297 1.2671 1.26434.0010 4.0001 4.0051 4.0052 3.9974 4.0010^3.9980 3.9960 4.0019 4.0076 4.0084 4.0127Structural Formula (Based on 8 Oxygens)K+ 0.0397 0.0462 0.0644 0.0383 0.0680 0.0491^0.0485 0.0567 0.0312 0.0444 0.0364 0.0457Na+ 0.6901 0.7232 0.6908 0.6688 0.6753 0.6862^0.6770 0.6889 0.6112 0.7002 0.6593 0.6706Ca2 + 0.2540 0.2103 0.2198 0.2457 0.2333 0.2367^0.2570 0.2306 0.3263 0.2190 0.2620 0.2414mg2+ 0.0000 0.0000 0.0000 0.0000 0.0000 0.0008^0.0000 0.0012 0.0000 0.0000 0.0003 0.0000Fe2 + 0.0059 0.0089 0.0077 0.0080 0.0065 0.0042^0.0056 0.0088 0.0044 0.0087 0.0068 0.00750.9897 0.9887 0.9826 0.9608 0.9831 0.9770^0.9880 0.9862 0.9731 0.9723 0.9649 0.9651End MembersK 0.0397 0.0462 0.0644 0.0383 0.0680 0.0491^0.0485 0.0567 0.0312 0.0444 0.0364 0.0457Na 0.6901 0.7232 0.6908 0.6688 0.6753 0.6862^0.6770 0.6889 0.6112 0.7002 0.6593 0.6706Ca 0.2540 0.2103 0.2198 0.2457 0.2333 0.2367^0.2570 0.2306 0.3263 0.2190 0.2620 0.24140.9837 0.9798 0.9750 0.9528 0.9766 0.9720^0.9825 0.9761 0.9687 0.9636 0.9577 0.9577Others 0.0059 0.0089 0.0077 0.0080 0.0065 0.0050^0.0056 0.0100 0.0044 0.0087 0.0072 0.0075113Table A4-2 Electron microprobe analyses of PLA amphiboles.Sample MM-293 MM-293 MM-299 MM-299 MM-29 MM-265 MM-265 MM-26Oxides 1 2 1 2 average 1 2 averageSi02 47.80 49.86 41.39 40.49 44.88 49.95 48.94 49.44Al203 5.99 5.21 12.96 13.23 9.35 5.33 5.93 5.63TiO2 0.97 0.90 3.05 2.98 1.98 1.02 0.97 1.00FeO 12.16 12.16 8.70 8.58 10.40 11.45 11.46 11.46MgO 16.08 16.24 15.94 15.88 16.03 17.00 16.54 16.77MnO 0.81 0.86 0.06 0.74 0.62 0.82 0.78 0.80Na20 1.41 1.35 2.91 2.89 2.14 1.30 1.43 1.36K2O 0.30 0.22 0.65 0.65 0.46 0.25 0.32 0.28CaO 11.41 11.43 12.33 12.32 11.87 11.28 11.28 11.28F 0.42 0.32 0.32 0.45 0.38 0.34 0.46 0.40Cl 0.13 0.06 0.02 0.02 0.06 0.08 0.13 0.10Total 97.47 98.59 98.34 97.58 97.99 98.88 98.25 98.66Structural Formula (Based on 24 (0, OH, F, and Cl))Si4 + 6.9192 7.2721 6.0533 5.8783 7.2746 7.0991 7.0991A13-1- (IV) 1.0213 0.7279 1.9467 2.1217 0.7254 0.9009 0.90097.9405 8.0000 8.0000 8.0000 8.0000 8.0000 8.0000mg2+ 3.4693 3.5300 3.4748 3.4369 3.6913 3.5754 3.5754Fe2 + 1.4717 1.4828 1.0635 1.0424 1.3945 1.3909 1.3909Mn2 + 0.0991 0.1064 0.0075 0.0091 0.1016 0.0961 0.0961A13 +(VI) 0.0000 0.1679 0.2865 0.1426 0.1899 0.1135 0.1135Ti4 + 0.1059 0.0986 0.3359 0.3256 0.1113 0.1063 0.1063Ca2 + 1.7700 1.7866 1.9322 1.9163 1.7605 1.7532 1.7532Na+(C) 0.0840 0.0000 0.0000 0.1272 0.0000 0.0000 0.00007.0000 7.1723 7.1003 7.0000 7.2491 7.0353 7.0353Na+(A) 0.3117 0.3812 0.8264 0.6875 0.3666 0.4019 0.4019K+ 0.0552 0.0399 0.1220 0.1207 0.0466 0.0591 0.05910.3669 0.4211 0.9484 0.8082 0.4131 0.4609 0.4609OH- 2.4421 1.3650 1.6215 2.3455 1.1547 1.6948 1.6948F- 0.1364 0.1029 0.1050 0.1465 0.1096 0.1500 0.1500Cl- 0.0265 0.0121 0.0032 0.0037 0.0151 0.0258 0.02582.6050 1.4800 1.7297 2.4956 1.2794 1.8707 1.87070 21.3950 22.5200 22.7297 21.5044 22.7206 22.1293 22.1293114Table A4-3 Electron microprobe analyses of PLA pyroxenes.Sample MM-395 MM-395 MM- MM-39 MM-291 MM-291 MM-29 MM-407 MM-407 MM-407 MM-40395Oxides^2^3^4^average 1^2^average 1^2^3^averageSiO2^54.27 54.50 54.09 54.28 52.07 50.16^51.12 54.20 54.48 52.80 53.82Al203^2.66^1.55^0.92^1.71^2.23^2.95^2.59^1.65^1.24^3.22^2.04FeO^12.61^13.54^15.06^13.74 5.79^5.58^5.68^13.82^13.98^13.16^13.66MgO^29.18^28.74 28.21 28.71 16.36^15.87^16.11 28.28^28.85^28.86^28.66MnO^0.32^0.41^0.39^0.37^0.14^0.14^0.14 0.43^0.37^0.34^0.38TiO2^0.20^0.18^0.20^0.19 0.49^0.572 0.53 0.20^0.06^0.22^0.16Na2O^0.07^0.02^0.03^0.04 0.31^0.34^0.32 0.02^0.03^0.02^0.02CaO^1.23^1.33^1.38^1.31^22.11^22.24^22.17 1.23^1.24^1.25^1.24Total^100.54 100.26 100.27 100.36 99.50 97.84^98.67 99.83^100.26 99.88^99.99Structural Formula (Based on 6 Oxygens)Si4 +^1.9240 1.9469 1.9482^1.9268 1.8931^1.9469 1.9502 1.8935A13 ±(IV) 0.0760 0.0531 0.0392 0.0732 0.1069 0.0531 0.0498 0.10652.0000 2.0000 1.9875^2.0000 2.0000^2.0000 2.0000 2.0000Ca2 +^0.0467 0.0508 0.0532^0.8765 0.8992^0.0474 0.0475 0.0482Mg2 +^1.5419 1.5303 1.5144 0.9020 0.8928 1.5141 1.5395 1.5428Fe2 +^0.3740 0.4047 0.4536^0.1793 0.1760^0.4153 0.4198 0.3946M3+ (VI) 0.0351 0.0120 0.0000 0.0240 0.0242 0.0169 0.0027 0.0298Mn2 +^0.0095 0.0123 0.0119^0.0045 0.0043^0.0131 0.0112 0.0103Na2 ÷^0.0050 0.0012 0.0023 0.0220 0.0252 0.0012 0.0020 0.0015Ti4 +^0.0054 0.0049 0.0053^0.0137 0.0162^0.0054 0.0015 0.00602.0176 2.0163 2.0406 2.0220 2.0378 2.0133 2.0231 2.0332End MembersMg^0.7710 0.7652 0.7572^0.4510 0.4464^0.7570 0.7698 0.7714Fe^0.1870 0.2023 0.2268 0.0896 0.0880 0.2077 0.2093 0.1973Ca^0.0233 0.0254 0.0266^0.4383 0.4496^0.0237 0.0238 0.0241 0.9813 0.9929 1.0106 0.9789 0.9840 0.9884 1.0029 0.9928Others^0.0275 0.0153 0.0097^0.0321 0.0349^0.0183 0.0087 0.0238115Table A4-4 Original electron microprobe analyses of PLA biotites.Sample MM-265 MM-265 MM-264 MM-264 MM-264 MM-264 MM-26Oxides core core rim rim core core averageSi02 36.87 36.96 37.44 37.04 37.30 37.32 37.16Al203 13.47 13.41 13.33 13.39 13.60 13.37 13.43FeO 16.76 17.65 18.02 17.16 17.48 17.15 17.37MgO 13.69 13.81 14.02 14.07 13.85 14.03 13.91MnO 0.38 0.40 0.38 0.31 0.34 0.40 0.37TiO2 4.16 4.21 4.32 4.36 4.33 4.23 4.27K20 8.80 9.07 9.05 8.83 8.62 8.68 8.84CaO 0.03 0.02 0.01 0.00 0.07 0.03 0.03Na20 0.54 0.53 0.56 0.57 0.61 0.58 0.56F 0.53 0.40 0.52 0.68 0.32 0.38 0.47Cl 0.13 0.25 0.16 0.18 0.18 0.14 0.17Total 95.37 96.62 97.80 96.59 96.69 96.30 96.56Structural Formula (Based on 24 (0, OH, F, and Cl))Si4 + 5.4779 5.6087 5.7512 5.5837 5.6353 5.6075Al3 + (IV) 2.3594 2.3913 2.2488 2.3794 2.3647 2.36797.8373 8.0000 8.0000 7.9631 8.0000 7.9754A0+(V1) 0.0000 0.0074 0.1640 0.0000 0.0565 0.0000Mg2 + 3.0334 3.1234 3.2104 3.1619 3.1195 3.1438Fe2 + 2.0819 2.2395 2.3152 2.1634 2.2082 2.1548Ti4 + 0.4650 0.4808 0.4988 0.4941 0.4922 0.4782NIn2-1- 0.0479 0.0513 0.0491 0.0401 0.0436 0.05035.6283 5.9024 6.2375 5.8595 5.9200 5.8266K+ 1.6689 1.7548 1.7723 1.6972 1.6621 1.6644Na+ 0.1545 0.1570 0.1674 0.1670 0.1785 0.1704Ca2 + 0.0053 0.0030 0.0018 0.0000 0.0120 0.00451.8288 1.9148 1.9416 1.8643 1.8526 1.8393OH- 4.5878 3.3286 2.2524 3.4327 3.3325 3.7084F- 0.1739 0.1368 0.1763 0.2265 0.1070 0.1259Cl- 0.0267 0.0525 0.0336 0.0381 0.0366 0.02904.7883 3.5178 2.4623 3.6972 3.4761 3.86340 19.2117 20.4822 21.5377 20.3028 20.5239 20.1366116Table A4.4 continued.SampleOxidesMM-297coreMM-297coreMM-297rimMM-297rimMM-294coreMM-294coreMM-294 MM-294 29rim^rim^averageSiO2 37.31 36.76 37.15 36.90 36.95 37.13 36.80 36.76 36.97Al2°3 13.23 13.34 13.43 13.42 13.29 13.57 13.32 13.32 13.36FeO 16.45 16.07 16.33 15.74 17.18 17.58 17.11 17.27 16.72MgO 14.50 14.48 14.55 14.71 13.81 13.75 13.72 14.01 14.19MnO 0.32 0.30 0.35 0.29 0.28 0.33 0.38 0.39 0.33TiO2 4.30 4.39 4.21 4.34 4.34 4.232 4.27 4.31 4.3K2O 8.58 8.25 8.83 8.65 8.55 8.64 8.59 8.38 8.56CaO 0.07 0.11 0.10 0.04 0.09 0.04 0.08 0.01 0.07Na2O 0.54 0.53 0.60 0.70 0.52 0.51 0.50 0.49 0.55F 0.73 0.72 0.60 0.53 0.48 0.64 0.62 0.64 0.62Cl 0.20 0.20 0.20 0.23 0.19 0.18 0.15 0.18 0.19Total 96.23 95.16 96.35 95.56 95.69 96.60 95.54 95.77 95.86Structural Formula (Based on 24 (0, OH, F, and Cl))Si4 + 5.5742 5.4162 5.5687 5.4685 5.5116 5.6000 5.4800 5.4865Al3 +(IV) 2.3302 2.3157 2.3719 2.3450 2.3370 2.4000 2.3373 2.34197.9044 7.7319 7.9407 7.8135 7.8486 8.0000 7.8173 7.8284Al3 +(VI) 0.0000 0.0000 0.0000 0.0000 0.0000 0.0114 0.0000 0.0000Mg2 + 3.2288 3.1806 3.2507 3.2505 3.0170 3.0914 3.0459 3.1154Fe2 + 2.0556 1.9805 2.0470 1.9514 2.1429 2.2176 2.1311 2.1553114+ 0.4827 0.4868 0.4744 0.4839 0.4873 0.4799 0.4778 0.4839Mn2 + 0.0412 0.0375 0.0448 0.0363 0.0352 0.0426 0.0485 0.04965.8053 5.6854 5.8186 5.7221 5.7363 5.8428 5.7034 5.8043K+ 1.6343 1.5508 1.6892 1.6354 1.6275 1.6625 1.6311 1.5960Na+ 0.1574 0.1525 0.1747 0.2018 0.1513 0.1482 0.1452 0.1428Ca2+ 0.0114 0.0177 0.0155 0.0069 0.0139 0.0070 0.0123 0.00181.8030 1.7210 1.8795 1.8441 1.7926 1.8177 1.7886 1.7406OH 3.7581 4.7513 3.6455 4.3821 4.2949 3.4162 4.4296 4.2077F 0.2410 0.2358 0.2015 0.1744 0.1604 0.2146 0.2076 0.2124Ct 0.0408 0.0400 0.0409 0.0472 0.0394 0.0368 0.0310 0.03793.0398 5.0270 3.8879 4.6037 4.4947 3.6676 4.6679 4.45810 19.9602 18.9730 20.1121 19.3963 19.5053 20.3324 19.3321 19.5419117Table A4.4 continued.SampleOxidesMM-3091MM-309 MM-3062^1MM-306 MM-302^averageMM-392coreMM-392 MM-392rim^rimMM-39averageSi02 36.29 36.02 36.74 37.24 36.57 36.34 36.62 37.13 36.69Al203 12.49 13.02 13.60 13.73 13.21 14.19 13.20 13.36 13.58FeO 15.99 16.44 16.94 17.03 16.60 16.85 17.14 17.16 17.05MgO 13.73 13.74 14.01 13.53 13.75 13.96 13.87 14.06 13.96MnO 0.37 0.37 0.36 0.48 0.40 0.25 0.51 0.29 0.35TiO2 4.05 4.18 4.22 4.40 4.21 4.32 4.27 4.35 4.31K20 7.93 7.90 8.79 8.67 8.32 8.68 9.02 8.94 8.88CaO 0.17 0.29 0.04 0.05 0.14 0.03 0.00 0.01 0.01Na20 0.59 0.59 0.71 0.70 0.65 0.58 0.64 0.58 0.60F 1.09 0.50 0.36 0.48 0.61 1.38 0.98 0.98 1.11Cl 0.16 0.15 0.24 0.19 0.19 0.16 0.16 0.14 0.15Total 92.85 93.30 96.01 99.50 95.39 96.73 96.40 96.99 96.71Structural Formula (Based on 24 (0, OH, F, and Cl))Si4 + 5.2084 5.2148 5.5111 5.6046 5.4610 5.5197 5.6166A13 +(IV) 2.1128 2.2215 2.4037 2.3954 2.5128 2.3456 2.38287.3212 7.4363 7.9149 8.0000 7.9738 7.8654 7.99940.0000 0.0000 0.0000 0.0395 0.0000 0.0000 0.0000A13 +(VI) 2.9372 2.9653 3.1314 3.0352 3.1278 3.1155 3.1693Mg2 + 1.9187 1.9904 2.1246 2.1437 2.1182 2.1612 2.1714Fe2+ 0.4368 0.4547 0.4764 0.4980 0.4883 0.4845 0.4947Ti4 + 0.0452 0.0454 0.0454 0.0612 0.0321 0.0648 0.0374Mn2 + 5.3379 5.4558 5.7779 5.7776 5.7663 5.8259 5.87271.4513 1.4597 1.6817 1.6637 1.6637 1.7341 1.7246K+ 0.1636 0.1657 0.2078 0.2030 0.1677 0.1842 0.1692Na+ 0.0265 0.0456 0.0065 0.0081 0.0045 0.0000 0.0011Ca2+ 1.6413 1.6710 1.8960 1.8748 1.8359 1.9213 1.89506.8534 6.5585 3.9879 3.5168 3.2789 3.6125 3.0390OH' 0.3486 0.1625 0.1214 0.1618 0.4611 0.3299 0.3301F- 0.0319 0.0297 0.0503 0.0398 0.0329 0.0335 0.0287Cl- 7.2339 6.7507 4.1597 3.7184 3.7728 3.9785 3.39780 16.7661 17.2493 19.8403 20.2816 20.2272 20.0241 20.6022118Table A4.4 continued.Sample MM-406 MM-406 MM-406 MM-406 MM-405 MM-405 MM-40Oxides rim rim core core core core averageSiO2 43.07 43.88 41.55 43.28 36.23 36.98 40.83Al203 10.97 10.49 13.09 11.94 13.55 13.29 12.22FeO 11.83 11.31 12.72 11.15 16.62 16.77 13.40MgO 14.07 14.64 13.22 14.15 13.92 14.19 14.03MnO 0.38 0.22 0.24 0.33 0.38 0.35 0.32TiO2 2.52 2.36 2.36 2.30 4.43 4.27 3.04K2O 0.26 0.28 0.51 0.57 8.76 8.80 3.20CaO 11.46 11.92 11.76 11.26 0.06 0.03 7.75Na2O 2.37 2.39 2.24 2.02 0.64 0.68 1.73F 0.17 0.23 0.32 0.13 0.75 0.46 0.34Cl 0.01 0.05 0.06 0.09 0.15 0.16 0.09Total 97.10 97.77 98.08 97.22 95.48 95.99 96.94Structural Formula (Based on 24 (0, OH, F, and Cl))Si4 + 6.2533 6.3973 6.1316 6.2623 5.3919 5.5376A13-1- (IV) 1.7467 1.6027 1.8684 1.7377 2.3771 2.34488.0000 8.0000 8.0000 8.0000 7.7690 7.8823A131- (VI) 0.1300 0.2004 0.4082 0.2983 0.0000 0.0000Mgt -1- 3.0440 3.1826 2.9073 3.0538 3.0876 3.1666Fe2 + 1.4362 1.3787 1.5699 1.3493 2.0686 2.0991Ti4 + 0.2757 0.2583 0.2614 0.2505 0.4957 0.4810Mn2 + 0.0472 0.0273 0.0294 0.0402 0.0474 0.04474.9331 5.0472 5.1761 4.9891 5.6992 5.7914K+ 0.0477 0.0513 0.0957 0.1056 1.6634 1.6815Nal- 0.6675 0.6768 0.6426 0.5680 0.1863 0.1984Ca2+ 1.7822 1.8623 1.8646 1.7458 0.0091 0.00472.4975 2.5904 2.6029 2.4194 1.8588 1.8846OH- 2.8042 2.1642 1.8914 2.6736 4.4892 4.0076F- 0.0546 0.0756 0.1046 0.0412 0.2476 0.1544Cl- 0.0020 0.0109 0.0120 0.0186 0.0310 0.03482.8608 2.2507 2.0080 2.7334 4.7679 4.19680 21.1392 21.7493 21.9920 21.2666 19.2321 19.80321193Parent Daughter Plag Opx Xenolith Relative_Error57.76 58.65 53.35 53.86 72.01 1.00.89 0.77 0.04 0.00 0.20 1.017.18 17.42 29.51 2.46 14.44 1.06.46 6.06 0.93 13.64 1.78 1.05.61 4.00 0.10 28.53 0.89 1.06.90 6.41 11.80 1.48 2.61 1.03.69 3.96 4.42 0.00 3.96 1.01.22 1.49 0.23 0.00 3.00 1.00.29 0.30 0.00 0.00 0.13 1.0An Assimilation Example from Literature.0.0059OxideSi02TiO2Al203FeOMgOCaONa2OK20P2050.90Input Data File In.DatDescription of Data File IN.DATLine 1:^(Characters) The title of the project.Line 2:^(Float) The interval of amounts of minerals to map the x 2 and R2 distributions aroundthe optimal solution.Line 3:^(Integers) The numbers of the oxide components (n) and mineral phases (m).Line 4:^(Character strings) Labels (such as "oxide"), the name of the parent rock, the name ofderivative rock, and the names of (m) mineral phases. Note: Each character stringmust not contain any spaces.Line 5 - (4+n): (One character string, and m+3 floats) The name of oxide component, thecompositions of the parent and derivative rocks, the compositions of the mineralphases, and the relative analytical error for each oxide.Line (5+n):^(Float) The confidence level used to calculate the confidence region.120Example of Ouput of Program MASSBAL.C. Filename = Out.dat.***********************************An Assimilation Example from Literature***********************************INPUT PARAMETRSnumber of oxide = 9,^number of mineral_phase = 3,^degrees_of freedom = 6Oxide Parent Daughter Plag Opx Xenolith Relative ErrorSiO2 57.76 58.65 53.35 53.86 72.01 1.00TiO2 0.89 0.77 0.04 0.00 0.20 1.00Al203 17.18 17.42 29.51 2.46 14.44 1.00FeO 6.46 6.06 0.93 13.64 1.78 1.00MgO 5.61 4.00 0.10 28.53 0.89 1.00CaO 6.90 6.41 1.80 1.48 2.61 1.00Na2O 3.69 3.96 4.42 0.00 3.96 1.00K2O 1.22 1.49 0.23 0.00 3.00 1.00P2O5 0.29 0.30 0.00 0.00 0.13 1.00The Parameters of the Compositional Difference Matrix# oxide ml-d m2-d m3-d p-d1 SiO2 (-5.30 -4.79 13.36) = -0.892 TiO2 (-0.73 -0.77 -0.57) = 0.123 Al203 (12.09 -14.96 -2.98) = -0.244 FeO (-5.13 7.58 -4.28) = 0.405 MgO (-3.90 24.53 -3.11) = 1.616 CaO (5.39 -4.93 -3.80) = 0.497 Na2O (0.46 -3.96 0.00) = -0.278 K2O (-1.26 -1.49 1.51) = -0.279 P2O5 (-0.30 -0.30 -0.17) = -0.01The Parameters of the Weighted [to 1/standard_deviation] Matrix# Oxide ml-d m2-d m3-d p-d Standard deviation1 SiO2 (-9.176 -8.293 23.130) = -1.541 0.5782 TiO2 (-82.022 -86.517 -64.045) = 13.483 0.0093 Al203 (70.373 -87.078 -17.346) = -1.397 0.1724 FeO (-79.412 117.337 -66.254) = 6.192 0.0655 MgO (-69.519 437.255 -55.437) = 28.699 0.0566 CaO (78.116 -71.449 -55.072) = 7.101 0.0697 Na2O (12.466 -107.317 0.000) = -7.317 0.0378 K2O (-103.279 -122.131 123.770) = -22.131 0.0129 P2O5 (-103.448 -103.448 -58.621) = -3.448 0.003121The x2 Solution.The Parameters of the Matrix U from SVDCMP(^0.020 0.052 -0.122^)(^0.133 0.388 0.457^)(^0.175 -0.286 0.144^)( -0.255 0.289 0.311^)( -0.860 0.110 -0.000^)(^0.138 -0.339 0.345^)(^0.207 -0.010 0.079^)(^0.243 0.557 -0.579^)(^0.162 0.494 0.446^)The Parameters of the Diagonal Matrix W518.23968506 219.75834656 173.74630737The Parameters of the Matrix V(^0.102 -0.993 -0.054 )( -0.986 -0.094 -0.134 )(^0.128 0.067 -0.989 )122The Minimized x2 Solution (out of the total 1.0 gram)mineral_phase weightplag^0.0137051 gramsOpx 0.0455283 gramsXeno^-0.1302265 gramsWeighted Residual on Each Oxide Equation = RHS - LHS# oxide residual1 Si02 1.9746202 TiO2 10.2058933 Al203 -0.6557984 FeO -6.6898865 MgO 2.5247086 CaO 2.1119357 Na20 -2.6019578 K20 0.9629239 P205 -4.954653Chi_Squared Value = 196.324387Residual on Each Oxide Equation = RHS - LHS (using unweighted compositions)# Oxide residual1 Si02 1.1405412 TiO2 0.0908323 Al203 -0.1126664 FeO -0.4321675 MgO 0.1416366 CaO 0.1457247 Na2O -0.0960128 K20 0.0117489 P205 -0.014368The R2 Associated with the x2 Solution= 1.559404Del-Chi-Squared Mapped against Amounts of PhasesPlag vs Opx (Xenolith = -0.1302 gramPlag \ Opx 0.026 0.031 0.036 0.041 0.046 0.051 0.056 0.061 0.066-0.006 107.33 65.87 37.54 22.32 20.21 31.22 55.36 92.60 142.97-0.001 102.94 60.37 30.92 14.59 11.37 21.27 44.29 80.42 129.670.004 101.08 57.40 26.83 9.38 5.05 13.84 35.74 70.76 118.900.009 101.74 56.95 25.27 6.71 1.26 8.94 29.72 63.63 110.650.014 104.93 59.03 26.23 6.56 0.00 6.56 26.23 59.03 104.930.019 110.65 63.63 29.72 8.94 1.26 6.71 25.27 56.95 101.740.024 118.90 70.76 35.74 13.84 5.05 9.38 26.83 57.40 101.080.029 129.67 80.42 44.29 21.27 11.37 14.59 30.92 60.37 102.940.034 142.97 92.60 55.36 31.22 20.21 22.31 37.54 65.87 107.33Plag vs Xenolith (Opx = 0.0455 gram )Plag Xeno -0.150 -0.145 -0.140 -0.135 -0.130 -0.125 -0.120 -0.115 -0.110-0.006 35.42 29.05 24.40 21.45 20.21 20.68 22.86 26.76 32.36-0.001 26.19 19.92 15.36 12.51 11.37 11.94 14.21 18.20 23.900.004 19.49 13.32 8.85 6.10 5.05 5.72 8.09 12.17 17.960.009 15.32 9.24 4.87 2.21 1.26 2.02 4.49 8.67 14.560.014 13.68 7.69 3.42 0.85 0.00 0.85 3.42 7.69 13.680.019 14.56 8.67 4.49 2.02 1.26 2.21 4.87 9.24 15.320.024 17.96 12.17 8.09 5.72 5.05 6.10 8.85 13.32 19.490.029 23.90 18.20 14.21 11.94 11.37 12.51 15.36 19.92 26.190.034 32.36 26.76 22.86 20.68 20.21 21.45 24.40 29.05 35.42Opx vs Xenolith (Plag = 0.0137 gram )Opx \ Xeno -0.150 -0.145 -0.140 -0.135 -0.130 -0.125 -0.120 -0.115 -0.1100.026 94.39 94.46 96.24 99.73 104.93 111.84 120.46 130.79 142.830.031 54.54 53.10 53.36 55.34 59.03 64.42 71.53 80.34 90.860.036 27.80 24.84 23.60 24.06 26.23 30.12 35.71 43.01 52.020.041 14.18 9.71 6.95 5.90 6.56 8.93 13.00 18.79 26.290.046 13.68 7.69 3.42 0.85 0.00 0.85 3.42 7.69 13.680.051 26.29 18.79 13.00 8.93 6.56 5.90 6.95 9.71 14.180.056 52.02 43.01 35.71 30.12 26.23 24.06 23.60 24.84 27.800.061 90.86 80.34 71.53 64.42 59.03 55.34 53.36 53.10 54.540.066 142.83 130.79 120.46 111.84 104.93 99.73 96.24 94.46 94.39****************************************************************124The R2 Solution.The Paramenters of the Matrix U from SVDCMP(^0.115 0.304 -0.822^)(^0.013 -0.181 -0.002^)(^0.546 0.405 0.405^)( -0.281 -0.606 0.095^)( -0.748 0.551 0.204^)(^0.188 0.051 0.308^)(^0.118 -0.155 -0.015^)(^0.033 -0.096 -0.120^)(^0.005 -0.068 -0.004^)The Parameters of the Diagonal Matrix W32.80237579 6.23699284 17.70904922The Parameters of the Matrix V(^0.346 0.758 0.552 )( -0.934 0.332 0.131 )(^0.084 0.561 -0.823 )The Minimized R2 Solution (out of the total 1.0 gram)Mineral_phase weightPlag^0.0636599 gramsOpx 0.0699627 gramsXeno^-0.0276482 gramsResidual on Each Oxide Equation = RHS - LHS# Oxide residual1 Si02 0.1518952 TiO2 0.2045844 FeO 0.0779245 MgO 0.0561036 CaO 0.3867277 Na20 -0.0222318 K20 -0.0437959 P205 0.025387The Least Residual_Squared = 0.228821125Del-Residual-Squared Mapped against Amounts of PhasesPlag vs Opx (Xenolith = -0.0276 grams )Plag\Opx 0.050 0.055 0.060 0.065 0.070 0.075 0.080 0.085 0.0900.044 0.226 0.123 0.067 0.059 0.099 0.186 0.320 0.502 0.7310.049 0.246 0.127 0.056 0.032 0.056 0.127 0.245 0.411 0.6250.054 0.278 0.143 0.056 0.017 0.025 0.080 0.183 0.333 0.530.059 0.323 0.172 0.069 0.014 0.006 0.046 0.133 0.267 0.4490.064 0.380 0.213 0.095 0.024 0.000 0.024 0.095 0.213 0.3800.069 0.449 0.267 0.133 0.046 0.006 0.014 0.069 0.172 0.3230.074 0.531 0.333 0.183 0.080 0.025 0.017 0.056 0.143 0.2780.079 0.625 0.411 0.245 0.127 0.056 0.032 0.056 0.127 0.2460.084 0.731 0.502 0.320 0.186 0.099 0.059 0.067 0.123 0.226Plag vs Xenolith (Opx = 0.0700 grams )Plag\Xeno -0.048 -0.043 -0.038 -0.033 -0.028 -0.023 -0.018 -0.013 -0.0080.044 0.116 0.094 0.084 0.086 0.099 0.124 0.160 0.208 0.2680.049 0.092 0.065 0.050 0.047 0.056 0.076 0.107 0.151 0.2050.054 0.080 0.049 0.029 0.021 0.025 0.040 0.067 0.105 0.1560.059 0.080 0.044 0.020 0.007 0.006 0.017 0.039 0.073 0.1180.064 0.093 0.052 0.023 0.006 0.000 0.006 0.023 0.052 0.0930.069 0.118 0.073 0.039 0.017 0.006 0.007 0.020 0.044 0.0800.074 0.156 0.105 0.067 0.040 0.025 0.021 0.029 0.049 0.0800.079 0.205 0.151 0.107 0.076 0.056 0.047 0.050 0.065 0.0920.084 0.268 0.208 0.160 0.124 0.099 0.086 0.084 0.094 0.116Opx vs Xenolith (Plag = 0.0637 grams )Opx\Xeno -0.048 -0.043 -0.038 -0.033 -0.028 -0.023 -0.018 -0.013 -0.0080.050 0.384 0.365 0.358 0.363 0.380 0.408 0.447 0.499 0.5610.055 0.240 0.216 0.203 0.203 0.213 0.236 0.270 0.316 0.3730.060 0.143 0.114 0.096 0.090 0.095 0.112 0.140 0.181 0.2320.065 0.094 0.059 0.036 0.024 0.024 0.035 0.058 0.093 0.1390.070 0.093 0.052 0.023 0.006 0.000 0.006 0.023 0.052 0.0930.075 0.139 0.093 0.058 0.035 0.024 0.024 0.036 0.059 0.0940.080 0.232 0.181 0.140 0.112 0.095 0.090 0.096 0.114 0.1430.085 0.373 0.316 0.270 0.236 0.213 0.203 0.203 0.216 0.2400.090 0.561 0.499 0.447 0.408 0.380 0.363 0.358 0.365 0.384126**************************************************************************Del Chi_Sq Value = 10.644501 less than which the Probability of the Chi_Sq Variable Is 0.90000**************************************************************************The Orientations and Lengths of Semi-Axes of the x2 Solution Confidence Region(0.102, -0.986, 0.128) semi-axis length = 0.0063(-0.993, -0.094, 0.067) semi-axis length = 0.0148(-0.054, -0.134, -0.989) semi-axis length = 0.0188127/*^ PROGRAM MASSBAL.0/* BY: PEIWEN KE/* May, 1992/* This program is written to calculate the amount of minerals fractionating from magma A to magma B.The principles involved are material balance, and the Chi-squared minimization fitting technique. */#include <stdio.h>#include <math.h>#include "ganunln.c"#include "gammq.c"#include "gcf.c"#include "gser.c"#include "malloc.h"#include "nrutil.c"#include "nrutil.h"#define TOL 1.0e-5/* THE MAIN PROGRAM MASSBAL.0^ */main(argc, argv)int argc;char **argv;{int i,j, k, il,j1, kl, I;int number_of oxide;int number_of mineral_phase;int degrees_of freedom;int number_of characters_in_title;int string_length();char title[81];char name_of component[10];char name_of oxide[12][12];char name_of_parent[20];char name_of daughter[20];char name_of mineral_phase[20][8];char name_of error[20];float step, alpha, r_2, df;float wmax, thresh, sum, sum_chi_to_r, chisq, residual_square, residual_chi_to_r;float mole_pt_of_p[12];float mole_pt_of d[12];float mole_pt_of m[12][8];float relative_error[12];float std_deviation[12];float diff of_p_and_d[12];float diff of m_and_d[12][8];float diff of m_and_d_save[12][8];float diff of_p_and_d_std[12];float diff of m_and_d_std[12][8];float diff of m_and_d_std_save[12][8];128float weight_of mineral_phase[8];float weight_of mineral_phase_r[8];float residual_std[12];float residual_std_chi_to_r[12];float residual_412];float w[8];float v[8][8];float weight_of mineral[8][16];float chisq_array[16][16][8];float chisq_diff array[16][16][8];float one_phase_chisq_diff array[16];float one_phase_residual_square_diff array[16];float w_r[8];float v_r[8][8];float weight_of mineral_r[8][16];float residual_square_array[16][16][8];float residual_square_diff array[16][16][8];float p;float one_minus_p;float half degrees_of freedom;float left_ side_ of_ root;.float nght side_of root;float fittedidelta_chisq;float left_side_of fitted_delta_chisq;float right_side_of fitted_delta_chisq;float motification_of theleft_and_right_sides;float delta_weight_of mineral_phase[8];float delta;float unit_length_of axis[8];float delta_chisq_length_of axis[8];FILE *fp;FILE *fo;FILE *fo2;fp = fopen(argv[1], "r");fo = fopen(argv[2], "w");fo2 = fopen("plot.out", "w");/* inputting data file *//* outputting data file *//* outputting data file for plotting only for threemineral phases*/if (argc != 3){printf("please check the opening data file\n");exit(0);else if ((fp = fopen(argv[1], "r")) = = NULL)printf("Can't open file %s\n", argv[1]);exit(0);}/* initiate the parameters used in calling the functions */for (i=0; i< =numberof mineral_phase; i+ +)){129{w[i]=0.00;for (j =0; j < =number— of mineral_phase; j + +)v[i][j]=0.00;}for (i = 0; i< = number_of mineral_phase; i = i +1)weight_of mineral_phase[i] = 0.0;/* input the raw data */fscanf (fp, "%[^ \n]", title);fscanf (fp, " %An", &step);fscanf (fp, "%d%d\n“, &number_of oxide, &number_of mineral_phase);fscanf (fp, " %s %s %s ", name_of component, name_of_parent, name_of_daughter);for (i=1; i< =number_of mineral_phase; i+ +)fscanf (fp, "%s", name_of mineral_phase[i]);fscanf (fp, "%s\n", name of error);number_of characters_in_title = string_length (title);degrees_of freedom = number_of_oxide - number_of_mineral_phase;half degrees_of freedom = degrees_of freedom/2.0;for (i =1; i< = 2 +number_of characters_intitle; i+ +)fprintf (fo, "*");fprintf (fo, "\n");fprintf (fo, " %s \n", title);for (i =1; i< =2 +number_of characters_in_title; i+ +)fprintf (fo, "*");fprintf (fo, "\n\n\n");fprintf (fo, "INPUT PARAMETRS\n");fprintf (fo, "^\n");fprintf (fo, "number_of oxide = %d, number_of_mineral_phase = %d,degrees_of freedom = %d\n", number_of_oxide, number_of_mineral_phase,degrees_of freedom);fprintf (fo, " % 10s % 10s % 10s", name_of component, name_of_parent, name_of_daughter);for (i=1; i< =number_of mineral_phase; i+ +)fprintf (fo, " % 10s", name_of mineral_phase[i]);fprintf (fo, " %18s\n",name_of error);for (i=1; i< =number_of oxide; i+ +){fscanf (fp, "%s%f%f", name_of oxide[i], &mole_pt_of_p[i], &mole_pt_of d[i]);fprintf (fo, " % 10s %10.2f%10.2f", name_of oxide[i], mole_pt_of_p[i], mole_pt_of d[i]);diff of_p_and_d[i] = mole_pt_of_p[i] - mole_pt_of d[i];for (j =1; j < =number_of mineral_phase; j+ +){fscanf (fp, "%f", &mole_pt_of m[i][j]);fprintf (fo, "%10.2f ", mole_pt_of m[i][j]);diff of m_and_d[i][j] = mole_pt_of m[i][j] - mole_pt_of d[i];diff of mand_dsave[i][j] = diff of mandd[i][j];}fscanf (fp, "%f", &relative_error[i]);130fprintf (fo, "%15.2f", relativeerror[i]);std_deviation[i] = 1.0;if (mole_pt_of_p[i] != 0.00)if (relative_error[i] != 0.00)std_deviation[i] = mole_pt_of_p[i]*relative_error[i]/100.0;diff of_p_and_d_std[i] = diff of_p_and_d[i]/std_deviation[i];for (j = 1; j < =number_of mineral_phase; j + +){diff of m_and_d_std[i][j] = diff of m_and_d[i][Wstd_deviation[i];diff of m_and_d_std_save[i][j]=diff of m_and_d_std[i][j];}fprintf (fo, " \n");}fscanf (fp, "%f", &p);one_minus_p = 1.0 - p;/* output the parameters of the compositional difference matrix */fprintf (fo, "\nThe Parameters of the Compositional Difference Matrix\n");fprintf (fo, "^ \n");fprintf (fo, "# oxide^");for (i= 1; i< =number_of mineral_phase; i=i +1)fprintf(fo, "m%i-d ", i);fprintf (fo, " p-d \n");for (i= 1; i< =number_of oxide; i=i +1){fprintf (fo, "%i %6s ( ", i, name_of oxide[i]);for (j =1; j < =number_of mineral_phase; j + +)fprintf (fo, " %5.2f ", diff of m_and_d[i][j]);fprintf (fo, ") = %5.2f\n", diff of_p_and_d[i]);}/* output the parameters of the weighted [to 1/std_deviation] matrix A */fprintf (fo, "\nThe Parameters of the weighted to [1/standard_deviation] Matrix\n");fprintf (fo, "^ \n");fprintf (fo, "# oxide^");for (i=1; i< =number_of mineral_phase; i+ +)fprintf(fo, "m%i-d^", i);fprintf (fo, "^p-d^standard deviation\n");for (i= 1; i < =number_of oxide; i + +){fprintf (fo, "%i %6s ( ", i, name_of oxide[i]);for (j = 1; j < =number_of mineral_phase; j + +)fprintf (fo, "%7.3f ", diff of m_and_d_std[i][j]);^fprintf (fo, ") = %7.3f^%7.3f\n", diff of_p_and_d_std[i], std_deviation[i]);}svdcmp(diff of m_and_d_std, number_of oxide, number_of mineral_phase, w, v);wmax =0.0;131for (j =1; j < = number_of mineral_phase; j+ +)if (w[j] > wmax ) wmax=w[j];thresh = TOL *wtnax;for (j =1; j < =number of mineral_phase; j + +)if (w[j] < thresh) w[j]=0.0;fprintf (fo, "\n\nThe Chi-Squared Solution. \n");fprintf (fo, "^ \n");1* output the parameters of the U matrix from svdcmp *1fprintf (fo, "\nThe Parameters of the Matrix U from SVDCMP\n");fprintf (fo, "^ \n");for (i=1; i< =number_of oxide; i= i + 1){fprintf (fo, "( ");for (j =1; j < =number of mineral_phase; j + +)fprintf (fo, " %7.3f ", diff of m_and_d_std[i][j]);fprintf (fo, ")\n");}1* output the parameters of the diagonal matrix W *1fprintf (fo, "\nThe Parameters of the Diagonal Matrix W\n");fprintf (fo, "^ \n");for (i=1; i< =numberof mineral_phase; i=i +1)fprintf (fo, " % 12.8f ", w[i]);fprintf (fo, "\n");1* output the parameters of the matrix V *1fprintf (fo, "\nThe Parameters of the Matrix V \n");fprintf (fo, "^\n");for (i =1; i< =numberof mineral_phase; i = i +1){fprintf (fo, "( ");for (j =1; j < =numberof mineral_phase; j=j + 1)fprintf (fo, " %7.3f ", v[i][j]);fprintf (fo, " )\n");}svbksb(diff of m_and_dstd, w, v, number_of oxide, number_of mineral_phase,diff of_p_and_d_std, weight_of mineral_phase);chisq = 0.0;residual_chi_tor = 0.0;for (i=1; i< =number_of oxide; i + +){sum = 0.0;sumchito_r = 0.0;for (j =1; j < =number_of mineral_phase; j + +){sum + = diff of m_and_d_std_save[i][j]*weight_of mineral_phase[j];sum_chi_to_r += diff of m_and_d_save[i][j]*weight_of mineral_phase[j];}residual_std[i] = diff of_p_and_d_std[i] - sum;132residual_std_chi_to_r[i] = diff of_p_and_d[i] - sum_chi_to_r;chisq + = residual_std[i]*residual_std[i];residual_chi_to_r + = residual_std_chi_to_r[i]*residual_std_chi_to_r[i];}/* output of the solution of least chi-square */fprintf (fo, " \nThe Minimized Chi-Squared Solution (out of the total 1.0 gram)\n");fprintf (fo, "^ \n");fprintf (fo, "^mineral_phase^weight \n");for (i= 1; i< =number_of mineral_phase; i=i + 1)fprintf(fo, " %20s % 10.7f grams \n", name_of mineral_phase[i],weight_of mineral_phase[i]);fprintf (fo, " \n");/* output of residual on each oxide equation */fprintf (fo, "Weighted Residual on Each Oxide Equation = RHS - LHS\n");fprintf (fo, "^ \n");fprintf (fo, "# oxide^residual \n");for (i= 1; i< =number_of oxide; i+ +)fprintf (fo, "%d %6s^%f \n", i, name_of oxide[i], residual_std[i]);fprintf (fo, " \n");/* output the least chi_square */fprintf (fo, "Chi_Squared Value = %fln", chisq );fprintf (fo, "Residual on Each Oxide Equation = RHS - LHS (using unweightedcompositions)\n");fprintf (fo, "^ \n");fprintf (fo, "# oxide^residual \n");for (i= 1; i< =number_of oxide; i+ +)fprintf (fo, "%d %6s^%f \n", i, name_of oxide[i], residual_std_chi_to_r[i]);fprintf (fo, "\n");fprintf (fo, "The Residual-Squared Associated with the Chi_Squared solution= %fin",residual_chi_to_r );for (i1=1; il < =number_of mineral_phase; il + +)for (i=1; i< =9; i++)weight_of mineral[il][i] = weight_of mineral_phase[il] + step*(i-5);/* output the distribution of the chi_square difference between the best solution and the chi_squarearound it */fprintf (fo, " \nDel-Chi-Sqaured Mapped against Amounts of Phases\n");fprintf (fo, "^ \n");1=0;for (i= 1; i< =number_of mineral_phase; i + +)if (number_of mineral_phase= =1){fprintf (fo, " % lOs ", name_of mineral_phase[i]);for (j=1; j<=9; j++)133fprintf (fo, " %6.3f ", weight_of mineral[1][j]);fprintf (fo, "\n");fprintf (fo, "delta chisq ");for (j=1; j<=9; j++){one_phase_chisq_diff array[j] = 0.0;{for (k =1; k< =number_of oxide; k+ +)onephase_chisq_diff array[j] + = ((diff ofp_and_d_std[k] -diff of m_and_d_std_save[k][1]*weight_of mineral[1][j])*(diff of_p_and_d_std[k] -diff of—m_and_d_std_save[k][1]*weight_of mineral[1][j]) - chisq);}fprintf (fo, " %6.3f ", one_phase_chisq_diff array[j]);}fprintf (fo, "\n");}else{for (j =i + 1; j < =numberof mineral_phase; j + +){I+ =1;fprintf (fo, " \n %s vs %s (", name_of mineral_phase[i], name_of mineral_phase[j]);for (k =1; k< =numberof mineral_phase; k+ +){if (k!=i)if (k! =j)fprintf (fo, "%s = %7.4f gram ", name of mineral_phase[k], weight_of mineralphase[k]);}fprintf (fo, ")\n");fprintf (fo, " %s\\ %s ", name_of mineral_phase[i], name_of mineral_phase[j]);for (i1=1; il< =9; il++)fprintf (fo, "%6.3f ", weight_of mineral[j][il]);fprintf (fo, "\n");for (il =1; il < =9; il+ +)fprintf (fo, " %6.3f ", weight_of mineral[i][il]);for (j1=1; jl<=9; j1++){chisq_array[il][jl][I]=0.0;for (k =1; k< =numberof oxide; k+ +) {sum = 0.0;for (kl =1; kl < =numberof mineral_phase; kl + +){if (kl= =i)sum += diff of m_and_d_std_save[k][kl]*weightof mineral[kl][il];else {if (kl= =j)sum += diff of m_and_d_std_save[k][k1]*weight_of mineral[kl][j1];elsesum + = diff of mand_d_stdsave[k][k1]*weight_of mineralphase[kl];}}residual_std[k] = diff of_p_and_d_std[k] - sum;chisq_array[i1][j1][I] + = residual_std[k]*residual_std[k];}chisq_diff array[i1] [j 1] [I] = chisq_array[il][jl][I]-chisq;fprintf (fo, " %6. 2f ", chisq_diff array[i 1 ][j MID;}134fprintf (fo, " \n");}}}/**//* output the calculation results based on the residual_square */svdcmp(diff of m_and_d, number_of oxide, number_of mineral_phase, w_r, v_r);wmax =0.0;for (j =1; j < =number_of mineral_phase; j+ +)if (w_r[j] > wmax ) wmax =w_r[j];thresh = TOL *wmax;for (j =1; j < =number_of mineral_phase; j + +)if (w_r[j] < thresh) w_r[j]=0.0;fprintf (fo,nn****************************************************************\ n\n ”) ;fprintf (fo, "The Residual_Squared Solution\n");fprintf (fo, "^ \n");/* output the parameters of the U matrix from svdcmp from r-square */fprintf (fo, " \nThe Paramenters of the Matrix U from SVDCMP\n");fpiintf (fo, "^ \n");for (i = 1; i< =number_of oxide; i=i + 1){fprintf (fo, "( ");for (j =1; j < =number_of tnineral_phase; j + +)fprintf (fo, " %7.3f ", diff of m_and_d[i][j]);fprintf (fo, ")\n");}/* output the parameters of the diagonal matrix W based on r_square */fprintf (fo, " \nThe Parameters of the Diagonal Matrix W\n");fprintf (fo, " ^\n");for (i=1; i< =number_of mineral_phase; i=i +1)fprintf (fo, " %12.8f ", w_r[i]);fprintf (fo, " \n");/* output the parameters of the matrix V based on r_square */fprintf (fo, " \nThe Parameters of the Matrix V \n");fprintf (fo, " ^\n");for (i = 1; i < =number_of mineral_phase; i=i + 1){fprintf (fo, "( ");for (j = 1; j < =number_of mineral_phase; j=j +1)fprintf (fo, " %7.3f ", v_r[i][j]);fprintf (fo, " )\n");135}svbksb(diff of m_and_d, w_r, v_r, number_of oxide, number_of mineral_phase,diff of_p_and_d, weight_of mineral_phase_r);residual_square = 0.0;for (i= 1; i< =number_of oxide; i+ +){sum = 0.0;for (j = 1; j < =number_of mineral_phase; j + +)sum + = diff of m_and_d_save[i][j]*weight_of mineral_phase_r[j];residual_r[i] = diff of_p_and_d[i] - sum;residual_square + = residual_r[i]*residual_r[i];}/* output of the solution of least residual-square */fprintf (fo, "\nThe Minimized Residual-Squared Solution (out of the total 1.0 gram)\n");fprintf (fo, "^ \n");fprintf (fo, "^Mineral_phase^weight \n");for (i = 1; i< =number_of mineral_phase; i= i + 1)fprintf(fo, " %20s^%10.7f grams \n", name_of mineral_phase[i],weight_of mineral_phase_r[i]);fprintf (fo, " \n");/* output of residual on each oxide equation */fprintf (fo, "Residual on Each Oxide Equation = RHS - LHS\n");fprintf (fo, "^ \n");fprintf (fo, "# oxide^residual \n");for (i= 1; i< =number_of oxide; i+ +)fprintf (fo, "%d %6s^%f \n", i, name_of oxide[i], residual_r[i]);fprintf (fo, " \n");/* output the least residual_square */fprintf (fo, "The Least Residual_Squares = %An", residual_square );for (i1=1; il < =number_of mineral_phase; il + +)for (i=1; i< =9; i++)weight_of mineral_r[il][i] = weight_of mineral_phase_r[il] + step*(i-5);/* output the distribution of the residual_square difference between the best solution and theresidual_square around it */fprintf (fo, "\nDel-Residual-Squared Mapped against Amounts of Phases\n");fprintf (fo, "^ \n");I=0;for (i = 1; i< =number_of mineral_phase; i+ +)if (number_of mineral_phase= =1){fprintf (fo, " % 10s ", name_of mineral_phase[i]);136for(j=1; j<=9; j++)fprintf (fo, " %6.3f ", weightof mineral_r[l][j]);fprintf (fo, "\n");fprintf (fo, "delta-resisq ");for(j=1; j<=9; j++){one_phase_residual_square_diff array[j] = 0.0;{for (k= 1; k< =number_of oxide; k+ +)one_phase_residual_square_diff array[j] + = ((diff of_p_and_d[k] -diff of m_and_d_save[k][1]*weight_of mineral_r[l][j])*(diff of_p_and_d[k] -diff of m_and_d_save[k][1]*weight_of mineral_r[l][j]) - residual_square);fprintf (fo, " %6.3f ", one_phase_residual_square_diff array[j]);}fprintf (fo, "\n");}else{for (j=i +1; j< =number_of mineral_phase; j + +){I+ =1;fprintf (fo, " \n %s vs %s (" , name_of mineral_phase[i], name_of mineral_phase[j]);for (k=1; k< =numberof mineral_phase; k+ +){if (k!=i)if (k!=j)fprintf (fo, "%s = %7.4f grams ", name_of mineral_phase[k],weight_of mineral_phase_r[k]);}fprintf (fo, ")\n");fprintf (fo, " %s\\ %s ",name_of mineral_phase[i],name_of mineral_phase[j]);for (i1=1; il< =9; il++)fprintf (fo, "% 6.3f ", weight_of mineral_r[j][i1]);fprintf (fo, "\n");for (i1=1; it < =9; il++) {fprintf (fo, " %6.3f ", weight_of mineral_r[i][il]);for (j1=1; jl<=9; j1++){residual_square_array[i1][j1][I] = 0.0;for (k=1; k< =numberof oxide; k+ +){sum = 0.0;for (k1=1; kl < =numberof mineral_phase; kl + +){if (kl= =i)sum += di ff of m_and_d_save[k][kl]*weight_of mineral_r[kl][i1];else {if (kl= =j)sum += di ff of m_and_d_save[k][kl]*weight_of mineral_r[kl] [j 1];elsesum + = diff of m_and_d_save[k][k1]*weight_of mineral_phase_r[kl];}}residual_r[k] = diff of_p_and_d[k] - sum;residual_square_array[i1][j1][I] + = residual_r[k]*residual_r[k];}residual_square_diff array[i1][j 1] [I] = residual_square_array[i1][j 1] [I]-residual_square;fprintf (fo, " %6.3f ", residual_square_diff array[il][jl][I]);137}fprintf (fo, " \n");}}}/* 4.//* calculation of the confidence level of the solutions around the best solution *//* call gammq(a,x) function to calculate the delta_chisq on the given significant level p */if (degrees_of freedom > 0){i = 1;motification_of theleft_and_right_sides = 1.0;while (motification_of the_left_and_right_sides > 0){left_side_of fitted_delta_chisq = i*0.001/2.0;left_side_of root = gammq(half degrees_of freedom, left_side_of fitted_delta_chisq) -one_minus_p;right_side_of fitted_delta_chisq = (i + 1)*O. 001 /2.0;right_side_of root = gammq(half degrees_of freedom, right_side_of fitted_delta_chisq) -one_minus_p;motification_of the_left_and_right_sides = left_side_of root * right_side_of root;i+ =1;}printf ("motification_of the_left_and_right_sides = %20.18fin",motification_of theleft_and_right_sides);printf ("left_side_of fitted_delta_chisq = %f \n", left_side_of fitted_delta_chisq);printf ("right_side_of fitted_delta_chisq = %f \n", right_side_of fitted_delta_chisq);fitted_delta_chisq = (left_side_of fitted_delta_chisq + right_side_of fitted_delta_chisq);printf ("fitted_delta_chisq = %f \n", fitted_delta_chisq);alpha = sqrt(fitted_delta_chisq);printf ("alpha = square root of fitted delta chisq = %nn", alpha);fprintf (fo, " \n\n*****************************************************\ n ") ;fprintf (fo, "Del_Chisq Value = %f Less than Which the Probability of Chi-Sq Variable Is %f\n",^fitted_delta_chisq, p);fprintf (f0, "*****************************************************\ n\n\e) ;/* outputting the parameters of the solution confidence region */fprintf (fo, " \n\nThe Orientations and Lengths of Semi-Axes of Chi-Squared Solution\n");fprintf (fo,"^ \n\n");for (i = 1; i< =number_of mineral_phase;i + +){fprintf (fo, "(");for (j = 1; j <number_of mineral_phase; j + +)fprintf (fo, " %5.3f, ", NUN);fprintf (fo, " %5.3f)^", v[number_of mineral_phase][i]);unit_length_of axis[i] = 1.0/w[i];delta_chisq_length_of axis[i] = alpha*unit_length_of axis[i];fprintf (fo, "Del-Chisq Length = %7.4f\n", delta_chisq_length_of axis[i]);138/* output the chi and residual square distribution for plotting profile into the data file plot.out *//* output the chi_square distribution around the best solution */if (numberof mineral_phase= =3){fprintf (fo2, "\nChi_squared mapped against amounts of phases\n");fprintf (fo2, "^ \n");I=0;for (i=1; i< =numberof mineral_phase; i+ +)for (j=i + 1; j < =number_of mineral_phase; j + +){I+ =1;fprintf (fo2, " \n%s vs %s (", name_of mineral_phase[i], name_of_mineral_phase[j]);for (k=1; k< =numberof mineral_phase; k+ +){if (k! =i)if (k!=j)fprintf (fo2, "%s = %7.4f gram ", name_of mineral_phase[k], weight_of mineral_phase[k]);}fprintf (fo2, ")\n");fprintf (fo2, " %s\\ %s ", name_of mineral_phase[i], name_of_mineral_phase[j]);for (i1=1; it < =9; il+ +)fprintf (fo2, " %6.3f ", weight_of mineral[j][il]);fprintf (fo2, "\n");for (il =1; il < =9; il+ +) {fprintf (fo2, " %6.3f^", weightof mineral[i][il]);for (j1=1; jl<=9; jl++)fprintf (fo2, " %6.2f ", chisq_array[il][jl][I]);fprintf (fo2, "\n");}}/* output the distribution of the chi_square difference between the best solution and the chi_squarearound it */fprintf (fo2, "\nDel-chi-squared against amounts of phases\n");fprintf (fo2, "^ \n");I=0;for (i =1; i< =number of mineral_phase; i+ +)for (j=i +1; j < =numberof mineral_phase; j + +){I+=1;fprintf (fo2, " \n%s vs %s (", nameof mineral_phase[i], name_of_mineral_phase[j]);for (k=1; k< =numberof mineral_phase; k+ +){if (k!=i)if (k! =j)fprintf (fo2, " %s = %7.4f gram ", name_of mineral_phase[k], weight_of mineral_phase[k]);}fprintf (fo2, ")\n");fprintf (fo2, " %s\\ %s ", nameof mineral_phase[i], name_of_mineral_phase[j]);for (i1=1; it < =9; il++)fprintf (fo2, " %6.3f ", weight_of mineral[j][il]);139fprintf (fo2, "\n");for (i1=1; il < =9; il+ +) {fprintf (fo2, " %6.3f^", weight_of mineral[i][il]);for (j1=1; jl<=9; jl++)fprintf (fo2, " %6.2f ", chisq_diff array[il][j1][11);fprintf (fo2, "\n");}}/* output the residual_square distribution around the best solution */fprintf (fo2, "\nR-squared mapped against amounts of phases\n");fprintf (fo2, "^ \n");I=0;for (i=1; i< =numberof mineral_phase; i+ +)for (j = i + 1; j < =numberof mineral_phase; j + +){I+ =1;fprintf (fo2, " \n%s vs %s (", name_of mineral_phase[i], name_of mineral_phase[j]);for (k =1; k< =numberof mineral_phase; k+ +){if (k!=i)if (k!=j)fprintf (fo2, "%s = %7.4f grams ", name_of mineral_phase[k],weight_of mineral_phase_r[k]);)fprintf (fo2, ")\n");fprintf (fo2, " %s\\ %s ",name_of mineral_phase[i],name_of mineral_phase[j]);for (i1=1; it < =9; il+ +)fprintf (fo2, " %6.3f ", weight_of mineral_r[j][i1]);fprintf (fo2, "\n");for (i1=1; it < =9; il+ +) {fprintf (fo2, " %6.3f ", weight_of mineral_r[i][il]);for (j1=1; jl<=9; jl++)fprintf (fo2, " % 6 . 3 f " , residual_square_array[il][jl][I]);fprintf (fo2, "\n");}}/* output the distribution of the residual_square difference between the best solution and theresidual_square around it */fprintf (fo2, "\nDelta-R-squared mapped against amounts of phases\n");fprintf (fo2, "^ \n");I=0;for (i =1; i < =number of mineral_phase; i+ +)for (j=i +1; j < =numberof mineral_phase; j + + )(I+ =1;fprintf (fo2, " \n %s vs %s (", name_of mineral_phase[i], name_of mineral_phase[j]);for (k =1; k< =numberof mineral_phase; k+ +){if (k!=i)if (k!=j)fprintf (fo2, "%s = %7.4f grams ", name_of mineral_phase[k],weight_of mineral_phase_r[k]);}fprintf (fo2, ")\n");140fprintf (fo2, " %s\ \ %s ",name_of mineral_phase[i],name_of mineral_phase[j]);for (i1=1; il< =9; il++)fprintf (fo2, " %6.3f ", weight_of mineral_r[j][il]);fprintf (fo2, "\n");for (i1=1; il< =9; il++) {fprintf (fo2, " %6.3f ", weight_of mineral_r[i][il]);for (j1=1; jl<=9;j1++)fprintf (fo2, " %6.3f ", residual_square_diff array[il][j1][I]);fprintf (fo2, "\n");}}}/* The SVDCMP.0 function called by main program.static float at,bt,ct;#define PYTHAG(a,b) ((at= fabs(a)) > (bt=fabs(b)) ? \(ct=bt/at,at*sqrt(1.0+ct*ct)) : (bt ? (ct=at/bt,bt*sqrt(1.0 +ct*ct)): 0.0))static float maxargl,maxarg2;#define MAX(a,b) (maxarg1=(a),maxarg2=(b),(maxargl) > (maxarg2) ?\(maxargl) : (maxarg2))#define SIGN(a,b) ((b) > 0.0 ? fabs(a) : -fabs(a))void svdcmp(a,m,n,w,v)float a[12][8], w[8], v[8][8];int m,n;{int flag,i,its,j,jj,k,l,nm;float c,f,h,s,x,y,z;float anorm=0.0,g=0.0,scale=0.0;float *rv1,*vector();void nrerror(),free_vector();if (m < n) nrerror("SVDCMP: You must augment A with extra zero rows");rvl=vector(1,n);for (i=1;i< =n;i+ +) {1=i+1;rvl[i]=scale*g;g=s=scale=0.0;if (i < = m)for (k=i;k < =m;k+ +) scale + = fabs(a[k][i]);if (scale) {for (k=i;k<=m;k+ +) {a[k][i] /= scale;s + = a[k][i]*a[k][i];f=a[i][i];g = -SIGN(sqrt(s),f);h=f*g-s;141a[i][i] = f-g;if (i != n) {for (j =1;j < =n;j + +) {for (s=0.0,k=i;k< =m;k+ +) s + =a[k][i]*a[k][j];f=s/h;for (k=i;k< = m;k + +) a[k][j] + = f*a[k][i];}}for (k=i;k< =m;k+ +) a[k][i] *= scale;}w[i] = scale*g;g=s=scale=0.0;if (i < = m && != n) {for (k=1;k< =n;k+ +) scale + = fabs(a[i][k]);if (scale) {for (k=1;k<=n;k+ +) {a[i][k] /= scale;s + = a[i][k]*a[i][k];}f=a[i][1];g = -SIGN(sqrt(s),0;h= f*g-s;a[i][1]=f-g;for (k=1;k< =n;k+ +) rvl[k]=a[i][k]/h;if (i != m) {for (j =1;j < =m;j+ +) {for (s=0.0,k=1;k< =n;k+ +) s +=for (k=1;k< =n;k+ +) a[j][k] + = s*rvl[k];}}for (k=1;k< =n;k+ +) a[i][k] *= scale;}}anorm=MAX(anorm,(fabs(w[i])+fabs(rvl[i])));}for (i=n;i > =1;i--) {if (i <n){if (g) {for (j=1;j< =n;j+ +)v[j][i] = (a[i][j]/a[i][1])/g;for (j =1;j < =n;j+ +) {for (s =0.0,k =1;k < =n;k+ +) s + = a[i][k]*v[k][j];for (k=1;k< =n;k+ +) v[k][j] + = s*v[k][i];}}for (j-----1;j < =n;j + +) v[i][j] =v[j][i] =0.0;}v[i][i] =1.0;g=--rvl[i];1=-i;}142for (i=n;i > =1;i--) {1=i+1;g=w[i];if (i < n)for (j =1;j < =n;j+ +) a[i][j]=0.0;if (g) {g= 1.0/g;if (i != n) {for (j =1;j < =n;j+ +) {for (s=0.0,k=1;k< = m;k+ +) s += a[k][i]*a[k][j];f=(s/a[i][i]rg;for (k=i;k< =m;k + +) a[k][j] + = f*a[k][i];}}for (j=i;j < =m;j+ +) a[j][i] *= g;} else {for (j=i;j < =m;j+ +) a[j][i]=0.0;}+ +a[i][i];}for (k=n;k > =1;k--) (for (its= 1;its < =30;its+ +) {flag= 1;for (1=k;1> =1;1--) {nm=1-1;if (fabs(rv1[1])+anorm = = anorm) {flag=0;break;}if (fabs(w[nm])+anorm == anorm) break;}if (flag) {c=0.0;s= 1.0;for (i=1;i<=k;i+ +) {f=s*rvl[i];if (fabs(f)+anorm != anorm) {g=w[i];h=PYTHAG(f,g);w[i]=h;h= 1.0/h;c=g*h;s=(-f*h);for (j=1;j<=m;j+ +) {Y=alllinin];z=a[j][i];a[j][nm]=y*c+z*s;a[j][i]=z*c-y*s;}}}}z=w[k];if (1 = = k) {143if (z < 0.0) {w[k] = -z;for (j=ld < =n;j + +) v[j][k]=(-v[j][k]);}break;}if (its = = 30) nrerror("No convergence in 30 SVDCMP iterations");x=w[1];nm=k-1;y =w[nm];g=rvl[nm];h=rvl[k];f=((y-z)*(y+z)+(g-h)*(g+h))/(2.0*h*y);g=PYTHAG(f,1.0);f=((x-z)*(x+z)+h*((y/(f+SIGN(g,f)))-h))/x;c=s=1.0;for (j =1;j < = mud + +)i=j+1;g=rvl[i];Y= wEil;h=s*g;g=c*g;z=PYTHAG(f,h);rvl[j]=z;c= f/z;s=h/z;f=x*c+g*s;g=g*c-x*s;h= y*s;y=y*c;for (jj=1;jj < =ndj+ +) {x=v[jj][j];z=v[jj][i];v[jj][j]= x*c + z*s;v[j][i]=z*c-x*s;}z=PYTHAG(f,h);w[j]=z;if (z) {z= 1.0/z;c=f*z;s=h*z;}f=(c*g)+(s*y);x=(c*y)-(s*g);for (jj=1;jj < =m;jj+ +) {Y=41.11U];z= ii[jj][i];a[jj][j]=y*c+z*s;a[jj][i]=z*c-y*s;}rv1[1] =0.0;rv l[k] = f;144w[k]=x;free_vector(rv1,1,n);}#undef SIGN#undef MAX#undef PYTHAG/* The SVBKSB.0 function called by main program.void svbksb(u,w,v,m,n,b,x)float u[12][8], w[8], v[8][8], b[S],x[8];int m,n;{int jj,j,i;float s,*tmp,*vector();void free vector();tmp = vector(1,n);for (j= 1;j < =n;j + +) {s =0.0;if (w{l]) {for (i = 1;i < =m;i+ +) s += u[i][j]*b[i];s /= with}tmp[j]=s;}for (j =1;j < =n;j + +) {s=0.0;for (jj =1; jj < =n;jj+ +) s + = v[j][jj]*tmp[jj];x[j]=s;}free_vector(tmp,l,n);}/* Function to count the number of characters in a string (using pointer).int string_length (string)char *string;{char *cptr = string;while ( *cptr )+ +cptr;return (cptr - string);}/* The Gammln.c function called by main program.^ *1float gammln(xx)float xx;{double x,tmp,ser;static double cof[6]= {76.18009173,-86.50532033,24.01409822,-1.231739516,0.120858003e-2,-0.536382e-5};int j;x =xx-1.0;tmp=x+ 5.5;tmp -= (x+0.5)*log(tmp);ser=1.0;for (j =0;j < =5;j+ +) {x + = 1.0;ser + = cof[j]/x;}return -tmp+log(2.50662827465*ser);}/* The Gammq.c function called by main program.^ *1float gammq(a,x)float a,x;{float gainser,gammcf,g1n;void gcf(),gser(),nrerror();if (x < 0.0 I I a < = 0.0) nrerror("Invalid arguments in routine GAMMQ");if (x < (a+1.0)) {gser(&gamser,a,x,&gln);return 1.0-gamser;} else {gcfi&gammcf,a,x,&gln);return gammcf;}}/* The Gcf.c function.^ *1#define ITMAX 100#define EPS 3.0e-7void gcfigammcf,a,x,g1n)float a,x,*gammcf,*gln;{int n;float gold=0.0,g,fac=1.0,b1=1.0;float b0=0.0,anf,ana,an,al,a0=1.0;float gammln();void nrerror();146*gln=gammln(a);al =x;for (n=1;n< =ITMAX;n+ +) {an= (float) n;ana = an-a;a0 = (al +a0*ana)*fac;b0 = (Ill +b0*ana)*fac;anf= an*fac;al =x*a0+anf*al;bl =x*b0+anf*bl;if (al) {fac =1.0/al ;g=bl*fac;if (fabs((g-gold)/g) < EPS) {*gammcf= exp(-x +a*log(x)-(*gln))*g;return;}gold= g;}}nrerror("a too large, ITMAX too small in routine GCF");}#undef ITMAX#undef EPS/* The Gser.c function.#define ITMAX 100#define EPS 3.0e-7void gser(gamser,a,x,g1n)float a,x,*gamser,*gln;{int n;float sum,del,ap;float gammln();void nrerror();*gln=gammln(a);if (x < = 0.0) {if (x < 0.0) nrerror("x less than 0 in routine GSER");*gamser = 0.0;return;} else {ap = a;del = sum= 1.0/a;for (n=1;n< =ITMAX;n+ +) {ap + = 1.0;del *= x/ap;sum + = del;if (fabs(del) < fabs(sum)*EPS) {147*gamser=sum*expex +a*log(x)-(*gln));return;})nrerror("a too large, ITMAX too small in routine GSER");return;}}#undef ITMAX#undef EPS/* The Malloc.h function^ .1* malloc.h: This file is referenced by the preprocessor directive#include <malloc.h> . You should either (i) copy this file to theappropriate directory so that it can be found, or (ii) replacethe statements #include <malloc.h> by different #include statements,compatible with your compiler's handling of malloc() and free().^*1char *malloc();/* The Nrutil.c function .void nrerror(error_text)char error_text[];{void exit();fprintfistderr,"Numerical Recipes run-time error... \n");fprintfistderr," %s\n",error_text);fprintfistderr,"...now exiting to system... \n");exit(1);}float *vector(nl,nh)int nl,nh;{float *v;v= (float *)mallocaunsigned) (nh-n1+1)*sizeof(float));if (!v) nrerror("allocation failure in vector()");return v-nl;}int *ivector(nl,nh)int nl,nh;{int *v;v =(int *)malloc((unsigned) (nh-n1+1)*sizeof(int));if (!v) nrerror("allocation failure in ivector()");return v-nl;)double *dvector(nl,nh)int nl,nh;{double *v;v=(double *)mallocaunsigned) (nh-n1+1)*sizeof(double));if (!v) nrerror("allocation failure in dvector()");return v-nl;)float **matrix(nrl,nrh,ncl,nch)int nrl,nrh,ncl,nch;{int i;float **m;m= (float **) malloc((unsigned) (nrh-nrl + 1)*sizeof(float*));if (!m) nrerror("allocation failure 1 in matrix()");m -= nrl;for(i=nr1;i < =nrh;i+ +) {m[i]=(float *) malloc((unsigned) (nch-nc1+1)*sizeof(float));if (!m[i]) nrerror("allocation failure 2 in matrix()");m[i] -= ncl;}return m;}double **dmatrix(nrl,nrh,ncl,nch)int nrl,nrh,ncl,nch;{int i;double **m;m=(double **) malloc((unsigned) (nrh-nrl +1)*sizeof(double*));if (!m) nrerror("allocation failure 1 in dmatrix()");m -= nrl;for(i=nrl;i< =nrh;i+ +) {m[i] = (double *) malloc((unsigned) (nch-ncl + 1)*sizeof(double));if (!m[i]) nrerror("allocation failure 2 in dmatrix()");m[i] -= ncl;}return m;}int **imatrix(nrl,nrh,ncl,nch)149int nrl,nrh,ncl,nch;{int i,**m;m=(int **)malloc((unsigned) (nrh-nrI+1)*sizeof(int*));if (!m) nrerror("allocation failure 1 in imatrix()");m -= nrl;for(i=nrI;i< =nrh;i+ +) {m[i]= (int *)malloc((unsigned) (nch-ncl + 1)*sizeof(int));if (!m[i]) nrerror("allocation failure 2 in imatrix0");m[i] -= ncl;}return m;}float **submatrix(a,oldrl,oldrh,oldcl,oldch,newrl,newc1)float **a;int oldrl,oldrh,oldcl,oldch,newrl,newcl;{int i,j;float **m;m= (float **) malloc((unsigned) (oldrh-oldr1+1)*sizeof(float*));if (!m) nrerror("allocation failure in submatrix()");m -= newrl;for(i=oldrl,j =newr1;i < =oldrh;i+ +,j + +) m[j]=a[i] +oldcl-newcl;return m;}void free_vector(v,n1,nh)float *v;int nl,nh;{free((char*) (v +n1));}void free_ivector(v,n1,nh)int *v,n1,nh;{free((char*) (v +n1));}void free_dvector(v,n1,nh)double *v;int nl,nh;{free((char*) (v +n1));150}void free_matrix(m,nrl,nrh,ncl,nch)float **m;int nrl,nrh,ncl,nch;{int i;for(i=nrh;i > =nrl;i--) free((char*) (m[i]+ncl));free((char*) (m+nr1));}void free_ dmatrix(m,nrl,nrh,ncl,nch)double **m;int nrl,nrh,ncl,nch;int i;for(i=nrh;i > =nrl;i--) free((char*) (m[i]+ncl));free((char*) (m+nr1));}void free_imatrix(m,nrl,nrh,ncl,nch)int **m;int nrl,nrh,ncl,nch;{int i;for(i=nrh;i > =nrl;i--) free((char*) (m[i]+ncl));free((char*) (m+nr1));}void freesubmatrix(b,nrl,nrh,ncl,nch)float **b;int nrl,nrh,ncl,nch;{free((char*) (b+nr1));}float **convert_matrix(a,nrl,nrh,ncl,nch)float *a;int nrl,nrh,ncl,nch;{int i,j,nrow,ncol;float **m;nrow=nrh-nr1 +1;ncol = nch-ncl +1;151m = (float **) malloc((unsigned) (nrow)*sizeof(float*));if (!m) nrerror("allocation failure in convert_matrix0");m -= nrl;for(i =0,j =nr1;i < = nrow-1 ;i + + ,j + + ) m[j] =a +ncol*i-ncl;return m;}void freeconvert_matrix(b,nrl,nrh,ncl ,nch)float **b;int nrl,nrh,ncl,nch;{free((char*) (b + nrl));}/* The Nrutil.h function.float *vector();float **matrix();float **convert_matrix();double *dvector();double **dmatrix();int *ivector();int **imatrix();float **submatrix();void free_vector();void free_dvector();void free_ivector();void free_matrix();void free_dmatrix();void freeimatrix();void free_submatrix();void free_convert_matrix();void nrerror();
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- A new approach to mass balance modeling: applications...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
A new approach to mass balance modeling: applications to igneous petrology Ke, Peiwan 1992
pdf
Page Metadata
Item Metadata
Title | A new approach to mass balance modeling: applications to igneous petrology |
Creator |
Ke, Peiwan |
Date Issued | 1992 |
Description | An overdetermined linear system of equations is proposed to investigate mass balance relationships in igneous systems. A chi square (x²) minimization technique is used to solve the set of linear equations. The x² technique is superior to R² techniques because it considers the analytical errors on measurement data. Singular Value Decomposition (SVD) is used to search for the x² optimal solution. The attributes of the x² technique are discussed. The x² minimization fitting technique, combined with SVD methods allows the confidence regions to be calculated and expressed by a geometric shape centered on the x² optimal solution. The geometric shape varies with the dimension of the linear system: one-dimensional cases yield a line fragment; two-dimensional cases yield an ellipse and the three dimensions require an ellipsoid. The size of the confidence regions is determined by the confidence limits and is related to the coefficient matrix of the linear problem. Model with smaller confidence regions are considered to be better than models associated with larger confidence regions where the same confidence limit is used. Therefore, the confidence region can be used to select the model that best explains chemical differentiation in igneous system. Synthetic data sets have been analyzed using the new program MASSBAL.0 which solves the mass balance problems using the principles of the x² technique and SVD. The x² technique presents a solution closer to the true one than does the traditional residual minimization fitting (R² technique). It overcomes the deficiencies that the R² technique has and appears to be a more robust and sensitive fitting technique. The x² technique has been applied to the apparent clinopyroxene paradox formed in lavas from Diamond Craters, Oregon. Three models are proposed to explain the chemical difference: a three phase model (01+P1+Px), a two phase model (O1+Pl) and a one phase model (01). The x² mass balance modeling results show that the three phase model (O1+Pl+Px) has the smallest minimum x² value and the most compact confidence region; thus it is preferred and clinopyroxene did participate in the fractionation process. The x² technique has also been applied to study the assimilation problem in Paricutin volcano and the fractionation problem in PLA group of Mount Meager volcanic complex. The mass balance calculations show that assimilation has play a major role in Paricutin volcano and about 15% xenolith has contaminated the magma. In PLA fractionation processes, plagioclase and amphibole are shown to have caused the chemical difference between lavas, although other phenocrysts ( biotite and pyroxene ) did not participate. |
Extent | 7491663 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Notes | 1 computer disk - University Archives - AW5 .B71 1993-0107 |
Date Available | 2008-09-30 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0052616 |
URI | http://hdl.handle.net/2429/2407 |
Degree |
Master of Science - MSc |
Program |
Geological Sciences |
Affiliation |
Science, Faculty of Earth, Ocean and Atmospheric Sciences, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 1993-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_1993_spring_ke_peiwan.pdf [ 7.14MB ]
- Metadata
- JSON: 831-1.0052616.json
- JSON-LD: 831-1.0052616-ld.json
- RDF/XML (Pretty): 831-1.0052616-rdf.xml
- RDF/JSON: 831-1.0052616-rdf.json
- Turtle: 831-1.0052616-turtle.txt
- N-Triples: 831-1.0052616-rdf-ntriples.txt
- Original Record: 831-1.0052616-source.json
- Full Text
- 831-1.0052616-fulltext.txt
- Citation
- 831-1.0052616.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0052616/manifest