Gems3k  3.1
GEMS3K standalone solver for geochemical equilibria
 All Classes Files Functions Variables Enumerations Enumerator
Public Member Functions | Protected Member Functions | Protected Attributes
TNode Class Reference

Implements a simple C/C++ interface between GEM IPM and FMT codes. More...

#include <node.h>

Inheritance diagram for TNode:
TNodeArray

List of all members.

Public Member Functions

 TNode ()
 Constructor of the class instance in memory for standalone GEMS3K or coupled program.
virtual ~TNode ()
 destructor
long int GEM_init (const char *ipmfiles_lst_name, const char *dbrfiles_lst_name=0, long int *nodeTypes=0, bool getNodT1=false)
 (1) Initialization of GEM IPM3 data structures in coupled programs that use GEMS3K module.
void GEM_restore_MT (long int &p_NodeHandle, long int &p_NodeStatusCH, double &p_TK, double &p_P, double &p_Vs, double &p_Ms, double *p_bIC, double *p_dul, double *p_dll, double *p_aPH)
 (6) Passes (copies) the GEMS3K input data from the work instance of DATABR structure.
void GEM_from_MT (long int p_NodeHandle, long int p_NodeStatusCH, double p_TK, double p_P, double p_Vs, double p_Ms, double *p_bIC, double *p_dul, double *p_dll, double *p_aPH)
 (8) Loads the GEMS3K input data for a given mass-transport node into the work instance of DATABR structure.
void GEM_from_MT (long int p_NodeHandle, long int p_NodeStatusCH, double p_TK, double p_P, double p_Vs, double p_Ms, double *p_bIC, double *p_dul, double *p_dll, double *p_aPH, double *p_xDC)
 (8a) Loads the GEMS3K input data for a given mass-transport node into the work instance of DATABR structure.
void GEM_from_MT (long int p_NodeHandle, long int p_NodeStatusCH, double p_TK, double p_P, double p_Vs, double p_Ms, double *p_bIC, double *p_dul, double *p_dll, double *p_aPH, double *p_xDC, double *p_gam)
 (8b) Loads the GEMS3K input data for a given mass-transport node into the work instance of DATABR structure.
