pygcc.pygcc_utils ================= .. py:module:: pygcc.pygcc_utils .. 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.pygcc_utils.eps pygcc.pygcc_utils.J_to_cal pygcc.pygcc_utils.IAPWS95_COEFFS Classes ------- .. autoapisummary:: pygcc.pygcc_utils.calcRxnlogK pygcc.pygcc_utils.write_database Functions --------- .. autoapisummary:: pygcc.pygcc_utils.calc_elem_count_molewt pygcc.pygcc_utils.importconfile pygcc.pygcc_utils.var_name pygcc.pygcc_utils.feval pygcc.pygcc_utils.roundup_tenth pygcc.pygcc_utils.roundup_hundredth pygcc.pygcc_utils.read_specific_lines pygcc.pygcc_utils.denormalize_phreeqc_species_charge pygcc.pygcc_utils.normalize_species_charges pygcc.pygcc_utils.normalize_phreeqc_species_charge pygcc.pygcc_utils.contains_missing_species pygcc.pygcc_utils.build_side_regex pygcc.pygcc_utils.derivative pygcc.pygcc_utils.info pygcc.pygcc_utils.drummondgamma pygcc.pygcc_utils.Henry_duan_sun pygcc.pygcc_utils.co2fugacity pygcc.pygcc_utils.gamma_correlation pygcc.pygcc_utils.Helgeson_activity pygcc.pygcc_utils.aw_correlation pygcc.pygcc_utils.outputfmt pygcc.pygcc_utils.main_function_name Module Contents --------------- .. py:data:: eps :value: 2.220446049250313e-16 .. py:data:: J_to_cal :value: 4.184 .. py:data:: IAPWS95_COEFFS .. py:function:: calc_elem_count_molewt(formula, **kwargs) This function calculates the molecular mass and the elemental composition of a substance given by its chemical formula. This is modified from https://github.com/cgohlke/molmass/blob/master/molmass/molmass.py :param formula: chemical formula :type formula: string :param Elementdic: dictionary, containing the atomic mass of element database :type Elementdic: dict :returns: * **elements** (*dict*) -- dictionary of elemental composition and their respective number of atoms * **molewt** (*float*) -- Calculated Molecular Weights [g/mol] Usage: ---------- [elements, molewt] = calc_elem_count_molewt(formula) Examples of valid formulas are "H2O", "[2H]2O", "CH3COOH", "EtOH", "CuSO4:5H2O", "(COOH)2", "AgCuRu4(H)2[CO]12{PPh3}2", "CGCGAATTCGCG", and, "MDRGEQGLLK" . .. rubric:: Examples >>> elements, molewt = calc_elem_count_molewt("CuSO4:5H2O") >>> elements, molewt {'O': 9.0, 'H': 10.0, 'S': 1, 'Cu': 1}, 249.6773 .. py:function:: importconfile(filename, *Rows) This function imports numeric data from a text file as column vectors. [Var1, Var2] = importconfile(filename) Reads data from text file filename for the default selection. [Var1, Var2] = importconfile(filename, StartRow, EndRow) Reads data from rows StartRow through EndRow of text file filename. Example: [Var1, Var2] = importconfile('100bar.con', 6, 13); .. py:function:: var_name(var) .. py:function:: feval(funcName, *args) .. py:function:: roundup_tenth(x) .. py:function:: roundup_hundredth(x) .. py:function:: read_specific_lines(file, lines_to_read) .. py:function:: denormalize_phreeqc_species_charge(j) Helper function to reverse the normalization of PHREEQC species charge notation. - Convert expanded charges back into compact notation (e.g., 'Ca++' -> 'Ca+2'). - Remove '(aq)' from neutral species. - Keep integers or floats unchanged. Examples: 'Ca++' -> 'Ca+2' 'Fe---' -> 'Fe-3' 'H2O(aq)' -> 'H2O' 'CH4(aq)' -> 'CH4' '-1' -> '-1' '3.5' -> '3.5' .. py:function:: normalize_species_charges(species) .. py:function:: normalize_phreeqc_species_charge(j) Helper function to normalize PHREEQC species charge notation. - Replace charge notation like +2, -3, +4 etc. with expanded '+' or '-' symbols. - Skip cases where the entire string is just an integer (e.g., '1', '-1', '+2'). - If the species has no explicit charge, append '(aq)'. Examples: 'Ca+2' -> 'Ca++' 'Fe-3' -> 'Fe---' '-1' -> '-1' 'H2O' -> 'H2O(aq)' 'CH4' -> 'CH4(aq)' .. py:function:: contains_missing_species(line, missing_spx) .. py:function:: build_side_regex(side) Build a regex to match one side of a chemical reaction. - Optional numeric coefficients (integer, decimal, scientific notation) - Flexible spacing around '+' signs - Species with charges, e.g., CO3-2, Eu+3, H+, e- .. py:function:: derivative(f, a, method='central', h=0.001) Compute the derivative of f, f'(a) with step size h. :param f: Vectorized function of one variable :type f: function :param a: Compute derivative at x = a :type a: number :param method: Difference formula: 'forward', 'backward' or 'central' :type method: string :param h: Step size in difference formula :type h: number :returns: Difference formula: central: f(a + h) - f(a - h))/2h forward: f(a + h) - f(a))/h backward: f(a) - f(a-h))/h :rtype: float .. py:function:: info(name, dic) This function checks for naming convention of the species in the direct-access or source database :param name: species name :type name: string :param dic: dictionary of species from direct-access or source database :type dic: dict :returns: **lst** -- resulting search of all species with the input name :rtype: list .. rubric:: Examples >>> from pygcc.pygcc_utils import db_reader >>> ps = db_reader() # utilizes the default direct-access database, speq21 >>> info('ss_', ps.dbaccessdic) ['ss_Anorthite', 'ss_Albite_high', 'ss_K-feldspar', 'ss_Ferrosilite', 'ss_Enstatite', 'ss_Clinoenstatite', 'ss_Hedenbergite', 'ss_Diopside', 'ss_Forsterite', 'ss_Fayalite', 'ss_Annite', 'ss_Phlogopite', 'ss_Anorthite', 'ss_Albite_high', 'ss_K-feldspar', 'ss_Ferrosilite', 'ss_Enstatite', 'ss_Clinoenstatite', 'ss_Hedenbergite', 'ss_Diopside', 'ss_Forsterite', 'ss_Fayalite', 'ss_Annite', 'ss_Phlogopite'] .. py:function:: drummondgamma(TK, I) This function models solubility of CO2 gas in brine using Drummond equation :param TK: Temperature [K] :type TK: float, vector :param I: Ionic strength :type I: float, vector Returns ---------- log10_gamma : float, vector co2 aqueous activity coefficients in log10 .. rubric:: Examples >>> log10_gamma = drummondgamma( 500, 0.5) >>> log10_gamma 0.0781512920184902 .. py:function:: Henry_duan_sun(TK, P, I) This function evaluates the solubility of CO2 gase in brine using Duan_Sun Formulation :param TK: Temperature [K] :type TK: float, vector :param P: Pressure [bar] :type P: float, vector :param I: Ionic strength :type I: float :returns: * **log10_co2_gamma** (*float, vector*) -- co2 aqueous activity coefficients in log10 * **mco2** (*float, vector*) -- co2 aqueous molalities Usage ---------- log10_co2_gamma, mco2 = Henry_duan_sun( TK, P, I) .. rubric:: Examples >>> log10_co2_gamma, mco2 = Henry_duan_sun( 500, 250, 0.5) >>> log10_co2_gamma, mco2 0.06202532, 1.6062557 .. rubric:: References (1) Duan, Z. and Sun, R., 2003. An improved model calculating CO2 solubility in pure water and aqueous NaCl solutions from 273 to 533 K and from 0 to 2000 bar. Chemical geology, 193(3-4), pp.257-271. .. py:function:: co2fugacity(TK, P, poy=True) This function computes the fugacity and density of CO2 by Duan and Sun 2003 Also Calculate the Saturation Index SI and partial pressure of CO2(g) at any given T, P using the Duan equation of state. A Poynting correction factor is also applied. :param TK: Temperature [K] :type TK: float, vector :param P: pressure [bar] :type P: float, vector :param poy: Option to apply a Poynting correction factor [True is the default] :type poy: bool :returns: * **fCO2** (*float, vector*) -- co2 fugacity [bar] * **denCO2** (*float, vector*) -- co2 density [g/cm3] * **SI** (*float, vector*) -- co2 saturation index * **pCO2** (*float, vector*) -- co2 partial pressure [bar] Usage ---------- [fCO2, denCO2, SI, pCO2] = co2fugacity( TK, P) .. rubric:: Examples >>> fCO2, denCO2, SI, pCO2 = co2fugacity( np.array([400,420]), np.array([150, 200])) >>> fCO2 array([114.50387423, 150.87672517]) >>> denCO2 array([2.6741242 , 3.33134397]) >>> SI array([1.73291553, 1.84597591]) >>> pCO2 array([54.78127571, 71.0710152]) .. py:function:: gamma_correlation(TC, P, method=None) This function calculates the CO2 activity correlation coefficients at given temperature T and pressure P :param TC: Temperature [°C] :type TC: float, vector :param P: Pressure [bar] :type P: float, vector :param method: specify the activity model [``'Duan_Sun'`` or ``'Drummond'``] :type method: string Returns ---------- cco2 : float, vector co2 correlation coefficients Usage ---------- cco2 = gamma_correlation( TC, P) .. rubric:: Examples >>> log10_co2_gamma, mco2 = Henry_duan_sun( 500, 250, 0.5) >>> log10_co2_gamma, mco2 0.06202532, 1.6062557 .. rubric:: References (1) Segal Edward Drummond, 1981, Boiling and Mixing of Hydrothermal Fluids: Chemical Effects on Mineral Precipitation, page 19 (2) Wolery, T. J., Lawrence Livermore National Laboratory, United States Dept. of Energy, 1992. EQ3/6: A software package for geochemical modeling of aqueous systems: package overview and installation guide (version 7.0) .. py:function:: Helgeson_activity(TC, P, I, Dielec_method=None, **rhoEDB) This function calculates the solute activity coefficient, solvent osmotic coefficient, and solvent activity at given temperature and pressure using equations 298, 190 and 106 in Helgeson, Kirkham and Flowers, 1981, A.J.S. p.1249-1516 Parameters: ---------- TC : Temperature [°C] P : pressure [bar] I : ionic strength Dielec_method : specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate dielectric constant (optional), if not specified default - 'JN91' rhoEDB : dictionary of water properties like density (rho), dielectric factor (E) and Debye–Hückel coefficients (optional) Returns: ---------- aw : solvent activity phi : solvent osmotic coefficient mean_act : solute activity coefficient Usage ------- [aw, phi, mean_act] = Helgeson_activity( TC, P, I) .. py:function:: aw_correlation(TC, P, Dielec_method=None, **rhoEDB) Calculates the water activity correlation coefficients at given temperature and pressure Parameters: ---------- TC : Temperature [°C] P : pressure [bar] Dielec_method : specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate dielectric constant (optional), if not specified default - 'JN91' rhoEDB : dictionary of water properties like density (rho), dielectric factor (E) and Debye-Hückel coefficients (optional) Returns: ---------- ch20 : water activity correlation coefficients Usage ------- [ch20] = aw_correlation( TC, P) .. py:class:: calcRxnlogK(**kwargs) This class implemetation calculates logK values for any reaction with the option for extrapolation where rho < 350kg/m3 :param T: Temperature [°C] :type T: float, vector :param P: Pressure [bar] :type P: float, vector :param Specie: specify the species for logK calculation, either the Product species of any reaction or solid solutions or clay like 'AnAb' or 'AbOr' or 'FoFa' or 'EnFe' or 'DiHedEnFe' or 'clay' :type Specie: string :param Specie_class: specify the class of species like 'aqueous', 'minerals', 'liquids', or 'gases' :type Specie_class: string, optional :param elem: list containing nine or ten parameters with clay names and elements compositions with the following format ['Montmorillonite_Lc_MgK', 'Si', 'Al', 'FeIII', 'FeII', 'Mg', 'K', 'Na', 'Ca', 'Li', 'H3O'] :type elem: list :param dbaccessdic: direct-acess database dictionary :type dbaccessdic: dict :param rhoEGextrap: dictionary of water properties like density (rho), dielectric factor (E) and Gibbs Energy for density region 350-550kg/m3 :type rhoEGextrap: dict :param group: specify the structural layering of the phyllosilicate, for layers composed of ``1 tetrahedral + 1 octahedral sheet (1:1 layer)`` - specify ``'7A'``, ``2 tetrahedral + 1 octahedral sheet (2:1 layer)`` - specify ``'10A'``, or the latter with ``a brucitic sheet in the interlayer (2:1:1 layer)`` - specify ``'14A'`` (optional), if not specified, default is '10A' for smectites, micas, et cetera :type group: string :param ClayMintype: specify either 'Smectite' or 'Chlorite' or 'Mica' as the clay type, if not specified default - 'Smectites' :type ClayMintype: string :param X: volume fractions of any (Anorthite, Albite, Forsterite, Enstatite) or mole fraction of Mg :type X: float :param cpx_Ca: number of moles of Ca in formula unit (=1 for Di, Hed), must be greater than zero :type cpx_Ca: float :param sourcedic: source database reactions dictionary :type sourcedic: dict :param specielist: source database species grouped into categories [element, basis, redox, aqueous, minerals, gases, oxides] :type specielist: list of list, optional :param Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate dielectric constant, default is 'JN91' :type Dielec_method: string :param heatcap_method: specify either 'SUPCRT' or 'Berman88' or 'HP11' or 'HF76' as the method to calculate thermodynamic properties of any mineral or gas, default is 'SUPCRT' :type heatcap_method: string :param ThermoInUnit: 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' :type ThermoInUnit: string :param Al_Si: 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' :type Al_Si: string :param Int_Mg_fract: specify the fraction of Mg to partition into Interlayer sheet and the remainder will be partitioned into Octahedral sheet, if not specified, default is 1 :type Int_Mg_fract: string :param Int_Li_fract: specify the fraction of Li to partition into Interlayer sheet and the remainder will be partitioned into Octahedral sheet, if not specified, default is 1 :type Int_Li_fract: string :param heatcap_approx: specify either 'Maier-Kelley' or 'constant' as the approximation method for clay minerals' specific heat capacity calculation, default is 'constant', based on ClayTherm's definition for specific cations (Octahedral sites: Li+, Mn2+, Cr3+, Ni2+, Co2+, Zn2+; Interlayer sites: Cs+, Rb+, Li+, Ba2+, Sr2+, Mg2+, Cu2+, Co2+, Zn2+, H3O+) :type heatcap_approx: string :param sourceformat: specify the source database format, either 'GWB', 'EQ36' or 'PHREEQC' :type sourceformat: string :param densityextrap: specify the extrapolation option for density-logK, 'Yes'/True or 'No'/False :type densityextrap: float, vector :returns: * **log_K** (*float, vector*) -- logarithmic K value(s) * **dGrxn** (*float, vector*) -- Total reaction Gibbs energy [cal/mol] Usage ---------- The general usage of water_dielec is as follows: (1) For water dielectric properties at any Temperature and Pressure: calclogK = calcRxnlogK(T = T, P = P, Dielec_method = 'JN91', **kwargs), where T is temperature in celsius and P is pressure in bar (2) For water dielectric properties at any Temperature and density : calclogK = calcRxnlogK(T = T, rho = rho, Dielec_method = 'JN91', **kwargs), 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: calclogK = calcRxnlogK(T = T, P = 'T', Dielec_method = 'JN91', **kwargs), where T is temperature in celsius, followed with a quoted character 'T' to reflect steam saturation pressure calclogK = calcRxnlogK(P = P, T = 'P', Dielec_method = 'JN91', **kwargs), where P is pressure in bar, followed with a quoted character 'P' to reflect steam saturation temperature .. rubric:: Examples >>> ps = db_reader(sourcedb = './default_db/thermo.com.dat', sourceformat = 'gwb', dbaccess = './default_db/speq21.dat') >>> calclogK = calcRxnlogK(T = 100, P = 50, Specie = 'H2S(aq)', Specie_class = 'aqueous', dbaccessdic = ps.dbaccessdic, sourcedic = ps.sourcedic, specielist = ps.specielist) >>> calclogK.logK, calclogK.dGrxn -6.4713, 11049.3168 .. py:attribute:: kwargs .. py:method:: __calc__(**kwargs) .. py:method:: AllRxnslogK(TC, P, rhoEG) This function calculates logK values of all reactions including solid-solution and clay minerals :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logK** (*logarithmic K value(s)*) * **dGrxn** (*Total reaction Gibbs energy [cal/mol]*) * **dGP** (*Product specie Gibbs energy [cal/mol]*) * **dGRs** (*Reactant species Gibbs energy [cal/mol]*) * **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 calcRxnlogK without the optional argument is as follows: (1) Not on steam saturation curve: [logK, dGrxn, dGP, dGRs] = AllRxnslogK(TC, P, Prod, dbaccessdic, sourcedic, specielist), where TC is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, dGrxn, dGP, dGRs] = AllRxnslogK(TC, 'T', Prod, dbaccessdic, sourcedic, specielist), where TC is temperature in celsius, followed with a quoted char 'T' [logK, dGrxn, dGP, dGRs] = AllRxnslogK(P, 'P', Prod, dbaccessdic, sourcedic, specielist), 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, dGrxn, dGP, dGRs] = AllRxnslogK( TC, P, Prod, dbaccessdic, sourcedic, specielist, Dielec_method = 'FGL97') .. py:method:: RxnlogK(TC, P, Prod, rhoEG) This function calculates logK values of any reaction :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param Prod: :type Prod: Product species of the reaction :param dbaccessdic: :type dbaccessdic: direct-acess database dictionary :param sourcedic: :type sourcedic: source database reactions dictionary :param specielist: [element, basis, redox, aqueous, minerals, gases, oxides] :type specielist: source database species grouped into :param Dielec_method: to calculate dielectric constant, default is 'JN91' :type Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method :param sourceformat: :type sourceformat: source database format, either 'GWB' or 'EQ36', default is 'GWB' :param heatcap_method: of any mineral, default is 'SUPCRT' :type heatcap_method: specify either 'SUPCRT' or 'Berman88' or 'HP11' or 'HF76' as the method to calculate thermodynamic properties :param rhoEG: dielectric factor (E) and Gibbs Energy (optional) :type rhoEG: dictionary of water properties like density (rho), :returns: * **logK** (*logarithmic K value(s)*) * **dGrxn** (*Total reaction Gibbs energy [cal/mol]*) * **dGP** (*Product specie Gibbs energy [cal/mol]*) * **dGRs** (*Reactant species Gibbs energy [cal/mol]*) Usage ------- The general usage of calcRxnlogK without the optional argument is as follows: (1) Not on steam saturation curve: [logK, dGrxn, dGP, dGRs] = RxnlogK(TC, P, Prod, dbaccessdic, sourcedic, specielist), where TC is temperature in celsius and P is pressure in bar; (2) On steam saturation curve: [logK, dGrxn, dGP, dGRs] = RxnlogK(TC, 'T', Prod, dbaccessdic, sourcedic, specielist), where TC is temperature in celsius, followed with a quoted char 'T' [logK, dGrxn, dGP, dGRs] = RxnlogK(P, 'P', Prod, dbaccessdic, sourcedic, specielist), 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, dGrxn, dGP, dGRs] = RxnlogK( TC, P, Prod, dbaccessdic, sourcedic, specielist, Dielec_method = 'FGL97') .. py:method:: densitylogKextrap(TC, P, rhoEGextrap) This function calculates logK values extrapolation for conditions where rho < 350kg/m3 :param TC: :type TC: temperature [°C] :param P: :type P: pressure [bar] :param rhoEGextrap: Gibbs Energy for density region 350-550kg/m3 :type rhoEGextrap: dictionary of water properties like density (rho), dielectric factor (E) and :returns: **logK** :rtype: extrapolated logarithmic K value(s) Usage ------- [logK] = densitylogKextrap(TC, P, 'H2S(aq)', dbaccessdic = dbaccessdic, sourcedic = sourcedic, specielist = specielist), .. py:function:: outputfmt(fid, logK, Rxn, *T, dataset=None, logK_form=None) This function writes logK and Rxn data to any file using GWB, EQ36, PHREEQC, Pflotran and ToughReact format :param fid: file ID :type fid: string :param logK: logarithmic K value(s) :type logK: float, vector :param Rxn: dictionary of reaction thermodynamic properties :type Rxn: dict :param T: Temperature value(s), optional, required when 'polycoeffs' is specified for logK_form :type T: float, vector :param dataset: specify the dataset format, either 'GWB', 'EQ36', 'PHREEQC', 'Pflotran' or 'ToughReact' :type dataset: string :param logK_form: specify the format of logK either as a set of eight values one for each of the dataset's principal temperatures, or blocks of polynomial coefficients, [values, polycoeffs], default is 'a set of eight values' (optional) :type logK_form: string :rtype: Output data to the file with filename described in fid with any format mentioned above. .. rubric:: Examples >>> fid = open('./logK_details.txt', 'w') >>> logKRxn = calcRxnlogK(T = 100, P = 'T', X = 0.634, Specie = 'Plagioclase', densityextrap = True) >>> # output in EQ36 format >>> outputfmt(fid, logKRxn.logK, logKRxn.Rxn, dataset = 'EQ36') >>> fid.close() .. py:class:: write_database(**kwargs) Class to write the new database for either GWB, EQ3/6, ToughReact, Pflotran or PHREEQC into a new folder called "output" \n :param T: Temperature [°C] \n :type T: float, vector :param P: Pressure [bar] \n :type P: float, vector :param cpx_Ca: number of moles of Ca in solid solution of clinopyroxene (optional) if it is ommitted solid solution of clinopyroxene will not be included, 0 < nCa >=1 \n :type cpx_Ca: string :param solid_solution: specify the inclusion of solid-solution [Yes/True or No/False], default is 'No' \n :type solid_solution: string, bool :param clay_thermo: specify the inclusion of clay thermodynamic properties [Yes/True or No/False], default is 'No' \n :type clay_thermo: float, vector :param logK_form: specify the format of logK either as a set of eight values one for each of the dataset’s principal temperatures, or blocks of polynomial coefficients, [values, polycoeffs] default is 'a set of eight values' \n :type logK_form: float, vector :param densityextrap: specify the utilization of density extrapolation [Yes/True or No/False], default is 'Yes' \n :type densityextrap: float, vector :param dbaccess: direct-access database filename and location (optional) \n :type dbaccess: string :param dbBerman_dir: filename and location of the Berman mineral database (optional) \n :type dbBerman_dir: string :param dbHP_dir: filename and location of the supcrtbl mineral and gas database, optional :type dbHP_dir: string :param dbaccessformat: specify the direct-access/sequential-access database format, either 'speq' or 'supcrtbl', default is 'speq' :type dbaccessformat: string, optional :param sourcedb: source database filename and location (optional) \n :type sourcedb: string :param sourceformat: source database format, either 'GWB', 'EQ36' or 'PHREEQC', default is 'GWB' :type sourceformat: string :param sourcedb_codecs: specify the name of the encoding used to decode or encode the sourcedb file, optional :type sourcedb_codecs: string :param objdb: new database filename and location (optional) \n :type objdb: string :param co2actmodel: co2 activity model equation [Duan_Sun or Drummond] (optional), if not specified, default is 'Drummond' \n :type co2actmodel: string :param Dielec_method: specify either 'FGL97' or 'JN91' or 'DEW' as the method to calculate dielectric constant, default is 'JN91' (optional) \n :type Dielec_method: string :param heatcap_method: specify either 'SUPCRT' or 'Berman88' or 'HP11' or 'HF76' the method to calculate thermodynamic properties of any mineral, default is 'SUPCRT' \n :type heatcap_method: string :param ThermoInUnit: 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' :type ThermoInUnit: string :param dataset: specify the dataset format, either 'GWB', 'EQ36', 'PHREEQC', 'Pflotran' or 'ToughReact', default is old GWB database ['GWB'] (optional) \n :type dataset: string :param print_msg: print debug message [True or False], default is False \n :type print_msg: string, bool :rtype: Output the new database to an ASCII file with filename described in 'objdb' if specified. Usage ------- With any Temperature and Pressure: (1) General format with default dielectric constant and CO2 activity model and exclusions of solid solutions for GWB \n write_database(T = T, P = P, cpx_Ca = nCa, dataset = 'GWB', sourceformat = 'GWB') \n (2) Inclusion of solid solutions and exclusion of solid solution of clinopyroxene and clay thermo \n write_database(T = T, P = P, solid_solution = 'Yes', clay_thermo = 'Yes', dataset = 'GWB', sourceformat = 'GWB') \n (3) Inclusion of all solid solutions and clay thermo with \\emph{'Duan_Sun'} CO2 activity model and 'FGL97' dielectric constant calculation \n write_database(T = T, P = P, cpx_Ca = nCa, solid_solution = 'Yes', clay_thermo = 'Yes', co2actmodel = 'Duan_Sun', Dielec_method = 'FGL97', dataset = 'GWB', sourceformat = 'GWB') \n With any Temperature or Pressure on the steam saturation curve: (4) General format with default dielectric constant and CO2 activity model \n write_database(T = T, P = 'T', cpx_Ca = nCa, dataset = 'GWB', sourceformat = 'GWB'), \n where T is temperature in celsius and P is assigned a quoted character 'T' to reflect steam saturation pressure \n write_database(P = P, T = 'P', cpx_Ca = nCa, dataset = 'GWB', sourceformat = 'GWB'), \n where P is pressure in bar and T is assigned a quoted character 'P' to reflect steam saturation temperature .. rubric:: Examples >>> write_database(T=np.array([ 0.010, 25., 60., 100., 150., 200., 250., 300.]), P=200, cpx_Ca = 0.5, sourceformat = 'GWB', solid_solution = True, clay_thermo = True, dataset = 'GWB') \n >>> write_database(T=np.array([ 0.010, 25., 60., 100., 150., 200., 250., 300.]), P=200, cpx_Ca = 0.5, sourceformat = 'GWB', solid_solution = True, clay_thermo = True, dataset = 'Pflotran', sourcedb = './pygcc/default_db/thermo.2021.dat') \n .. py:attribute:: kwargs .. py:attribute:: T :value: None .. py:attribute:: P :value: None .. py:attribute:: dbaccessformat .. py:attribute:: dbBerman_dir .. py:attribute:: dbHP_dir .. py:attribute:: dataset .. py:attribute:: objdb .. py:attribute:: co2actmodel .. py:attribute:: ThermoInUnit .. py:attribute:: sourcedb_codecs .. py:attribute:: Dielec_method .. py:attribute:: cpx_Ca .. py:attribute:: logK_form .. py:attribute:: heatcap_method .. py:attribute:: dbaccess .. py:attribute:: sourceformat .. py:method:: __calc__(**kwargs) .. py:method:: _resolve_sourcedb(sourceformat, sourcedb) Resolve default source database based on format and user input. .. py:method:: _flag_to_yesno(value, default='No') Convert boolean/None/custom flag to Yes/No string. .. py:method:: _process_PT() Process temperature and pressure arrays safely. Accepted string sentinels: P == 'T' -> P is computed along the saturation curve for the given T T == 'P' -> T is computed along the saturation curve for the given P Any other string raises ValueError. All `==` against the sentinels is gated on ``isinstance(..., str)`` so an ndarray is never compared to a string (which produces NumPy FutureWarning / breaks under NumPy 2.x). Saturation properties are computed directly via ``Dummy().vcalcsatpropT`` / ``Dummy().vcalcsatpropP`` instead of going through a full ``iapws95(...)`` call. .. py:method:: _expand_temperature(T) Expand temperature input into proper array. .. py:method:: _build_message() Construct final success message. .. py:method:: write_GWBdb(T, P) This function writes the new GWB database into a new folder called "output" \n :param T: :type T: temperature [°C] \n :param P: :type P: pressure [bar] \n :rtype: Outputs the new database to an ASCII file with filename described in 'objdb'. \n Usage ------- Example: (1) General format with default dielectric constant and CO2 activity model and exclusions of solid solutions \n write_GWBdb(T, P) \n (2) Inclusion of solid solutions and clay thermo and exclusion of solid solution of clinopyroxene \n write_GWBdb(T, P) \n (3) Inclusion of all solid solutions and clay thermo with \\emph{'Duan_Sun'} CO2 activity model and 'JN91' dielectric constant calculation \n write_GWBdb(T, P) \n .. py:method:: write_EQ36db(T, P) This function writes the new EQ3/6 database into a new folder called "output" \n :param T: :type T: temperature [°C] \n :param P: :type P: pressure [bar] \n :rtype: Outputs the new database to an ASCII file with filename described in 'objdb'. \n Usage ------- Example: (1) General format with default dielectric constant and CO2 activity model and exclusions of solid solutions \n write_EQ36db(T, P ) \n (2) Inclusion of solid solutions and clay thermo and exclusion of solid solution of clinopyroxene \n write_EQ36db(T, P ) \n (3) Inclusion of all solid solutions and clay thermo with \\emph{'Duan_Sun'} CO2 activity model and 'FGL97' dielectric constant calculation \n write_EQ36db(T, P ) \n .. py:method:: write_PHREEQCdb(T, P) This function writes the new PHREEQC database into a new folder called "output" \n :param T: :type T: temperature [°C] \n :param P: :type P: pressure [bar] \n :rtype: Outputs the new database to an ASCII file with filename described in 'objdb'. \n Usage ------- Example: (1) General format with default dielectric constant and CO2 activity model and exclusions of solid solutions \n write_PHREEQCdb(T, P ) \n (2) Inclusion of solid solutions and clay thermo and exclusion of solid solution of clinopyroxene \n write_PHREEQCdb(T, P ) \n (3) Inclusion of all solid solutions and clay thermo with \\emph{'Duan_Sun'} CO2 activity model and 'FGL97' dielectric constant calculation \n write_PHREEQCdb(T, P ) \n .. py:method:: write_pflotrandb(T, P) This function writes the new pflotran database into a new folder called "output" \n :param T: :type T: temperature [°C] \n :param P: :type P: pressure [bar] \n :rtype: Outputs the new database to an ASCII file with filename described in 'objdb'. \n Usage ------- Example: (1) General format with default dielectric constant and CO2 activity model and exclusions of solid solutions \n write_pflotrandb(T, P ) \n (2) Inclusion of solid solutions and clay thermo and exclusion of solid solution of clinopyroxene \n write_pflotrandb(T, P ) \n (3) Inclusion of all solid solutions and clay thermo with \\emph{'Duan_Sun'} CO2 activity model and 'FGL97' dielectric constant calculation \n write_pflotrandb(T, P ) \n .. py:method:: write_ToughReactdb(T, P) This function writes the new ToughReact database into a new folder called "output" \n :param T: :type T: temperature [°C] \n :param P: :type P: pressure [bar] \n :rtype: Outputs the new database to an ASCII file with filename described in 'objdb'. \n Usage ------- Example: (1) General format with default dielectric constant and CO2 activity model and exclusions of solid solutions \n write_ToughReactdb(T, P, dbaccess = 'location', sourcedb = 'location', objdb = 'location', sourceformat = 'GWB') \n (2) Inclusion of solid solutions and clay thermo and exclusion of solid solution of clinopyroxene \n write_ToughReactdb(T, P, solid_solution = 'Yes', clay_thermo = 'Yes', dbaccess = 'location', sourcedb = 'location', objdb = 'location', sourceformat = 'GWB') \n (3) Inclusion of all solid solutions and clay thermo with \\emph{'Duan_Sun'} CO2 activity model and 'FGL97' dielectric constant calculation \n write_ToughReactdb(T, P, nCa, solid_solution = 'Yes', clay_thermo = 'Yes', dbaccess = 'location', sourcedb = 'location', objdb = 'location', co2actmodel = 'Duan_Sun', Dielec_method = 'FGL97', sourceformat = 'GWB') \n .. py:function:: main_function_name(module) Return the names of essential public functions and classes in the given module. Filters out imported third-party names, stdlib helpers, internal constants, and private utility functions so that only the user-facing API is listed.