Cpx-Pl-Qz: Perchuk, 2003

Ab ↔ Jd+Qz.

Использовать Проверить

source('recalc_cpx_sturm.r');
source('recalc_pl_tpf_realc0.r');

sensors.p.cpx_pl_perchuk <- function(Cpx, Pl, Temp=600, usePFU = FALSE, ...){

Name <- paste0(Cpx$Name,'-',Pl$Name);

if(!usePFU)
{
Cpx <- recalc.cpx_sturm(Cpx);
Pl <- recalc.pl_tpf_realc0(Pl);
}

Cpx_pfu <- selectNames(Cpx,c('Mg_pfu','Fe_p2_pfu','Na_pfu','Fe_p3_pfu'));
Pl_pfu <- selectNames(Pl,c('Na_pfu','Ca_pfu','K_pfu'));

X_Mg = Cpx_pfu$Mg_pfu/(Cpx_pfu$Mg_pfu+Cpx_pfu$Fe_p2_pfu);
X_Jd = Cpx_pfu$Na_pfu - Cpx_pfu$Fe_p3_pfu;
X_Ab = Pl_pfu$Na_pfu/(Pl_pfu$Na_pfu+Pl_pfu$K_pfu+Pl_pfu$Ca_pfu);

X_Di = X_Mg*(1-X_Jd)
X_Hd = 1-X_Jd-X_Di

T_K = Temp+273

W31 = 13270-9.3*T_K
W13 = 16570-13.5*T_K
W21 = 15090-8*T_K
W12 = 3670-9*T_K

G1 = W12*(2*X_Hd^3-X_Hd^2)+W21*2*(X_Hd^2-X_Hd^3)+W13*(2*X_Di^3-X_Di^2)+W31*2*(X_Di^2-X_Di^3)+0.5*X_Jd*X_Di*X_Hd*(W12+W21+W13+W31)
X_ord = ifelse(T_K<1173,sqrt(1-T_K/1138)/2,0)
Tc = ifelse(X_Jd<0,5.2*X_Jd*1138,2*(1-X_Jd)*1138)
Z = ifelse(T_K<1138,(4*X_Jd*(1-X_Jd)*1138-T_K)/1138,0)
y2 = ifelse(Z>0,SQRT(Z),0)
G2 = X_Mg*(-11.3/3*y2*(2*1138*(1-X_Jd)*(4*(1-X_Jd)-1)-T_K))

a0 = 5.479
b = 6854
aod0 = 41.62
bod = -9301
cod = 43600
Tcod = 824.1 
Tc2 = 1251
Q = -1.9*(10^(-6))*T_K^2+0.0021*T_K+0.3795
Qod = -0.0012*T_K+1.5127
d0 = -2.17
d1 = -3.043
d2 = -1.569*10^(-3)
d3 = 2.109*10^(-6)

F12 = ifelse(T_K<1290,
             0.5*a0*(T_K-Tc2)*Q^2+
             0.25*b*Q^4+
             0.5*aod0*(T_K-Tcod)*Qod^2+
             0.25*bod*Qod^4+
             1/6*cod*Qod^6+
             (d0+d1*T_K+d2*T_K^2+d3*T_K^3)*Q*Qod,
             0)

GG = G1+G2

Ge_Ab = 4.18*(1-X_Ab)^2*((6860-3.874*T_K)-2*X_Ab*(4880-2.348*T_K))

Pcalc = 10^-3*(1.827+0.036*X_Ab*Qod-
    (
        (1.827+0.036*X_Ab*Qod)^2-
        1.24*10^-5*(-9360+44*T_K+8.314*T_K*log(X_Jd/X_Ab)+GG-X_Ab*F12-Ge_Ab)
    )^0.5)/(6.2*10^-6)

return (data.frame(Name=Name, T_seq=Temp, P_calc_cpx_pl_perchuk=Pcalc, X_Jd = X_Jd, X_Ab = X_Ab))
}

sensors.p.cpx_pl_perchuk.args <- function() {
    ars <- list()
    ars$minerals = c('Cpx','Pl')
    ars$ranges = c('Temp')
    return(ars)
}