void GEM_set_MT (double p_Tm, double p_dt)
 (9) Optional, for passing the current mass transport iteration information into the work DATABR structure (e.g.
long int GEM_read_dbr (const char *fname, bool binary_f=false)
 (5) Reads another DBR file (with input system composition, T,P etc.) \ .
long int GEM_run (bool uPrimalSol)
 (2) Main call for GEM IPM calculations using the input bulk composition, temperature, pressure and metastability constraints provided in the work instance of DATABR structure.
double GEM_CalcTime ()
 Returns GEMIPM2 calculation time in seconds elapsed during the last call of GEM_run() - can be used for monitoring the performance of calculations.
long int GEM_Iterations (long int &PrecLoops, long int &NumIterFIA, long int &NumIterIPM)
 To obtain the number of GEM IPM2 iterations performed during the last call of GEM_run() e.g.
void GEM_write_dbr (const char *fname, bool binary_f=false, bool with_comments=true, bool brief_mode=false)
 (3) Writes the contents of the work instance of the DATABR structure into a disk file with path name fname.
void GEM_print_ipm (const char *fname)
 (4) Produces a formatted text file with detailed contents (scalars and arrays) of the GEM IPM work structure.
void GEM_to_MT (long int &p_NodeHandle, long int &p_NodeStatusCH, long int &p_IterDone, double &p_Vs, double &p_Ms, double &p_Gs, double &p_Hs, double &p_IC, double &p_pH, double &p_pe, double &p_Eh, double *p_rMB, double *p_uIC, double *p_xDC, double *p_gam, double *p_xPH, double *p_vPS, double *p_mPS, double *p_bPS, double *p_xPA, double *p_aPH, double *p_bSP)
 (7) Retrieves the GEMIPM2 chemical speciation calculation results from the work DATABR structure instance into memory provided by the mass transport part.
DATACHpCSD () const
 Get the pointer to chemical system definition data structure.
DATABRpCNode () const
 Get pointer to work node data structure usage on the level of TNodearray is not recommended !
double cTC () const
 Get current node Temperature T, Celsius.
double cTK () const
 Get current node Temperature T, Kelvin.
double cP () const
 Get current node Pressure P, Pa.
double cMs () const
 Get current node mass in kg (reactive part)
double cVs () const
 Get current node volume in m3 (reactive part)
void setNodeHandle (long int jj)
 Set current node identification handle to value of.
long int IC_name_to_xCH (const char *Name)
 Returns DCH index of IC given the IC Name string (null-terminated) or -1 if no such name was found in the DATACH IC name list.
long int DC_name_to_xCH (const char *Name)
 Returns DCH index of DC given the DC Name string or -1 if no such name was found in the DATACH DC name list.
long int Ph_name_to_xCH (const char *Name)
 Returns DCH index of Phase given the Phase Name string or -1 if no such name was found in the DATACH Phase name list.
long int IC_name_to_xDB (const char *Name)
 Returns DBR index of IC given the IC Name string or -1 if no such name was found in the DATACH IC name list.
long int DC_name_to_xDB (const char *Name)
 Returns DBR index of DC given the DC Name string or -1 if no such name was found in the DATACH DC name list.
long int Ph_name_to_xDB (const char *Name)
 Returns DBR index of Phase given the Phase Name string or -1 if no such name was found in the DATACH Phase name list.
long int IC_xCH_to_xDB (const long int xCH)
 Converts the IC DCH index into the IC DBR index or returns -1 if this IC is not used in the data bridge.
long int DC_xCH_to_xDB (const long int xCH)
 Converts the DC DCH index into the DC DBR index or returns -1 if this DC is not used in the data bridge.
long int Ph_xCH_to_xDB (const long int xCH)
 Converts the Phase DCH index into the Phase DBR index or returns -1 if this Phase is not used in the data bridge.
long int IC_xDB_to_xCH (const long int xBR) const
 Converts the IC DBR index into the IC DCH index.
long int DC_xDB_to_xCH (const long int xBR) const
 Converts the DC DBR index into the DC DCH index.
long int Ph_xDB_to_xCH (const long int xBR) const
 Converts the Phase DBR index into the Phase DCH index.
long int Phx_to_DCx (const long int Phx)
 Returns the DCH index of the first DC belonging to the phase with DCH index Phx.
long int PhtoDC_DCH (const long int Phx, long int &nDCinPh)
 Returns the DCH index of the first DC belonging to the phase with DCH index Phx, plus returns through the nDCinPh (reference) parameter the number of DCs included into this phase.
long int DCtoPh_DCH (const long int xCH)
 Returns the DCH index of the Phase to which the Dependent Component with index xCH belongs.
long int PhtoDC_DBR (const long int Phx, long int &nDCinPh)
 Returns the DBR index of the first DC belonging to the phase with DBR index Phx, plus returns through the nDCinPh (reference) parameter the number of DCs included into DBR for this phase.
long int DCtoPh_DBR (const long int xBR)
 Returns the DBR index of the Phase to which the Dependent Component with index xBR belongs.
void packDataBr ()
 Packs GEMIPM calculation results into work node structure.
void unpackDataBr (bool uPrimalSol)
 unpacks work DATABR content into GEMIPM data structure
bool check_TP (double T, double P)
 Checks if given temperature T (K) and pressure P (Pa) fit within the interpolation intervals of the DATACH lookup arrays (returns true) or not (returns false)
long int check_grid_T (double TK)
 Tests TK as a grid point for the interpolation of thermodynamic data.
long int check_grid_P (double P)
 Tests P as a grid point for the interpolation of thermodynamic data.
long int check_grid_TP (double T, double P)
 Tests T (K) and P (Pa) as a grid point for the interpolation of thermodynamic data using DATACH lookup arrays.
long int gridTP () const
 Returns number of temperature and pressure grid points for one dependent component.
double get_Ppa_sat (double Tk)
 Returns 1 if a Psat value corresponding to the temperature of interest was found in the GEMS3K input file.
long int get_grid_index_Ppa_sat (double Tk)
 Returns index of Tk point - Psat point pair.
double DC_G0 (const long int xCH, const double P, const double TK, bool norm=true)
 Retrieves (interpolated) molar Gibbs energy G0(P,TK) value for Dependent Component from the DATACH structure.
double DC_V0 (const long int xCH, const double P, const double TK)
 Retrieves (interpolated, if necessary) molar volume V0(P,TK) value for Dependent Component (in J/Pa) from the DATACH structure.
double DC_H0 (const long int xCH, const double P, const double TK)
 Retrieves (interpolated) molar enthalpy H0(P,TK) value for Dependent Component (in J/mol) from the DATACH structure.
double DC_S0 (const long int xCH, const double P, const double TK)
 Retrieves (interpolated) absolute molar enropy S0(P,TK) value for Dependent Component (in J/K/mol) from the DATACH structure.
double DC_Cp0 (const long int xCH, const double P, const double TK)
 Retrieves (interpolated) constant-pressure heat capacity Cp0(P,TK) value for Dependent Component (in J/K/mol) from the DATACH structure.
double DC_A0 (const long int xCH, const double P, const double TK)
 Retrieves (interpolated) Helmholtz energy of Dependent Component (in J/mol) from the DATACH structure.
double DC_U0 (const long int xCH, const double P, const double TK)
 Retrieves (interpolated) Internal energy of Dependent Component (in J/mol) from the DATACH structure.
void DensArrayH2Ow (const double P, const double TK, vector< double > &DensAW)
 Retrieves (interpolated) density and its derivatives of liquid water at (P,TK) from the DATACH structure or 0.0, if TK (temperature, Kelvin) or P (pressure, Pa) parameters go beyond the valid lookup array intervals or tolerances.
void EpsArrayH2Ow (const double P, const double TK, vector< double > &EpsAW)
 Retrieves (interpolated) dielectric constant and its derivatives of liquid water at (P,TK) from the DATACH structure or 0.0, if TK (temperature, Kelvin) or P (pressure, Pa) parameters go beyond the valid lookup array intervals or tolerances.
double EpsH2Ow (const double P, const double TK)
 Retrieves (interpolated) dielectric constant of liquid water at (P,TK) from the DATACH structure.
double DenH2Ow (const double P, const double TK)
 Retrieves (interpolated) density of liquid water (in kg/m3) at (P,TK) from the DATACH structure.
double VisH2Ow (const double P, const double TK)
 Retrieves (interpolated) viscosity of liquid water (in kg/m3) at (P,TK) from the DATACH structure.
double EpsH2Og (const double P, const double TK)
 Retrieves (interpolated) dielectric constant of H2O vapor at (P,TK) from the DATACH structure.
double DenH2Og (const double P, const double TK)
 Retrieves (interpolated) density of H2O vapor (in kg/m3) at (P,TK) from the DATACH structure.
double Ph_Volume (const long int xBR)
 Retrieves the current phase volume in m3 in the reactive sub-system.
double Ph_Mass (const long int xBR)
 Retrieves the phase mass in kg.
double Ph_SatInd (const long int xph)
 Retrieves the phase saturation index.
double * Ph_BC (const long int xph, double *ARout=0)
 Retrieval of the phase bulk composition into memory indicated by ARout.
double Get_mIC (const long xic)
 Retrieves total dissolved aqueous molality of Independent Component with DBR index xic.
void Set_bIC (const long int xic, const double bIC)
 Sets the amount of IC in the bIC input vector of the work DATABR structure.
double Get_bIC (const long int xic) const
 Retrieves the current amount of Independent Component.
void Set_dll (const long int xdc, const double dll)
 Sets the metastability constraint from below to the amount of DC in the dll vector of the work DATABR structure.
void Set_dul (const long int xdc, const double dul)
 Sets the metastability constraint from above to the amount of DC in the dul vector of the work DATABR structure.
void Set_nDC (const long int xdc, const double nDC)
 Sets the amount of DC in the xDC vector of the work DATABR structure.
double Get_nDC (const long int xdc) const
 Retrieves the current mole amount of Dependent Component.
double Get_muDC (const long int xDC, bool norm=true)
 Retrieval of (dual-thermodynamic) chemical potential of the DC.
double Get_aDC (const long int xdc, bool scale=true)
 Retrieval of (dual-thermodynamic) activity of the DC.
double Get_cDC (const long int xdc)
 Retrieves concentration of Dependent Component in its phase in the respective concentration scale.
double Get_gDC (const long int xdc) const
 Retrieves the activity coefficient of Dependent Component in its phase in the respective scale.
double DCmm (const long int xdc) const
 Retrieves the molar mass of Dependent Component in kg/mol.
double ICmm (const long int xic) const
 Retrieves the molar mass of Independent Component in kg/mol.
double DCaJI (const long int xdc, const long int xic) const
 Retrieves the stoichiometry coefficient a[xdc][xic] of IC in the formula of DC.
void Set_IC_b (const double b_val, const long int xCH)
 Sets the total amount of Independent Component.
double IC_b (const long int xCH) const
 Retrieves the current total amount of Independent Component.
double DC_n (const long int xCH) const
 Retrieves the current mole amount of DC directly from GEM IPM work structure.
double DC_a (const long int xCH)
 Retrieves the current (dual-thermodynamic) activity of DC directly from GEM IPM work structure.
void Get_IPc_IPx_DCc_indices (long &index_phase_aIPx, long &index_phase_aIPc, long &index_phase_aDCc, const long &index_phase)
 Functions for retrieveing and setting values needed for activity coefficient calculation by TSolMod class Sets values of LsMod and LsMdc array.
void Get_NPar_NPcoef_MaxOrd_NComp_NP_DC (long &NPar, long &NPcoef, long &MaxOrd, long &NComp, long &NP_DC, const long &index_phase)
void Get_aIPc (vector< double > &aIPc, const long &index_phase_aIPc, const long &index_phase)
void Get_aIPx (vector< long > &aIPx, const long &index_phase_aIPx, const long &index_phase)
void Get_aDCc (vector< double > &aDCc, const long &index_phase_aDCc, const long &index_phase)
void Set_aIPc (const vector< double > aIPc, const long &index_phase_aIPc, const long &index_phase)
void Set_aDCc (const vector< double > aDCc, const long &index_phase_aDCc, const long &index_phase)
void Set_Tk (double &T_k)
 These methods set contents of fields in the work node structure.
void Set_Pb (double &P_b)
double DC_c (const long int xCH)
 Retrieves the current concentration of Dependent Component in its phase directly from GEM IPM work structure.
double DC_g (const long int xCH) const
 Retrieves the current activity coefficient of DC in its plase directly from GEM IPM work structure.
double DC_lng (const long int xCH) const
 Retrieves the natural logarithm of the internal activity coefficient of species at DCH index xCH.
double DC_mu (const long int xCH, bool norm=true)
 Retrieves the current (dual-thermodynamic) chemical potential of DC directly from GEM IPM work structure.
double DC_mu0 (const long int xCH, bool norm=true)
 Retrieves the standard chemical potential of DC directly from GEM IPM work structure at current pressure and temperature.
virtual void databr_to_vtk (fstream &ff, const char *name, double time, long int cycle, long int nFilds, long int(*Flds)[2])
 Writes work node (DATABR structure) to a text VTK file.

Protected Member Functions

void allocMemory ()
void freeMemory ()
void datach_realloc ()
void datach_free ()
void datach_reset ()
void databr_realloc ()
void databr_reset (DATABR *CNode, long int level=0)
DATABRdatabr_free (DATABR *data_BR_=0)
 Deletes fields of DATABR structure indicated by data_BR_ and sets the pointer data_BR_ to NULL.
void datach_to_file (GemDataStream &ff)
 Writes CSD (DATACH structure) to a binary DCH file.
void datach_from_file (GemDataStream &ff)
 Reads CSD (DATACH structure) from a binary DCH file.
void databr_to_file (GemDataStream &ff)
 Writes node (work DATABR structure) to a binary DBR file.
void databr_from_file (GemDataStream &ff)
 Reads node (work DATABR structure) from a binary DBR file.
void datach_to_text_file (fstream &ff, bool with_comments=true, bool brief_mode=false, const char *path=" ")
 Writes CSD (DATACH structure) to a text DCH file.
void datach_from_text_file (fstream &ff)
 Reads CSD (DATACH structure) from a text DCH file.
void databr_to_text_file (fstream &ff, bool with_comments=true, bool brief_mode=false, const char *path=" ")
 Writes work node (DATABR structure) to a text DBR file.
void databr_from_text_file (fstream &ff)
 Reads work node (DATABR structure) from a text DBR file.
void databr_element_to_vtk (fstream &ff, DATABR *CNode_, long int nfild, long int ndx)
void databr_name_to_vtk (fstream &ff, long int nfild, long int ndx, long int ndx2=0)
void databr_size_to_vtk (long int nfild, long int &nel, long int &nel2)
void databr_head_to_vtk (fstream &ff, const char *name, double time, long cycle, long nx=1, long ny=1, long nz=1)
virtual void InitNodeArray (const char *, long int *, bool, bool)
virtual void setNodeArray (long int, long int *)
virtual long int nNodes () const

Protected Attributes

MULTIpmm
 Pointer to GEM IPM work data structure (ms_multi.h)
TMultimulti
TProfilprofil
DATACHCSD
 Pointer to chemical system data structure CSD (DATACH)
DATABRCNode
 Pointer to a work node data bridge structure (node)
double CalcTime
 GEMIPM2 calculation time, s.
long int PrecLoops
 Number of performed IPM-2 precision refinement loops.
long int NumIterFIA
 Total Number of performed FIA entry iterations.
long int NumIterIPM
 Total Number of performed IPM main iterations.

Detailed Description

Implements a simple C/C++ interface between GEM IPM and FMT codes.

) Implements a simple C/C++ interface between GEM IPM and FMT codes.

