A fundamental problem in biochemical thermodynamics pertains to the treatment of complex biochemical equilibria in aqueous solution. The problem is interesting because it involves the coupling of several different types of equilibria with ligands such as
and
, as well as the effects of temperature and ionic strength. The effects caused by changes in these variables can dramatically affect the position of equilibrium. The problem is important because the reactions that are dealt with form the basis of metabolism and energy utilization in living systems [1].
Traditional methods [2-4] for performing these calculations involve first the derivation of analytic expressions for the fractions of the various species in solution. Since these fractions must be introduced into the computer codes, some programming effort is required for each new problem. However, a more general approach has been implemented by Alberty [5] and by Alberty and Krambeck [6] in which the chemical equilibrium problem is formulated in terms of the conservation matrix for the system of reactions [7, 8]. A Newton--Raphson minimization of the Gibbs free energy is then used to solve for the molalities of the species.
The package described here builds upon these earlier efforts by including the effects of ionic strength in the computation and by allowing for the input of the data in a convenient chemical format. Thus, the transpose of the stoichiometric number matrix
is constructed directly from the chemical equations. The package allows one to calculate the molalities and mole fractions of the species in the solution, the activity coefficients of the species and the activity of the solvent (
), values of apparent equilibrium constants
, standard transformed Gibbs free energies
and standard transformed enthalpies of reaction
, calorimetric enthalpies of reaction
, and changes in the binding
of ligands such as
and
. The package allows for the calculation of these thermodynamic quantities as functions of pX (e.g., pH and pMg), ionic strength, and temperature. It is also capable of treating systems involving several overall biochemical reactions at a specified temperature, ionic strength, and pX. Thus, the package allows for the convenient calculation of the quantities that are of primary interest to biochemical thermodynamics. Although designed for biochemical applications, this package can also be used for the study of other complex systems at equilibrium.
The reader wishing to gain additional background on the thermodynamics of biochemical reactions is referred to [9] and the references cited therein. The terminology, symbols, and units used for the physical quantities in the present article adhere
generally to IUPAC recommendations [9, 10] (see glossary).
Mathematica [11] has several features that are useful for the implementation of a computer program suitable for the above purposes. Firstly, the capabilities to conveniently perform symbolic manipulations and matrix calculations are particularly important. Also, the Newton-Raphson minimization used by Krambeck and Alberty [6] has been implemented in a very brief but elegant Mathematica module. This module has also proven to be extremely powerful and rugged in its ability to provide reliable values of the molalities of the species in solution at equilibrium. The capability of Mathematica to carry out equilibrium calculations with arbitrary precision will also almost certainly prove valuable in dealing with some classes of problems. Finally, the capability to conveniently integrate graphics with calculations in a notebook format is a significant convenience to the user.
Treatment of Activity Coefficients
In these calculations, the activity coefficients of the species are represented by the extended
Debye-Hückel equation:
The calculated values of the activity coefficients are used to adjust the values of the (standard) equilibrium constants, which are input as data and pertain to
to the actual ionic strength of the solution. Values of the standard molar enthalpies of reaction
at
are also adjusted to the appropriate ionic strength. The activity of water
is calculated by using a
Gibbs-Duhem integration [4]. Following an initial guess for the ionic strength I, this quantity is calculated at each cycle of the equilibrium calculations until it has converged to within 1 part in
. The standard state used in the thermodynamic calculations is that of the hypothetical ideal solution, namely
as
. As a consequence of the
Gibbs-Duhem equation,
also
as
.
The package consists of 35 modules, each of which performs a specific function. However, not all of the modules need to be directly accessed. Documentation for use of BioEqCalc is primarily by way of example, solving generic problems that users can copy and adapt for their own use. In these examples, the names of functions will be denoted in the text in bold with [] following the name of the function; the names of variables or parameters used by the functions will also be in bold. These functions are also defined with usage statements in the package BioEqCalc. The package also uses several global variables which are summarized below.