pygcc.solid_solution ==================== .. py:module:: pygcc.solid_solution .. autoapi-nested-parse:: Created on Wed Mar 17 16:02:22 2021 @author: adedapo.awolayo and Ben Tutolo, University of Calgary Copyright (c) 2020 - 2021, Adedapo Awolayo and Ben Tutolo, University of Calgary This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program. If not, see . Attributes ---------- .. autoapisummary:: pygcc.solid_solution.MW pygcc.solid_solution.J_to_cal Classes ------- .. autoapisummary:: pygcc.solid_solution.solidsolution_thermo Functions --------- .. autoapisummary:: pygcc.solid_solution.roundup_tenth pygcc.solid_solution.Molecularweight Module Contents --------------- .. py:function:: roundup_tenth(x) .. py:function:: Molecularweight() This function stores the Molecular weight of all elements .. py:data:: MW .. py:data:: J_to_cal :value: 4.184 .. py:class:: solidsolution_thermo(**kwargs) Class to calculate Ideal mixing model for solid-solutions, supporting multisite ideal formalism :param X: End member mineral volume fraction or mole fraction of Mg in clinopyroxene :type X: float :param cpx_Ca: number of moles of Ca in clinopyroxene formula unit (=1 for Di, Hed) :type cpx_Ca: float :param T: Temperature [°C] :type T: float, vector :param P: Pressure [bar] :type P: float, vector :param dbaccessdic: dictionary of species from direct-access database, optional, default is speq23 :type dbaccessdic: dict :param solidsolution_type: specify either 'All' or 'plagioclase' or 'olivine' or 'pyroxene' or 'cpx' or 'alk-feldspar' or 'biotite' to carry-out all or any solid-solution calculations, default is 'All' :type solidsolution_type: string :param Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: string ThermoInUnit : string specify either 'cal' or 'KJ' as the input units for species properties (optional), particularly used to covert KJ data to cal by supcrtaq function if not specified default - 'cal' Al_Si : string specify either 'pygcc' or 'Arnórsson_Stefánsson' as the input to express Al and Si species in solid solution (optional), 'Arnórsson_Stefánsson' expresses them as 'Al(OH)4-' and 'H4SiO4(aq)', respectively while pygcc uses 'Al3+' and 'SiO2(aq)', respectively if not specified default - 'pygcc' rhoEG : dict dictionary of water properties like density (rho), dielectric factor (E) and Gibbs Energy (optional) :returns: * *if 'All' is specified for solidsolution_type, each case has the follwowing prefix - AnAb_, FoFa_, EnFe_, cpx_, AlkFed_, Biotite_* * **logK** (*float, vector*) -- logarithmic K values * **Rxn** (*dict*) -- The calculated dictionary of reaction thermodynamic properties has the following properties: * type: solid-solution mineral type, [-] * name: solid-solution name, [K] * formula: solid-solution mineral formula, [-] * MW: Molecular weight, [g/mol] * min: solid-solution mineral properties, ['formula', 'source date', dG[cal/ml], dH[cal/mol], S[cal/mol-K], V[cm3/mol], a[cal/mol-K], b[10^3 cal/mol/K^2], c[10^-5 cal/mol/K]] * spec: list of species, [-] * coeff: list of corresponding coefficients of species above, [-] * nSpec: Total number of species, [-] * V: Molar volume, [cm3/mol] * source: Source of thermo data, [kJ/kg·K] * elements: list of elements and their total numbers, [-] Usage ---------- The general usage of water_dielec is as follows: (1) For water dielectric properties at any Temperature and Pressure: ss = solidsolution_thermo(T = T, P = P, X = X, solidsolution_type = 'All', Dielec_method = 'JN91'), where T is temperature in celsius and P is pressure in bar (2) For water dielectric properties at any Temperature and density : ss = solidsolution_thermo(T = T, rho = rho, X = X, solidsolution_type = 'All', Dielec_method = 'JN91'), where T is temperature in celsius and rho is density in kg/m³ (3) For water dielectric properties at any Temperature and Pressure on steam saturation curve: ss = solidsolution_thermo(T = T, P = 'T', X = X, solidsolution_type = 'All', Dielec_method = 'JN91'), where T is temperature in celsius and P is assigned a quoted character 'T' to reflect steam saturation pressure ss = solidsolution_thermo(P = P, T = 'P', X = X, solidsolution_type = 'All', Dielec_method = 'JN91'), where P is pressure in bar and T is assigned a quoted character 'P' to reflect steam saturation temperature .. rubric:: Examples >>> ss = solidsolution_thermo(cpx_Ca = 0.5, X = 0.85, T = 30, P = 50, solidsolution_type = 'All') >>> ss.AnAb_logK, ss.FoFa_logK, ss.EnFe_logK, ss.cpx_logK, ss.AlkFed_logK, ss.Biotite_logK 21.750590, 26.570902, 10.793372, 20.577152, 3.756792, 32.96363 >>> ss = solidsolution_thermo(cpx_Ca = 1, X = 0.85, T = 50, P = 100, solidsolution_type = 'cpx') >>> ss.logK 18.81769186 >>> ss.Rxn {'type': 'cpx', 'name': 'Di85Hed15', 'formula': 'Ca1.00Mg0.85Fe0.15Si2O6', 'MW': 221.2817, 'min': ['Ca1.00Mg0.85Fe0.15Si2O6', ' R&H95, Stef2001', -711392.7391800061, nan, 37.032850185405415, 66.369, 106.66383843212236, -0.019588551625239002, -16326.4818355, -1052.9517208413004, 5.714746653919693e-06], 'spec': ['H+', 'Ca++', 'Mg++', 'Fe++', 'SiO2(aq)', 'H2O'], 'coeff': [-4, 1.0, 0.85, 0.15, 2, 2], 'nSpec': 6, 'V': 66.369, 'source': ' R&H95, Stef2001', 'elements': ['1.0000', 'Ca', '0.8500', 'Mg', '0.1500', 'Fe', '2.0000', 'Si', '6.0000', 'O']} >>> ss = solidsolution_thermo(X = 0.5, T = 50, P = 100, solidsolution_type = 'plagioclase') >>> ss.logK 12.62488394 >>> ss.Rxn {'type': 'plag', 'name': 'An50', 'formula': 'Ca0.50Na0.50Al1.50Si2.50O8', 'MW': 270.215145, 'min': ['Ca0.50Na0.50Al1.50Si2.50O8', ' R&H95, Stef2001', -919965.222804649, nan, 47.33081919017781, 100.43, 131.53908221797, -0.022148661567686426, 32265.761759082205, -1316.0836615678777, 7.719885114722753e-06], 'spec': ['H+', 'Al+++', 'Na+', 'Ca++', 'SiO2(aq)', 'H2O'], 'coeff': [-6.0, 1.5, 0.5, 0.5, 2.5, 3.0], 'nSpec': 6, 'V': 100.43, 'source': ' R&H95, Stef2001', 'elements': ['0.5000', 'Ca', '0.5000', 'Na', '1.5000', 'Al', '2.5000', 'Si', '8.0000', 'O']} .. py:attribute:: kwargs .. py:method:: __calc__(**kwargs) .. py:method:: __checker__(**kwargs) .. py:method:: calclogKAnAb(XAn, TC, P, dbaccessdic, rhoEG, Dielec_method=None) This function calculates thermodynamic properties of solid solution of Plagioclase minerals :param XAn: :type XAn: volume fraction of Anorthite :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param dbaccessdic: :type dbaccessdic: dictionary of species from direct-access database :param Dielec_method: dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logKplag** (*logarithmic K values*) * **Rxn** (*dictionary of reaction thermodynamic properties*) Usage ------- The general usage of calclogKAnAb without the optional argument is as follows: (1) Not on steam saturation curve: [logK, Rxn] = calclogKAnAb(XAn, TC, P, dbaccessdic), where T is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, Rxn] = calclogKAnAb(XAn, TC, 'T', dbaccessdic), where T is temperature in celsius, followed with a quoted char 'T' [logK, Rxn] = calclogKAnAb(XAn, P, 'P', dbaccessdic), where P is pressure in bar, followed with a quoted char 'P'. (3) Meanwhile, usage with any specific dielectric constant method ('FGL97') for condition not on steam saturation curve is as follows. Default method is 'JN91' [logK, Rxn] = calclogKAnAb(XAn, TC, P, dbaccessdic, Dielec_method = 'FGL97') .. py:method:: calclogKFoFa(XFo, TC, P, dbaccessdic, rhoEG, Dielec_method=None) This function calculates thermodynamic properties of solid solution of olivine minerals :param XFo: :type XFo: volume fraction of Forsterite :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param dbaccessdic: :type dbaccessdic: dictionary of species from direct-access database :param Dielec_method: dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logK_ol** (*logarithmic K values*) * **Rxn** (*dictionary of reaction thermodynamic properties*) Usage ------- The general usage of calclogKFoFa without the optional argument is as follows: (1) Not on steam saturation curve: [logK, Rxn] = calclogKFoFa(XFo, TC, P, dbaccessdic), where T is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, Rxn] = calclogKFoFa(XFo, TC, 'T', dbaccessdic), where T is temperature in celsius, followed with a quoted char 'T' [logK, Rxn] = calclogKFoFa(XFo, P, 'P', dbaccessdic), where P is pressure in bar, followed with a quoted char 'P'. (3) Meanwhile, usage with any specific dielectric constant method ('FGL97') for condition not on steam saturation curve is as follows. Default method is 'JN91' [logK, Rxn] = calclogKFoFa(XFo, TC, P, dbaccessdic, Dielec_method = 'FGL97') .. py:method:: calclogKEnFe(XEn, TC, P, dbaccessdic, rhoEG, Dielec_method=None) This function calculates thermodynamic properties of solid solution of pyroxene minerals :param XEn: :type XEn: volume fraction of Enstatite :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param dbaccessdic: :type dbaccessdic: dictionary of species from direct-access database :param Dielec_method: dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logK_opx** (*logarithmic K values*) * **Rxn** (*dictionary of reaction thermodynamic properties*) Usage ------- The general usage of calclogKEnFe without the optional argument is as follows: (1) Not on steam saturation curve: [logK, Rxn] = calclogKEnFe(XEn, TC, P, dbaccessdic), where T is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, Rxn] = calclogKEnFe(XEn, TC, 'T', dbaccessdic), where T is temperature in celsius, followed with a quoted char 'T' [logK, Rxn] = calclogKEnFe(XEn, P, 'P', dbaccessdic), where P is pressure in bar, followed with a quoted char 'P'. (3) Meanwhile, usage with any specific dielectric constant method ('FGL97') for condition not on steam saturation curve is as follows. Default method is 'JN91' [logK, Rxn] = calclogKEnFe(XEn, TC, P, dbaccessdic, Dielec_method = 'FGL97') .. py:method:: calclogKDiHedEnFe(nCa, XMg, TC, P, dbaccessdic, rhoEG, Dielec_method=None) This function calculates thermodynamic properties of solid solution of clinopyroxene minerals (Di, Hed, En and Fe) :param nCa: and must be greater than zero :type nCa: number of moles of Ca in formula unit (=1 for Di, Hed) :param XMg: XMg = (nMg/(nFe + nMg)) :type XMg: mole fraction of Mg :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param dbaccessdic: :type dbaccessdic: dictionary of species from direct-access database :param Dielec_method: dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logK_cpx** (*logarithmic K values*) * **Rxn** (*dictionary of reaction thermodynamic properties*) Usage ------- The general usage of calclogKDiHedEnFe without the optional argument is as follows: (1) Not on steam saturation curve: [logK, Rxn] = calclogKDiHedEnFe(nCa, XMg, TC, P, dbaccessdic), where T is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, Rxn] = calclogKDiHedEnFe(nCa, XMg, TC, 'T', dbaccessdic), where T is temperature in celsius, followed with a quoted char 'T' [logK, Rxn] = calclogKDiHedEnFe(nCa, XMg, P, 'P', dbaccessdic), where P is pressure in bar, followed with a quoted char 'P'. (3) Meanwhile, usage with any specific dielectric constant method ('FGL97') for condition not on steam saturation curve is as follows. Default method is 'JN91' [logK, Rxn] = calclogKDiHedEnFe(nCa, XMg, TC, P, dbaccessdic, Dielec_method = 'FGL97') .. py:method:: calclogKAbOr(XAb, TC, P, dbaccessdic, rhoEG, Dielec_method=None) This function calculates thermodynamic properties of solid solution of Alkaline-Feldspar minerals :param XAb: :type XAb: volume fraction of Albite :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param dbaccessdic: :type dbaccessdic: dictionary of species from direct-access database :param Dielec_method: dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logKalkfeld** (*logarithmic K values*) * **Rxn** (*dictionary of reaction thermodynamic properties*) Usage ------- The general usage of calclogKAnAb without the optional argument is as follows: (1) Not on steam saturation curve: [logK, Rxn] = calclogKAbOr(XAb, TC, P, dbaccessdic), where T is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, Rxn] = calclogKAbOr(XAb, TC, 'T', dbaccessdic), where T is temperature in celsius, followed with a quoted char 'T' [logK, Rxn] = calclogKAbOr(XAb, P, 'P', dbaccessdic), where P is pressure in bar, followed with a quoted char 'P'. (3) Meanwhile, usage with any specific dielectric constant method ('FGL97') for condition not on steam saturation curve is as follows. Default method is 'JN91' [logK, Rxn] = calclogKAbOr(XAb, TC, P, dbaccessdic, Dielec_method = 'FGL97') .. py:method:: calclogKBiotite(XPh, TC, P, dbaccessdic, rhoEG, Dielec_method=None) This function calculates thermodynamic properties of solid solution of Biotite minerals :param XAb: :type XAb: volume fraction of Albite :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param dbaccessdic: :type dbaccessdic: dictionary of species from direct-access database :param Dielec_method: dielectric constant (optional), if not specified default - 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logKalkfeld** (*logarithmic K values*) * **Rxn** (*dictionary of reaction thermodynamic properties*) Usage ------- The general usage of calclogKAnAb without the optional argument is as follows: (1) Not on steam saturation curve: [logK, Rxn] = calclogKAbOr(XAb, TC, P, dbaccessdic), where T is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, Rxn] = calclogKAbOr(XAb, TC, 'T', dbaccessdic), where T is temperature in celsius, followed with a quoted char 'T' [logK, Rxn] = calclogKAbOr(XAb, P, 'P', dbaccessdic), where P is pressure in bar, followed with a quoted char 'P'. (3) Meanwhile, usage with any specific dielectric constant method ('FGL97') for condition not on steam saturation curve is as follows. Default method is 'JN91' [logK, Rxn] = calclogKAbOr(XAb, TC, P, dbaccessdic, Dielec_method = 'FGL97')