Works with DATACH and work DATABR structures without using the TNodearray class.


Member Function Documentation

long int TNode::check_grid_P ( double  P)

Tests P as a grid point for the interpolation of thermodynamic data.

Returns:
index in the lookup grid array or -1 if it is not a grid point
long int TNode::check_grid_T ( double  TK)

Tests TK as a grid point for the interpolation of thermodynamic data.

Returns:
index in the lookup grid array or -1 if it is not a grid point
long int TNode::check_grid_TP ( double  T,
double  P 
)

Tests T (K) and P (Pa) as a grid point for the interpolation of thermodynamic data using DATACH lookup arrays.

Returns:
-1L if interpolation is needed, or 1D index of the lookup array element if TK and P fit within the respective tolerances.
void TNode::databr_to_text_file ( fstream &  ff,
bool  with_comments = true,
bool  brief_mode = false,
const char *  path = " " 
) [protected]

Writes work node (DATABR structure) to a text DBR file.

Parameters:
brief_mode- Do not write data items that contain only default values
with_comments- Write files with comments for all data entries
void TNode::datach_to_text_file ( fstream &  ff,
bool  with_comments = true,
bool  brief_mode = false,
const char *  path = " " 
) [protected]

Writes CSD (DATACH structure) to a text DCH file.

Parameters:
brief_mode- Do not write data items that contain only default values
with_comments- Write files with comments for all data entries
double TNode::DC_a ( const long int  xCH)

