R. Feistel et al.: Oceanographic application and numerical implementation of TEOS-IO: Part 1
661
www.ocean-sci.net/6/633/2010/
Ocean Sci., 6, 633-677, 2010
be unstable if thermodynamic stability criteria such as posi
tive compressibility are violated. (Independent of any simul
taneous existence of other phases, a negative compressibil
ity would amplify any pressure-density fluctuation, causing
the fluid to collapse.) The boundary between metastable and
unstable existence is regarded as the spinodal line beyond
which the phase can no longer stably and homogeneously
exist.
Since it is practically impossible to measure thermody
namic properties on or beyond the spinodal, its location
in the phase diagram is not exactly known from empirical
equations of state. Even though the IAPWS-95 formulation
extrapolates well into the metastable regions (Wagner and
Pilli), 2002; Feistel et al., 2008), with increasing distance
from the saturation line the values computed from the
Helmholtz function outside its validity range will become
unreliable. This turns out not to be a major issue for de
termination of the seawater Gibbs function since we require
only relatively small excursions into the metastable regions
to deal with the effects of shifted phase transition boundaries
in the presence of sea salt.
Panels a and b of Fig. Al indicate the initializations used in
the library for our iterative solutions for the liquid and vapour
phases, respectively. Figure A2 provides additional infor
mation regarding the industrial formulation IF-97 (IAPWS,
2007; Wagner and Kretzschmar, 2008) referred to in Fig. Al.
Starting in either fluid state at a point near CP located
along the saturation line joining TP and CP, we may circum
scribe the critical point along a closed T—P path. Along
any such curve, the properties change only gradually; noth
ing like a transition point between liquid and vapour is en
countered. Since, for numerical purposes, we distinguish be
tween the Gibbs functions of liquid, Eq. (Al), and vapour,
Eq. (A2), we need to specify such a transition point for tech
nical rather than for physical reasons. Here we define the
Gibbs functions Eqs. (Al) and (A2) to be different at sub-
critical conditions (T<Tq and P<Pq) and to be identical at
supercritical conditions (T>Tc or P>Pc)■ The critical tem
perature of water is 7c=647.096K and the critical density
is pc=322kgm -3 (IAPWS, 2009a); the critical pressure fol
lows from Eq. (4.1) to be Pc=22.064MPa. In order to cover
metastable states of liquid water as required in the regions of
vapour-pressure lowering or freezing-point lowering caused
by the presence of dissolved sea salt, the Gibbs function for
liquid water is also available for T and P in regions extend
ing somewhat beyond the saturation curve and beyond the
melting curve.
According to our numerical definition of liquid, vapour
and fluid states, the initial values required for the iteration of
Eq. (A3) can be chosen identically for the fluid density in the
supercritical region (T>Tc or P>Pc), as shown in Fig. Al,
but must be different for liquid and vapour in the subcriticai
quarter (T<Tq and P<Pq). In the subcriticai region, sepa
rate Gibbs functions are available from the industrial formu
lation IF-97 (IAPWS, 2007; Wagner and Kretzschmar, 2008)
Fig. A2. Thermodynamic relations available from the Industrial
Formulation IF-97 (IAPWS, 2007; Wagner and Kretzschmar, 2008)
defined in different temperature-pressure regions, derived with re
duced accuracy from the IAPWS-95 Helmholtz function for differ
ent independent variables. Here we make use of the Gibbs func
tions “g(p, T)" available in region 1 (liquid/fluid) and in region 2
(vapour/fluid) as shown in Fig. Al. In region 3, separate equations
for the specific volume are available for various subregions which
are not used here. Region 4 is the saturation curve. Graphics repro
duced from IAPWS (2007), with permission of IAPWS.
in regions 1 and 2 (Fig. A2) as defined therein which pro
vide excellent starting values for the liquid and the vapour
state. These Gibbs functions can also be used for the fluid,
region 1 below 623.15 K and 100 MPa, and region 2 between
273.15 K and 1073.15 K, and below 16.529 MPa. In the sub
limation region and in the supercritical region, the ideal-gas
density, p=P/(RifjT), provides a sufficient starting estimate
below 273.15 K and above 650 K, and a constant value of
p=1000kgm -3 can be used below 650K, Fig. Al. These
latter choices are sufficient to ensure numerical convergence
but do not necessarily optimize the speed of the code. Addi
tional considerations apply to the immediate neighbourhood
of the critical point as discussed below. Note that all of these
rules are built into the library routines discussed in Part 2
(Wright et al., 2010a) so that the user can make use of the
routines without dealing with (or even being fully aware of)
the details.
The critical region is defined here as the T — P rectangle
623.15-650K and 16.529-35 MPa, Fig. Al. The coefficients
of an auxiliary cubic polynomial equation of state
have been determined by regression to IAPWS-95 data
points in the stable liquid, vapour and fluid region with an
r.m.s. deviation of 1% or less for each phase; the resulting
coefficients are given in Table Al. The cubic polynomial
used for Eq. (A4) permits analytical inversion to determine
p (T. p) for both the liquid and the vapour branches. The
critical point of Eq. (A4) was chosen to be identical with the