Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

MAGEMin_C 1.6.5 #63

Merged
merged 2 commits into from
Dec 3, 2024
Merged

MAGEMin_C 1.6.5 #63

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
5bfc368
MAGEMin_C 1.6.5
NicolasRiel Dec 2, 2024
a1ef05c
waiting for MAGEMin_jll 1.6.1
NicolasRiel Dec 3, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions Project.toml
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "MAGEMin_C"
uuid = "e5d170eb-415a-4524-987b-12f1bce1ddab"
authors = ["Nicolas Riel <[email protected]> & Boris Kaus <[email protected]>"]
version = "1.6.4"
authors = ["Nicolas Riel <[email protected]> & Boris Kaus <[email protected]>"]
version = "1.6.5"

[deps]
CEnum = "fa961155-64e5-5f13-b03f-caf6b980ea82"
Expand All @@ -17,6 +17,6 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
CEnum = "0.4"
CSV = "0.10.14"
DataFrames = "1.6.1"
MAGEMin_jll = "1.6.0 - 1.6.0"
MAGEMin_jll = "1.6.1 - 1.6.1"
ProgressMeter = "1"
julia = "1.6"
9 changes: 6 additions & 3 deletions gen/magemin_library.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1177,10 +1177,13 @@ struct stb_systems
rho_F::Cdouble
bulk_S_wt::Ptr{Cdouble}
frac_S_wt::Cdouble
frac_S_vol::Cdouble
bulk_M_wt::Ptr{Cdouble}
frac_M_wt::Cdouble
frac_M_vol::Cdouble
bulk_F_wt::Ptr{Cdouble}
frac_F_wt::Cdouble
frac_F_vol::Cdouble
n_ph::Cint
n_PP::Cint
n_SS::Cint
Expand Down Expand Up @@ -1323,7 +1326,7 @@ mutable struct metapelite_datasets
n_pp::Cint
n_ss::Cint
ox::NTuple{11, NTuple{20, Cchar}}
PP::NTuple{23, NTuple{20, Cchar}}
PP::NTuple{24, NTuple{20, Cchar}}
SS::NTuple{18, NTuple{20, Cchar}}
verifyPC::NTuple{18, Cint}
n_SS_PC::NTuple{18, Cint}
Expand All @@ -1349,7 +1352,7 @@ mutable struct metabasite_datasets
n_pp::Cint
n_ss::Cint
ox::NTuple{10, NTuple{20, Cchar}}
PP::NTuple{23, NTuple{20, Cchar}}
PP::NTuple{24, NTuple{20, Cchar}}
SS::NTuple{17, NTuple{20, Cchar}}
verifyPC::NTuple{17, Cint}
n_SS_PC::NTuple{17, Cint}
Expand Down Expand Up @@ -1505,7 +1508,7 @@ mutable struct metapelite_datasets_ext
n_pp::Cint
n_ss::Cint
ox::NTuple{13, NTuple{20, Cchar}}
PP::NTuple{25, NTuple{20, Cchar}}
PP::NTuple{26, NTuple{20, Cchar}}
SS::NTuple{23, NTuple{20, Cchar}}
verifyPC::NTuple{23, Cint}
n_SS_PC::NTuple{23, Cint}
Expand Down
20 changes: 14 additions & 6 deletions julia/MAGEMin_wrappers.jl
View file
Open in desktop

Large diffs are not rendered by default.

2 changes: 1 addition & 1 deletion src/MAGEMin.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -914,7 +914,7 @@ global_variable SetupDatabase( global_variable gv,

}