Retrieves the current (dual-thermodynamic) activity of DC directly from GEM IPM work structure.

Also activity of a DC not included into DATABR list can be retrieved. If DC has zero amount, its dual-thermodynamic activity is returned anyway. For single condensed phase component, this value has a meaning of the saturation index, also in the presence of metastability constraint(s). These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis DC DCH index
double TNode::DC_A0 ( const long int  xCH,
const double  P,
const double  TK 
)

Retrieves (interpolated) Helmholtz energy of Dependent Component (in J/mol) from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
Returns:
Helmholtz energy (in J/mol) or 7777777., if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DC_c ( const long int  xCH)

Retrieves the current concentration of Dependent Component in its phase directly from GEM IPM work structure.

Also activity of a DC not included into DATABR list can be retrieved. For aqueous species, molality is returned; for gas species, partial pressure; for surface complexes - density in mol/m2; for species in other phases - mole fraction. If DC has zero amount, the function returns 0.0. These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis DC DCH index
double TNode::DC_Cp0 ( const long int  xCH,
const double  P,
const double  TK 
)

Retrieves (interpolated) constant-pressure heat capacity Cp0(P,TK) value for Dependent Component (in J/K/mol) from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
Returns:
Cp0(P,TK) (in J/K/mol) or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DC_g ( const long int  xCH) const [inline]

Retrieves the current activity coefficient of DC in its plase directly from GEM IPM work structure.

Also activity coefficient of a DC not included into DATABR list can be retrieved. If DC has zero amount, this function returns 1.0. These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis DC DCH index
double TNode::DC_G0 ( const long int  xCH,
const double  P,
const double  TK,
bool  norm = true 
)

Retrieves (interpolated) molar Gibbs energy G0(P,TK) value for Dependent Component from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
normdefines in wnich units the value is returned: false - in J/mol; true (default) - in mol/mol
Returns:
G0(P,TK) or 7777777., if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DC_H0 ( const long int  xCH,
const double  P,
const double  TK 
)

Retrieves (interpolated) molar enthalpy H0(P,TK) value for Dependent Component (in J/mol) from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
Returns:
H0(P,TK) (in J/mol) or 7777777., if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DC_lng ( const long int  xCH) const [inline]

Retrieves the natural logarithm of the internal activity coefficient of species at DCH index xCH.

Parameters:
xCHindex of species DCH
double TNode::DC_mu ( const long int  xCH,
bool  norm = true 
)

Retrieves the current (dual-thermodynamic) chemical potential of DC directly from GEM IPM work structure.

Also for any DC not included into DATABR or having zero amount. These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis DC DCH index
normdefines in wnich units the chemical potential value is returned: false - in J/mol; true (default) - in mol/mol
double TNode::DC_mu0 ( const long int  xCH,
bool  norm = true 
)

Retrieves the standard chemical potential of DC directly from GEM IPM work structure at current pressure and temperature.

