MAGEMin_C v1.4.8 (#42)

- added the option to use aH2O, aO2, aTiO2, aMgO, aFeO and aAl2O3 as buffers
- added trace-element partitioning routine using O. Laurent Kd's (PhD 2012)
- added several zircon saturation models (Watson & Harrison, 1983; Boehnke et al., 2013; Crisp and Berry 2022)
- added routine to save minimization results as CSV
- corrected normalization when hitting liquidus

Project.toml

@@ -1,17 +1,22 @@
 name = "MAGEMin_C"
 uuid = "e5d170eb-415a-4524-987b-12f1bce1ddab"
 authors = ["Boris Kaus <[email protected]> & Nicolas Riel <[email protected]>"]
-version = "1.4.7"
+version = "1.4.8"
 
 [deps]
 CEnum = "fa961155-64e5-5f13-b03f-caf6b980ea82"
+CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
+DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
+Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"
 Libdl = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
 MAGEMin_jll = "763ebaa8-b0d2-5f6b-90ef-4fc23b5db1c4"
 ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [compat]
 CEnum = "0.4"
-MAGEMin_jll = "1.4.4 - 1.4.4"
+DataFrames = "1.6.1"
+CSV = "0.10.14"
+MAGEMin_jll = "1.4.5 - 1.4.5"
 ProgressMeter = "1"
 julia = "1.6"

gen/magemin_library.jl

@@ -582,6 +582,7 @@ mutable struct global_variables
     numPoint::Cint
     global_ite::Cint
     H2O_id::Cint
+    Al2O3_id::Cint
     K2O_id::Cint
     O_id::Cint
     TiO2_id::Cint
@@ -1412,7 +1413,7 @@ mutable struct metapelite_datasets
     n_pp::Cint
     n_ss::Cint
     ox::NTuple{11, NTuple{20, Cchar}}
-    PP::NTuple{17, NTuple{20, Cchar}}
+    PP::NTuple{23, NTuple{20, Cchar}}
     SS::NTuple{16, NTuple{20, Cchar}}
     verifyPC::NTuple{16, Cint}
     n_SS_PC::NTuple{16, Cint}
@@ -1438,7 +1439,7 @@ mutable struct metabasite_datasets
     n_pp::Cint
     n_ss::Cint
     ox::NTuple{10, NTuple{20, Cchar}}
-    PP::NTuple{18, NTuple{20, Cchar}}
+    PP::NTuple{24, NTuple{20, Cchar}}
     SS1::NTuple{14, NTuple{20, Cchar}}
     verifyPC1::NTuple{14, Cint}
     n_SS_PC1::NTuple{14, Cint}
@@ -1468,7 +1469,7 @@ mutable struct igneous_datasets
     n_pp::Cint
     n_ss::Cint
     ox::NTuple{11, NTuple{20, Cchar}}
-    PP::NTuple{17, NTuple{20, Cchar}}
+    PP::NTuple{23, NTuple{20, Cchar}}
     SS::NTuple{15, NTuple{20, Cchar}}
     verifyPC::NTuple{15, Cint}
     n_SS_PC::NTuple{15, Cint}
@@ -1494,7 +1495,7 @@ mutable struct ultramafic_datasets
     n_pp::Cint
     n_ss::Cint
     ox::NTuple{7, NTuple{20, Cchar}}
-    PP::NTuple{15, NTuple{20, Cchar}}
+    PP::NTuple{21, NTuple{20, Cchar}}
     SS::NTuple{12, NTuple{20, Cchar}}
     verifyPC::NTuple{12, Cint}
     n_SS_PC::NTuple{12, Cint}

julia/TE_Zr_tests.jl

@@ -0,0 +1,30 @@
+using MAGEMin_C
+
+# Initialize database  - new way
+dtb          =  "ig"
+data        =   Initialize_MAGEMin(dtb, verbose=false);
+test        =   0         #KLB1
+data        =   use_predefined_bulk_rock(data, test);
+
+# Call optimization routine for given P & T & bulk_rock
+P           =   8.0
+T           =   1400.0
+out         =   point_wise_minimization(P,T, data);
+
+elem_TE     = ["SiO2","TiO2","Al2O3","O","FeO","MnO","MgO","CaO","K2O","Na2O","H2O","Li","Be","B","Sc","V","Cr","Ni","Cu","Zn","Rb","Sr","Y","Zr","Nb","Cs","Ba","La","Ce","Pr","Nd","Sm","Eu","Gd","Tb","Dy","Ho","Er","Tm","Yb","Lu","Hf","Ta","Pb","Th","U"]
+bulk_TE     = [50.60777553,0.953497243,13.70435413,0.19,11.28130762,0.202560796,8.496024312,9.502380068,0.700881685,2.07434927,4,29.14258603,0.434482763,29.69200003,38.23663423,257.4346716,529.333066,208.2057375,88.87615683,91.7592182,16.60777308,163.4533209,20.74016207,66.90677472,3.808354064,1.529226981,122.8449739,6.938172601,16.04827796,2.253943183,10.18276823,3.3471043,0.915941652,3.28230146,1.417695298,3.851230952,0.914966282,2.20425,0.343734976,2.136202593,0.323405135,1.841502082,0.330971265,5.452969044,1.074692058,0.290233271]
+
+KDs_dtb     = get_OL_KDs_database();
+
+C0          = adjust_chemical_system(    KDs_dtb,bulk_TE,elem_TE);
+
+out_TE      = TE_prediction(C0,KDs_dtb, "CB",out,dtb);
+
+
+if ~isnothing(out_TE.bulk_cor_wt)      #then we can recompute the equilibrium after removing the SiO2 entering zircon
+    sys_in      = "wt"    
+    out         = single_point_minimization(P, T, data, X=out_TE.bulk_cor_wt, Xoxides=out.oxides, sys_in=sys_in)
+    out_TE      = TE_prediction(C0,KDs_dtb, "CB",out,dtb)
+end
+
+# add parallel test

julia/Zircon_saturation.jl

@@ -1,11 +1,7 @@
 # Routine to compute zircon saturation and adjust bulk-rock composition when zircon crystallizes
 # NR 12/04/2023
 
-# include("TE_routine.jl")
-
-# using MAGEMin_C  
-
-function zirconium_saturation(  out     :: gmin_struct{Float64, Int64}; 
+function zirconium_saturation(  out     :: MAGEMin_C.gmin_struct{Float64, Int64}; 
                                 model   :: String = "WH"    )
 
     if out.frac_M > 0.0                            
@@ -63,16 +59,18 @@ function zirconium_saturation(  out     :: gmin_struct{Float64, Int64};
     return C_zr_liq 
 end
 
+function adjust_bulk_4_zircon(  zr_liq  ::  Float64,
+                                sat_liq ::  Float64 )
 
+    SiO2_wt         = 0.0   
+    O2_wt           = 0.0        
+    zircon_wt       = 0.0             
+    zircon_excess   = (zr_liq - sat_liq)/1e4
 
-# dtb         = "ig"
-# data        =  Initialize_MAGEMin(dtb, verbose=true);
-# test        =  0 
-# data        =  use_predefined_bulk_rock(data, test);
-# P           =  2.0
-# T           =  1800.0
-# out         =  point_wise_minimization(P,T, data);
+    zircon_wt   = zircon_excess*0.497644
+    SiO2_wt     = (zircon_wt *0.327765)
+    O2_wt       = (zircon_wt *0.174570)
 
 
-# C_zr_liq    = zirconium_saturation( out; 
-#                                     model="CB")
+    return zircon_wt ,SiO2_wt ,O2_wt 
+end