if (gv.verbose == 1){
if (gv.verbose == 2){
printf("\n");
printf("--verbose : verbose = %i \n", gv.verbose );
printf("--rg : research group = %s \n", gv.research_group );
Expand Down
8 changes: 4 additions & 4 deletions src/MAGEMin.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -769,10 +769,10 @@ typedef struct stb_systems {
double *bulk_M; double frac_M; double rho_M; /* Melt system informations */
double *bulk_F; double frac_F; double rho_F; /* Fluid system informations */

double *bulk_S_wt; double frac_S_wt; /* Solid system informations */
double *bulk_M_wt; double frac_M_wt; /* Melt system informations */
double *bulk_F_wt; double frac_F_wt; /* Fluid system informations */
double *bulk_S_wt; double frac_S_wt; double frac_S_vol; /* Solid system informations */
double *bulk_M_wt; double frac_M_wt; double frac_M_vol; /* Melt system informations */
double *bulk_F_wt; double frac_F_wt; double frac_F_vol; /* Fluid system informations */

int n_ph; /* number of predicted stable phases */
int n_PP; /* number of predicted stable pure phases */
int n_SS; /* number of predicted stable solution phases */
Expand Down
2 changes: 1 addition & 1 deletion src/TC_database/NLopt_opt_function.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7718,7 +7718,7 @@ SS_ref NLopt_opt_aq17_function(global_variable gv, SS_ref SS_ref_db){
SS_ref_db.ub[i] = SS_ref_db.bounds[i][1];
}

SS_ref_db.opt = nlopt_create(NLOPT_LN_COBYLA, (n));
SS_ref_db.opt = nlopt_create(NLOPT_LD_SLSQP, (n));
nlopt_set_lower_bounds(SS_ref_db.opt, SS_ref_db.lb);
nlopt_set_upper_bounds(SS_ref_db.opt, SS_ref_db.ub);
nlopt_set_min_objective(SS_ref_db.opt, obj_aq17, &SS_ref_db);
Expand Down
4 changes: 2 additions & 2 deletions src/TC_database/TC_gem_function.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1035,11 +1035,11 @@ PP_ref G_FS_function( int len_ox,
PP_ref_db.factor = factor;
PP_ref_db.phase_shearModulus = 0.0;

// printf(" %4s %+10f | factor: %+10f\n",name,PP_ref_db.gbase,PP_ref_db.factor);
// printf(" %4s %+10f %+10f | factor: %+10f\n",name,G,PP_ref_db.gbase,PP_ref_db.factor);
// for (i = 0; i < len_ox; i++){
// printf("%+10f",PP_ref_db.Comp[i]*PP_ref_db.factor);
// }
// printf("\n");
// printf("\n\n");

return (PP_ref_db);
}
12 changes: 6 additions & 6 deletions src/TC_database/TC_init_database.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -35,10 +35,10 @@ oxide_data oxide_info = {
metapelite_dataset metapelite_db = {
62, /* Endmember default dataset number */
11, /* number of oxides */
23, /* number of pure phases */
24, /* number of pure phases */
17, /* number of solution phases */
{"SiO2" ,"Al2O3","CaO" ,"MgO" ,"FeO" ,"K2O" ,"Na2O" ,"TiO2" ,"O" ,"MnO" ,"H2O" },
{"q" ,"crst" ,"trd" ,"coe" ,"stv" ,"ky" ,"sill" ,"and" ,"ru" ,"sph" ,"O2" ,"H2O" ,
{"q" ,"crst" ,"trd" ,"coe" ,"stv" ,"ky" ,"sill" ,"and" ,"ru" ,"sph" ,"O2" ,"H2O" ,"zo" ,
"qfm" ,"qif" ,"nno" ,"hm" ,"cco" ,"aH2O" , "aO2" ,"aMgO" ,"aFeO" ,"aAl2O3" ,"aTiO2" },
{"liq" ,"fsp" ,"bi" ,"g" ,"ep" ,"ma" ,"mu" ,"opx" ,"sa" ,"cd" ,"st" ,"chl" ,"ctd" ,"sp" ,"mt" ,"ilm" ,"ilmm" ,"aq17" },

Expand All @@ -63,10 +63,10 @@ metapelite_dataset metapelite_db = {
metabasite_dataset metabasite_db = {
62, /* Endmember default dataset number */
10, /* number of oxides */
23, /* number of pure phases */
24, /* number of pure phases */
17, /* number of solution phases */
{"SiO2" ,"Al2O3","CaO" ,"MgO" ,"FeO" ,"K2O" ,"Na2O" ,"TiO2" ,"O" ,"H2O" },
{"q" ,"crst" ,"trd" ,"coe" ,"law" ,"ky" ,"sill" ,"and" ,"ru" ,"sph" ,"ab" ,"H2O" ,
{"q" ,"crst" ,"trd" ,"coe" ,"law" ,"ky" ,"sill" ,"and" ,"ru" ,"sph" ,"ab" ,"H2O" ,"zo" ,
"qfm" ,"qif" ,"nno" ,"hm" ,"cco" ,"aH2O" ,"aO2" ,"aMgO" ,"aFeO" ,"aAl2O3" ,"aTiO2" },
{"sp" ,"opx" ,"fsp" ,"liq" ,"mu" ,"ilmm" ,"ilm" ,"ol" ,"hb" ,"ep" ,"g" ,"chl" ,"bi" ,"dio" ,"aug" ,"abc" ,"spn" },

Expand Down Expand Up @@ -232,10 +232,10 @@ mantle_dataset mantle_db = {
metapelite_dataset_ext metapelite_ext_db = {
62, /* Endmember default dataset number */
13, /* number of oxides */
25, /* number of pure phases */
26, /* number of pure phases */
23, /* number of solution phases */
{"SiO2" ,"Al2O3","CaO" ,"MgO" ,"FeO" ,"K2O" ,"Na2O" ,"TiO2" ,"O" ,"MnO" ,"H2O" ,"CO2" ,"S" },
{"q" ,"crst" ,"trd" ,"coe" ,"stv" ,"ky" ,"sill" ,"and" ,"ru" ,"sph" ,"O2" ,"pyr" ,"gph" ,"law" ,
{"q" ,"crst" ,"trd" ,"coe" ,"stv" ,"ky" ,"sill" ,"and" ,"ru" ,"sph" ,"O2" ,"pyr" ,"gph" ,"law" ,"zo" ,
"qfm" ,"qif" ,"nno" ,"hm" ,"cco" ,"aH2O" , "aO2" ,"aMgO" ,"aFeO" ,"aAl2O3" ,"aTiO2" },
{"liq" ,"fsp" ,"bi" ,"g" ,"ep" ,"ma" ,"mu" ,"opx" ,"sa" ,"cd" ,"st" ,"chl" ,"ctd" ,"sp" ,"mt" ,"ilm" ,"ilmm" ,"occm" ,"fl" ,"po" ,"dio" ,"aug" ,"hb" },

Expand Down
6 changes: 3 additions & 3 deletions src/TC_database/TC_init_database.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@
int n_pp;
int n_ss;
char ox[11][20];
char PP[23][20];
char PP[24][20];
char SS[18][20];

int verifyPC[18];
Expand Down Expand Up @@ -54,7 +54,7 @@
int n_pp;
int n_ss;
char ox[10][20];
char PP[23][20];
char PP[24][20];

char SS[17][20];
int verifyPC[17];
Expand Down Expand Up @@ -246,7 +246,7 @@
int n_pp;
int n_ss;
char ox[13][20];
char PP[25][20];
char PP[26][20];
char SS[23][20];

int verifyPC[23];
Expand Down
59 changes: 59 additions & 0 deletions src/TC_database/aqueous.md
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,59 @@
# Aqueous species formulation, using Miron et al., 2017 database.


```math
G_{aq} = \sum_{i=1}^{n^{em}} n_i \mu_i
```

where $n_i$ the molar fraction and $\mu_i$ is the chemical potential of the species $i$.

```math
\mu_i = G_i^0 + RT log(\gamma_i m_i)
```

where $m_i$ is the molal fraction (!?)

And then develops into:

```math
\mu_i = G_i^0 + RT( log(\gamma_i) + log(m_i))
```

---

For solute (aqueous species, excluding water), and assuming $\gamma_i$ = 1 (for now)

```math
\mu_i = G_i^0 + RT log(m_i)
```

where

```math
m_i = \frac{n_i}{m_{H2O}}
```

so here $n_i$ is again the molar fraction and $m_{H2O}$ is mass of water? Should it be water density?

I could find the following definition:

```math
m_i = \frac{n_i}{n_{H2O} M_{H2O}}
```
and here $n_i$ is the molarity and $M_{H2O}$ is the Molar mass of water (18.015). However, here density of water like we discussed is not appearing, is that right?

Should it simply be:

```math
m_i = \frac{n_i}{\rho_{H2O}}
```

---

For solvant (water):

```math
\mu_i = G_i^0 + RT log(\frac{1}{M_{H2O}})
```

!?
5 changes: 2 additions & 3 deletions src/TC_database/objective_functions.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5645,7 +5645,6 @@ double obj_mp_sp(unsigned n, const double *x, double *grad, void *SS_ref_db){
sf[3] = 1.0 - x[0];
sf[4] = x[0];


mu[0] = R*T*creal(clog(sf[0]*sf[4])) + gb[0] + mu_Gex[0];
mu[1] = R*T*creal(clog(sf[0]*sf[3])) + gb[1] + mu_Gex[1];
mu[2] = R*T*creal(clog(sf[4]*sf[1] + d_em[2])) + gb[2] + mu_Gex[2];
Expand Down Expand Up @@ -14523,7 +14522,7 @@ double obj_aq17(unsigned n, const double *x, double *grad, void *SS_ref_db) {
d->g,
d->epsilon,
xiw ); //fraction of water
// cor = HSC_to_SUPCRT(ElH, Comp[i], len_ox);
cor = HSC_to_SUPCRT(ElH, Comp[i], len_ox);
mu[i] = gb[i] + (log(pow(10.0,loggamma)) + log(1000.0/18.0153) + log(x[i]/Xw) - log(xiw/Xw) - xiw/Xw + 1.0 )/1000.0;
m_all += x[i];
if (charge[i] != 0.0){
Expand Down Expand Up @@ -14552,7 +14551,7 @@ double obj_aq17(unsigned n, const double *x, double *grad, void *SS_ref_db) {
m_all );

/* set chemical potential of water */
// cor = HSC_to_SUPCRT(ElH, Comp[0], len_ox);
cor = HSC_to_SUPCRT(ElH, Comp[0], len_ox);
mu[0] = gb[0] + ( log(logawater) + log(xiw/Xw) - Xw/xiw - xiw/Xw + 2.0)/1000.0;

px_aq17(SS_ref_db,x);
Expand Down
24 changes: 12 additions & 12 deletions src/TC_database/tc_gss_function.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4246,12 +4246,12 @@ SS_ref G_SS_mp_sp_function(SS_ref SS_ref_db, char* research_group, int EM_datase
};
//int si = sizeof(SF_tmp)/sizeof(SF_tmp[0]); printf("sz_tmp: %d n_sf %d\n",si,SS_ref_db.n_sf);

SS_ref_db.W[0] = 16.;
SS_ref_db.W[1] = 2.;
SS_ref_db.W[2] = 20.;
SS_ref_db.W[3] = 18.;
SS_ref_db.W[4] = 36.;
SS_ref_db.W[5] = 30.;
SS_ref_db.W[0] = 0.0;
SS_ref_db.W[1] = 18.5;
SS_ref_db.W[2] = 27.0;
SS_ref_db.W[3] = 40.0;
SS_ref_db.W[4] = 30.0;
SS_ref_db.W[5] = 0.0;


em_data herc_eq = get_em_data( research_group, EM_dataset,
Expand Down Expand Up @@ -13401,12 +13401,12 @@ SS_ref G_SS_mpe_sp_function(SS_ref SS_ref_db, char* research_group, int EM_datas
};
//int si = sizeof(SF_tmp)/sizeof(SF_tmp[0]); printf("sz_tmp: %d n_sf %d\n",si,SS_ref_db.n_sf);

SS_ref_db.W[0] = 16.;
SS_ref_db.W[1] = 2.;
SS_ref_db.W[2] = 20.;
SS_ref_db.W[3] = 18.;
SS_ref_db.W[4] = 36.;
SS_ref_db.W[5] = 30.;
SS_ref_db.W[0] = 0.0;
SS_ref_db.W[1] = 18.5;
SS_ref_db.W[2] = 27.0;
SS_ref_db.W[3] = 40.0;
SS_ref_db.W[4] = 30.0;
SS_ref_db.W[5] = 0.0;


em_data herc_eq = get_em_data( research_group, EM_dataset,
Expand Down
2 changes: 1 addition & 1 deletion src/TC_database/tc_gss_init_function.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1730,7 +1730,7 @@ SS_ref G_SS_mpe_opx_init_function(SS_ref SS_ref_db, global_variable gv){
SS_ref_db.n_cat = 0;
SS_ref_db.is_liq = 0;
SS_ref_db.symmetry = 0;
SS_ref_db.n_sf = 10;
SS_ref_db.n_sf = 11;
SS_ref_db.n_em = 7;
SS_ref_db.n_v = 7;
SS_ref_db.n_w = 21;
Expand Down
32 changes: 32 additions & 0 deletions src/dump_function.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -173,6 +173,9 @@ void reset_output_struct( global_variable gv,
sp[0].frac_S_wt = 0.0;
sp[0].frac_M_wt = 0.0;
sp[0].frac_F_wt = 0.0;
sp[0].frac_S_vol = 0.0;
sp[0].frac_M_vol = 0.0;
sp[0].frac_F_vol = 0.0;
sp[0].rho_S = 0.0;
sp[0].rho_M = 0.0;
sp[0].rho_F = 0.0;
Expand Down Expand Up @@ -570,6 +573,7 @@ void fill_output_struct( global_variable gv,
if (gv.n_phase == 1){
sp[0].frac_M = 1.0;
sp[0].frac_M_wt = 1.0;
sp[0].frac_M_vol = 1.0;
}
else{
sp[0].frac_M = cp[i].ss_n_mol;
Expand Down Expand Up @@ -686,6 +690,19 @@ void fill_output_struct( global_variable gv,
for (int i = 0; i < gv.len_cp; i++){
if ( cp[i].ss_flags[1] == 1){
sp[0].ph_frac_vol[n] = sp[0].ph_frac_wt[n] / sp[0].SS[n].rho;

if (strcmp( cp[i].name, "liq") == 0 || strcmp( cp[i].name, "fl") == 0 ){
if (strcmp( cp[i].name, "liq") == 0){
sp[0].frac_M_vol = sp[0].ph_frac_vol[n];
}
else{
sp[0].frac_F_vol = sp[0].ph_frac_vol[n];
}
}
else{
sp[0].frac_S_vol += sp[0].ph_frac_vol[n];
}

sum_vol += sp[0].ph_frac_vol[n];
sum_mol += sp[0].ph_frac[n];
sum_wt += sp[0].ph_frac_wt[n];
Expand All @@ -696,6 +713,14 @@ void fill_output_struct( global_variable gv,
for (int i = 0; i < gv.len_pp; i++){
if (gv.pp_flags[i][1] == 1 && gv.pp_flags[i][4] == 0){
sp[0].ph_frac_vol[n] = sp[0].ph_frac_wt[n] / sp[0].PP[m].rho;

if (strcmp( gv.PP_list[i], "H2O") != 0){
sp[0].frac_S_vol += sp[0].ph_frac_vol[n];
}
if (strcmp( gv.PP_list[i], "H2O") == 0){
sp[0].frac_F_vol = sp[0].ph_frac_vol[n];
}

sum_vol += sp[0].ph_frac_vol[n];
sum_mol += sp[0].ph_frac[n];
sum_wt += sp[0].ph_frac_wt[n];
Expand All @@ -709,6 +734,7 @@ void fill_output_struct( global_variable gv,
sp[0].ph_frac_wt[i] /= sum_wt;
}


/* The following section normalizes the entries for S (solid), M (melt) and F (fluid) which are entries useful for geodynamic coupling */
// normalize rho_S and bulk_S
if (sp[0].frac_S > 0.0){
Expand Down Expand Up @@ -763,6 +789,12 @@ void fill_output_struct( global_variable gv,
sp[0].frac_M_wt /= sum;
sp[0].frac_S_wt /= sum;

sum = sp[0].frac_F_vol + sp[0].frac_M_vol + sp[0].frac_S_vol;

sp[0].frac_F_vol /= sum;
sp[0].frac_M_vol /= sum;
sp[0].frac_S_vol /= sum;

/* compute cp as J/K/kg for given bulk-rock composition */
sp[0].s_cp = sp[0].cp_wt/mass_bulk*1e6;
if (strcmp(gv.research_group, "tc") == 0 ){
Expand Down
Loading
Loading