Also for any DC not included into DATABR or having zero amount. These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis DC DCH index
normdefines in which units the chemical potential value is returned: false - in J/mol; true (default) - in mol/mol
double TNode::DC_n ( const long int  xCH) const [inline]

Retrieves the current mole amount of DC directly from GEM IPM work structure.

Also amount of DCs not included into DATABR list can be retrieved. Internal re-scaling to mass of the system is applied. These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis DC DCH index
double TNode::DC_S0 ( const long int  xCH,
const double  P,
const double  TK 
)

Retrieves (interpolated) absolute molar enropy S0(P,TK) value for Dependent Component (in J/K/mol) from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
Returns:
S0(P,TK) (in J/K/mol) or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DC_U0 ( const long int  xCH,
const double  P,
const double  TK 
)

Retrieves (interpolated) Internal energy of Dependent Component (in J/mol) from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
Returns:
Internal energy (in J/mol) or 7777777., if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DC_V0 ( const long int  xCH,
const double  P,
const double  TK 
)

Retrieves (interpolated, if necessary) molar volume V0(P,TK) value for Dependent Component (in J/Pa) from the DATACH structure.

Parameters:
xCHis the DC DCH index
Ppressure, Pa
TKtemperature, Kelvin
Returns:
V0(P,TK) (in J/Pa) or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DCaJI ( const long int  xdc,
const long int  xic 
) const [inline]

Retrieves the stoichiometry coefficient a[xdc][xic] of IC in the formula of DC.

Parameters:
xdcis DC DBR index
xicis IC DBR index
double TNode::DCmm ( const long int  xdc) const [inline]

Retrieves the molar mass of Dependent Component in kg/mol.

Parameters:
xdcis DC DBR index
double TNode::DenH2Og ( const double  P,
const double  TK 
)

Retrieves (interpolated) density of H2O vapor (in kg/m3) at (P,TK) from the DATACH structure.

Parameters:
Ppressure, Pa
TKtemperature, Kelvin
Returns:
density of H2O vapor or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::DenH2Ow ( const double  P,
const double  TK 
)

Retrieves (interpolated) density of liquid water (in kg/m3) at (P,TK) from the DATACH structure.

Parameters:
Ppressure, Pa
TKtemperature, Kelvin
Returns:
density or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
void TNode::DensArrayH2Ow ( const double  P,
const double  TK,
vector< double > &  DensAW 
)

Retrieves (interpolated) density and its derivatives of liquid water at (P,TK) from the DATACH structure or 0.0, if TK (temperature, Kelvin) or P (pressure, Pa) parameters go beyond the valid lookup array intervals or tolerances.

Parameters:
Prefers to the pressure in Pascal
TKrefers to the temperature in Kelvin
DensAWcontains the density of water (at P and Tk) and its temperature and pressure derivatives
void TNode::EpsArrayH2Ow ( const double  P,
const double  TK,
vector< double > &  EpsAW 
)

Retrieves (interpolated) dielectric constant and its derivatives of liquid water at (P,TK) from the DATACH structure or 0.0, if TK (temperature, Kelvin) or P (pressure, Pa) parameters go beyond the valid lookup array intervals or tolerances.

Parameters:
Prefers to the pressure in Pascal
TKrefers to the temperature in Kelvin
DensAWcontains the permittivity of water (at P and Tk) and its temperature and pressure derivatives
double TNode::EpsH2Og ( const double  P,
const double  TK 
)

Retrieves (interpolated) dielectric constant of H2O vapor at (P,TK) from the DATACH structure.

Parameters:
Ppressure, Pa
TKtemperature, Kelvin
Returns:
dielectric constant of H2O vapor or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::EpsH2Ow ( const double  P,
const double  TK 
)

Retrieves (interpolated) dielectric constant of liquid water at (P,TK) from the DATACH structure.

Parameters:
Ppressure, Pa
TKtemperature, Kelvin
Returns:
dielectric constant or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.
double TNode::GEM_CalcTime ( )

Returns GEMIPM2 calculation time in seconds elapsed during the last call of GEM_run() - can be used for monitoring the performance of calculations.

Returns:
double number, may contain 0.0 if the calculation time is less than the internal time resolution of C/C++ function
void TNode::GEM_from_MT ( long int  p_NodeHandle,
long int  p_NodeStatusCH,
double  p_TK,
double  p_P,
double  p_Vs,
double  p_Ms,
double *  p_bIC,
double *  p_dul,
double *  p_dll,
double *  p_aPH 
)

(8) Loads the GEMS3K input data for a given mass-transport node into the work instance of DATABR structure.

This call is usually preceeding the GEM_run() call

Parameters:
p_NodeHandleNode identification handle
p_NodeStatusCHNode status code (NEED_GEM_SIA or NEED_GEM_AIA)
p_TKTemperature T, Kelvin + - -
p_PPressure P, Pa + - -
p_VsVolume V of reactive subsystem, m3 - - +
p_MsMass of reactive subsystem, kg - - +
p_bICBulk mole amounts of IC [nICb] + - -
p_dulUpper restrictions to amounts of DC [nDCb] + - -
p_dllLower restrictions to amounts of DC [nDCb] + - -
p_aPHSpecific surface areas of phases, m2/kg [nPHb] + - -
void TNode::GEM_from_MT ( long int  p_NodeHandle,
long int  p_NodeStatusCH,
double  p_TK,
double  p_P,
double  p_Vs,
double  p_Ms,
double *  p_bIC,
double *  p_dul,
double *  p_dll,
double *  p_aPH,
double *  p_xDC 
)

(8a) Loads the GEMS3K input data for a given mass-transport node into the work instance of DATABR structure.

This overloaded variant uses the xDC speciation vector for setting the new bulk chemical composition to be used in the next GEM_run() calculation.

Parameters:
p_NodeHandleNode identification handle
p_NodeStatusCHNode status code (NEED_GEM_SIA or NEED_GEM_AIA)
p_TKTemperature T, Kelvin + - -
p_PPressure P, Pa + - -
p_VsVolume V of reactive subsystem, m3 - - +
p_MsMass of reactive subsystem, kg - - +
p_bICBulk mole amounts of IC [nICb] + - -
p_dulUpper restrictions to amounts of DC [nDCb] + - -
p_dllLower restrictions to amounts of DC [nDCb] + - -
p_aPHSpecific surface areas of phases, m2/kg [nPHb] + - -
p_xDCMole amounts of DCs [nDCb] - will be convoluted and added to the bIC GEM input vector (if full speciation and not just increments then p_bIC vector must be zeroed off - it will be calculated from p_xDC and stoichiometry matrix A
void TNode::GEM_from_MT ( long int  p_NodeHandle,
long int  p_NodeStatusCH,
double  p_TK,
double  p_P,
double  p_Vs,
double  p_Ms,
double *  p_bIC,
double *  p_dul,
double *  p_dll,
double *  p_aPH,
double *  p_xDC,
double *  p_gam 
)

(8b) Loads the GEMS3K input data for a given mass-transport node into the work instance of DATABR structure.

In addition, provides access to speciation vector p_xDC and DC activity coefficients p_gam that will be used in GEM "smart initial approximation" SIA mode if dBR->NodeStatusCH == NEED_GEM_SIA (5) and uPrimalSol = true are set for the GEM_run() call (see Section 2) . This works only when the DATACH structure contains a full list of Dependent Components used in GEM IPM2 calculations.

Parameters:
p_NodeHandleNode identification handle
p_NodeStatusCHNode status code (NEED_GEM_SIA or NEED_GEM_AIA)
p_TKTemperature T, Kelvin + - -
p_PPressure P, Pa + - -
p_VsVolume V of reactive subsystem, m3 - - +
p_MsMass of reactive subsystem, kg - - +
p_bICBulk mole amounts of IC [nICb] + - -
p_dulUpper restrictions to amounts of DC [nDCb] + - -
p_dllLower restrictions to amounts of DC [nDCb] + - -
p_aPHSpecific surface areas of phases, m2/kg [nPHb] + - -
p_xDCMole amounts of DCs [nDCb] - old primal soln. + - -
p_gamDC activity coefficients [nDCb] - old primal s. + - -
long int TNode::GEM_init ( const char *  ipmfiles_lst_name,
const char *  dbrfiles_lst_name = 0,
long int *  nodeTypes = 0,
bool  getNodT1 = false 
)

(1) Initialization of GEM IPM3 data structures in coupled programs that use GEMS3K module.

Also reads in the IPM, DCH and one or many DBR text input files.

Parameters:
ipmfiles_lst_namepointer to a null-terminated C string with a path to a text file containing the list of names of GEMS3K input files. Example: file "test.lst" with a content: -t "dch.dat" "ipm.dat" "dbr-0.dat" (-t tells that input files are in text format)
dbrfiles_lst_nameoptional parameter (used only at the TNodeArray level) - a pointer to a null-terminated C string with a path to a text file containing the list of comma-separated names of DBR input files. Example: file "test-dbr.lst" with a content: "dbr-0.dat" , "dbr-1.dat" , "dbr-2.dat"
nodeTypesoptional parameter used only at the TNodeArray level, the initial node contents from DATABR files will be distributed among nodes in array according to the distribution index list nodeTypes
getNodT1optional parameter used only when reading multiple DBR files after the modeling task interruption in GEM-Selektor
Returns:
0 if successful; 1 if input file(s) were not found or corrupt; -1 if internal memory allocation error occurred.
long int TNode::GEM_Iterations ( long int &  PrecLoops,
long int &  NumIterFIA,
long int &  NumIterIPM 
)

To obtain the number of GEM IPM2 iterations performed during the last call of GEM_run() e.g.

for monitoring the performance of GEMS3K in AIA or SIA modes, or for problem diagnostics. Parameters: long int variables per reference (must be allocated before calling GEM_Iterations(), previous values will be lost. See Return values.

Returns:
Function Total number of EFD + IPM iterations from the last call to GEM_run()
Parameters:
PrecLoopsNumber of performed IPM-2 precision refinement loops
NumIterFIATotal number of performed MBR() iterations to obtain a feasible initial approximation for the IPM algorithm.
NumIterIPMTotal number of performed IPM main descent algorithm iterations.
void TNode::GEM_print_ipm ( const char *  fname)

(4) Produces a formatted text file with detailed contents (scalars and arrays) of the GEM IPM work structure.

This call is useful when GEM_run() returns with a NodeStatusCH value indicating a GEM calculation error (see above). Another use is for a detailed comparison of a test system calculation after the version upgrade of GEMS3K.

Parameters:
fnamenull-terminated (C) string containing a full path to the disk file to be written. NULL - the disk file name path stored in the dbr_file_name field of the TNode class instance will be used, extended with ".dump.out". Usually the dbr_file_name field contains the path to the last input DBR file.
long int TNode::GEM_read_dbr ( const char *  fname,
bool  binary_f = false 
)

(5) Reads another DBR file (with input system composition, T,P etc.) \ .

The DBR file must be compatible with the currently loaded IPM and DCH files (see description of GEM_init() function call).

Parameters:
fnameNull-terminated (C) string containing a full path to the input DBR disk file.
binary_fFlag defining whether the file specified in fname is in text fromat (false or 0, default) or in binary format (true or 1)
Returns:
0 if successful; 1 if input file(s) has not found been or is corrupt; -1 if internal memory allocation error occurred.
void TNode::GEM_restore_MT ( long int &  p_NodeHandle,
long int &  p_NodeStatusCH,
double &  p_TK,
double &  p_P,
double &  p_Vs,
double &  p_Ms,
double *  p_bIC,
double *  p_dul,
double *  p_dll,
double *  p_aPH 
)

(6) Passes (copies) the GEMS3K input data from the work instance of DATABR structure.

This call is useful after the GEM_init() (1) and GEM_run() (2) calls to initialize the arrays which keep the chemical data for all nodes used in the mass-transport model.

Parameters:
p_NodeHandleNode identification handle
p_NodeStatusCHNode status code; see typedef NODECODECH
p_TKTemperature T, Kelvin + - -
p_PPressure P, Pa + - -
p_VsVolume V of reactive subsystem, m3 (+) - +
p_MsMass of reactive subsystem, kg - - +
p_bICBulk mole amounts of IC [nICb] + - -
p_dulUpper restrictions to amounts of DC [nDCb] + - -
p_dllLower restrictions to amounts of DC [nDCb] + - -
p_aPHSpecific surface areas of phases,m2/kg[nPHb] + - -
long int TNode::GEM_run ( bool  uPrimalSol)

(2) Main call for GEM IPM calculations using the input bulk composition, temperature, pressure and metastability constraints provided in the work instance of DATABR structure.

Actual calculation will be performed only when dBR->NodeStatusCH == NEED_GEM_SIA (5) or dBR->NodeStatusCH = NEED_GEM_AIA (1). By other values of NodeStatusCH, no calculation will be performed and the status will remain unchanged. In "smart initial approximation" (SIA) mode, the program can automatically switch into the "automatic initial approximation" (AIA) mode and return OK_GEM_AIA instead of OK_GEM_SIA.

Parameters:
uPrimalSolflag to define the mode of GEM smart initial approximation (only if dBR->NodeStatusCH = NEED_GEM_SIA has been set before GEM_run() call). false (0) - use speciation and activity coefficients from previous GEM_run() calculation true (1) - use speciation provided in the DATABR memory structure (e.g. after reading the DBR file)
Returns:
NodeStatusCH (the same as set in dBR->NodeStatusCH). Possible values (see "databr.h" file for the full list)
void TNode::GEM_set_MT ( double  p_Tm,
double  p_dt 
)

(9) Optional, for passing the current mass transport iteration information into the work DATABR structure (e.g.

\ for using it in tracing/debugging or in writing DBR files for nodes)

Parameters:
p_TmActual total simulation time, s + - -
p_dtActual time step, s + - -
void TNode::GEM_to_MT ( long int &  p_NodeHandle,
long int &  p_NodeStatusCH,
long int &  p_IterDone,
double &  p_Vs,
double &  p_Ms,
double &  p_Gs,
double &  p_Hs,
double &  p_IC,
double &  p_pH,
double &  p_pe,
double &  p_Eh,
double *  p_rMB,
double *  p_uIC,
double *  p_xDC,
double *  p_gam,
double *  p_xPH,
double *  p_vPS,
double *  p_mPS,
double *  p_bPS,
double *  p_xPA,
double *  p_aPH,
double *  p_bSP 
)

(7) Retrieves the GEMIPM2 chemical speciation calculation results from the work DATABR structure instance into memory provided by the mass transport part.

Dimensions and order of elements in the arrays must correspond to those in currently existing DATACH memory structure.

Parameters:
p_NodeHandleNode identification handle
p_NodeStatusCHNode status code (changed after GEM calculation); see typedef NODECODECH
p_IterDoneNumber of iterations performed in the last GEM IPM calculation
p_VsTotal volume V of reactive subsystem at given P,T, m3 - - + +
p_MsTotal mass of the reactive subsystem, kg - - + +
p_GsTotal Gibbs energy of the reactive subsystem, J - - + +
p_HsTotal enthalpy of reactive subsystem, J (reserved) - - + +
p_ICEffective aqueous ionic strength, molal - - + +
p_pHpH of aqueous solution - - + +
p_pepe of aqueous solution - - + +
p_EhEh of aqueous solution, V - - + +
p_rMBMole balance residuals for Independent Components [nICb] - - + +
p_uICDual solution: IC chemical potentials, mol/mol [nICb] - - + +
p_xDCPrimal solution: DC mole amounts [nDCb] - - + +
p_gamExternal activity coefficients of DC [nDCb] - - + +
p_xPHTotal mole amounts of all phases [nPHb] - - + +
p_vPSTotal volumes of multicomponent phases, m3 [nPSb] - - + +
p_mPSTotal mass of multicomponent phase (carrier),kg [nPSb] - - + +
p_bPSBulk compositions of multicomponent phases [nPSb][nICb] - - + +
p_xPAAmount of carrier in a multicomponent asymmetric phase[nPSb]- - + +
p_aPHCalculated surface areas of phases (m2) [nPHb] - - + +
p_bSPBulk composition of all solids, moles [nICb] - - + +
void TNode::GEM_write_dbr ( const char *  fname,
bool  binary_f = false,
bool  with_comments = true,
bool  brief_mode = false 
)

(3) Writes the contents of the work instance of the DATABR structure into a disk file with path name fname.

Parameters:
fnamenull-terminated (C) string containing a full path to the DBR disk file to be written. NULL - the disk file name path stored in the dbr_file_name field of the TNode class instance will be used, extended with ".out". Usually the dbr_file_name field contains the path to the last input DBR file.
binary_fdefines if the file is to be written in binary format (true or 1, good for interruption of coupled modeling task if called in the loop for each node), or in text format (false or 0, default).
with_comments(text format only): defines the mode of output of comments written before each data tag and content in the DBR file. If set to true (1), the comments will be written for all data entries (default). If false (0), comments will not be written.
brief_modeif true, tells that do not write data items, that contain only default values in text format
double TNode::Get_aDC ( const long int  xdc,
bool  scale = true 
)

Retrieval of (dual-thermodynamic) activity of the DC.

Parameters:
xdcis DC DBR index
scaleif true then activity is returned, if false then log10(activity)
double TNode::Get_bIC ( const long int  xic) const [inline]

Retrieves the current amount of Independent Component.

Parameters:
xicis IC DBR index
double TNode::Get_cDC ( const long int  xdc)

Retrieves concentration of Dependent Component in its phase in the respective concentration scale.

For aqueous species, molality is returned; for gas species, mole fraction not partial pressure; for surface complexes - molality; for species in other phases - mole fraction.

Parameters:
xdcis DC DBR index
Returns:
0.0, if DC has zero amount.
double TNode::Get_gDC ( const long int  xdc) const [inline]

Retrieves the activity coefficient of Dependent Component in its phase in the respective scale.

Parameters:
xdcis DC DBR index
Returns:
1.0, if DC has zero amount.
double TNode::Get_mIC ( const long  xic)

Retrieves total dissolved aqueous molality of Independent Component with DBR index xic.

Parameters:
xicis IC DBR index
Returns:
total dissolved aqueous molality or 0.0, if there is no water in the node or no aqueous phase in DATACH.
double TNode::Get_muDC ( const long int  xDC,
bool  norm = true 
)

Retrieval of (dual-thermodynamic) chemical potential of the DC.

Parameters:
xdcis DC DBR index
normdefines the scale: if true (1) then in mol/mol, otherwise in J/mol
double TNode::Get_nDC ( const long int  xdc) const [inline]

Retrieves the current mole amount of Dependent Component.

Parameters:
xdcis DC DBR index
double TNode::IC_b ( const long int  xCH) const [inline]

Retrieves the current total amount of Independent Component.

Also amount of ICs not included into DATABR list can be retrieved. Internal re-scaling to mass of the system is applied These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis IC DCH index
double TNode::ICmm ( const long int  xic) const [inline]

Retrieves the molar mass of Independent Component in kg/mol.

Parameters:
xicis IC DBR index
double * TNode::Ph_BC ( const long int  xph,
double *  ARout = 0 
)

Retrieval of the phase bulk composition into memory indicated by ARout.

This function works for multicomponent and for single-component phases

Parameters:
xphis DBR phase index
ARoutis array of at least [dCH->nICb elements] or ARout = NULL
Returns:
pointer to ARout which may also be allocated inside of Ph_BC() in the case if parameter ARout = NULL is specified; to avoid a memory leak, you will have to free this memory wherever appropriate.
double TNode::Ph_Mass ( const long int  xBR)

Retrieves the phase mass in kg.

Works for multicomponent and for single-component phases.

Parameters:
xphis DBR phase index
Returns:
the phase mass in kg or 0.0, if the phase mole amount is zero.
double TNode::Ph_SatInd ( const long int  xph)

Retrieves the phase saturation index.

Works for multicomponent and for single-component phases.

Parameters:
xphis DBR phase index
Returns:
the phase saturation index or 0.0, if the phase mole amount is zero.
double TNode::Ph_Volume ( const long int  xBR)

Retrieves the current phase volume in m3 in the reactive sub-system.

Works both for multicomponent and for single-component phases.

Parameters:
xphis DBR phase index
Returns:
the current phase volume in m3 or 0.0, if the phase mole amount is zero.
void TNode::Set_bIC ( const long int  xic,
const double  bIC 
) [inline]

Sets the amount of IC in the bIC input vector of the work DATABR structure.

Parameters:
xicis IC DBR index
bICis amount of IC
void TNode::Set_dll ( const long int  xdc,
const double  dll 
) [inline]

Sets the metastability constraint from below to the amount of DC in the dll vector of the work DATABR structure.

Parameters:
xdcis DC DBR index
void TNode::Set_dul ( const long int  xdc,
const double  dul 
) [inline]

Sets the metastability constraint from above to the amount of DC in the dul vector of the work DATABR structure.

Parameters:
xdcis DC DBR index
void TNode::Set_IC_b ( const double  b_val,
const long int  xCH 
) [inline]

Sets the total amount of Independent Component.

Also amount of ICs not included into DATABR list can be retrieved. Internal re-scaling to mass of the system is applied. These methods can only be used for the current work node (direct access to GEM IPM data)

Parameters:
xCHis IC DCH index
void TNode::Set_nDC ( const long int  xdc,
const double  nDC 
) [inline]

Sets the amount of DC in the xDC vector of the work DATABR structure.

Parameters:
xdcis DC DBR index
void TNode::setNodeHandle ( long int  jj) [inline]

Set current node identification handle to value of.

Parameters:
jj
double TNode::VisH2Ow ( const double  P,
const double  TK 
)

Retrieves (interpolated) viscosity of liquid water (in kg/m3) at (P,TK) from the DATACH structure.

Parameters:
Ppressure, Pa
TKtemperature, Kelvin
Returns:
viscosity or 0.0, if TK or P go beyond the valid lookup array intervals or tolerances.

Member Data Documentation

DATABR* TNode::CNode [protected]

Pointer to a work node data bridge structure (node)

used for exchanging input data and results between FMT and GEM IPM


The documentation for this class was generated from the following files: