eos80

seawater.eos80.adtg(s, t, p)

Calculates adiabatic temperature gradient as per UNESCO 1983 routines.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature [℃ (ITS-90)]

parray_like

pressure [db]

Returns:
adtgarray_like

adiabatic temperature gradient [℃ db -1]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

[2]

Bryden, H. 1973. New Polynomials for thermal expansion, adiabatic temperature gradient and potential temperature of sea water. Deep-Sea Res. Vol20,401-408. doi:10.1016/0011-7471(73)90063-6

Examples

>>> # Data from UNESCO 1983 p45.
>>> import seawater as sw
>>> from seawater.library import T90conv
>>> t = T90conv(
...     [
...         [0, 0, 0, 0, 0, 0],
...         [10, 10, 10, 10, 10, 10],
...         [20, 20, 20, 20, 20, 20],
...         [30, 30, 30, 30, 30, 30],
...         [40, 40, 40, 40, 40, 40],
...     ]
... )
>>> s = [
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
... ]
>>> p = [0, 5000, 10000, 0, 5000, 10000]
>>> sw.adtg(s, t, p)
array([[1.6871000e-05, 1.0470000e-04, 1.6942600e-04, 3.5803000e-05,
        1.1795650e-04, 1.7700700e-04],
       [1.0019458e-04, 1.6095905e-04, 2.0687417e-04, 1.1488728e-04,
        1.7136420e-04, 2.1299177e-04],
       [1.7381984e-04, 2.1353400e-04, 2.4448376e-04, 1.8427324e-04,
        2.2108780e-04, 2.4913796e-04],
       [2.4172046e-04, 2.6476410e-04, 2.8295959e-04, 2.4793456e-04,
        2.6946655e-04, 2.8615039e-04],
       [3.0787012e-04, 3.1698860e-04, 3.2300648e-04, 3.0984492e-04,
        3.1883970e-04, 3.2473388e-04]])
seawater.eos80.alpha(s, t, p, *, pt=False)

Calculate the thermal expansion coefficient.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature or potential temperature [℃ (ITS-90)]

parray_like

pressure [db].

ptbool

True if temperature is potential, default is False

Returns:
alphaarray_like

thermal expansion coeff \(\\alpha\) [℃ -1]

References

[1]

McDougall, Trevor J., 1987: Neutral Surfaces. J. Phys. Oceanogr., 17, 1950-1964. doi: 10.1175/1520-0485(1987)017<1950:NS>2.0.CO;2

Examples

>>> # Data from McDougall 1987
>>> import seawater as sw
>>> s, t, p = 40, 10, 4000
>>> sw.alpha(s, t, p, pt=True)
0.00025061316481624323
seawater.eos80.aonb(s, t, p, *, pt=False)

Calculate \(\\alpha/\\beta\).

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature or potential temperature [℃ (ITS-90)]

parray_like

pressure [db].

ptbool

True if temperature is potential, default is False

Returns:
aonbarray_like

\(\\alpha/\\beta\) [psu ℃ -1]

References

[1]

McDougall, Trevor J., 1987: Neutral Surfaces. J. Phys. Oceanogr., 17, 1950-1964. doi: 10.1175/1520-0485(1987)017<1950:NS>2.0.CO;2

Examples

>>> # Data from McDougall 1987.
>>> import seawater as sw
>>> s, t, p = 40, 10, 4000
>>> sw.aonb(s, t, p, pt=True)
0.347650567047807
seawater.eos80.beta(s, t, p, *, pt=False)

Calculate the saline contraction coefficient \(\\beta\) as defined by T.J. McDougall.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature or potential temperature [℃ (ITS-90)]

parray_like

pressure [db].

ptbool

True if temperature is potential, default is False

Returns:
betaarray_like

saline Contraction Coefficient [psu -1]

References

[1]

McDougall, Trevor J., 1987: Neutral Surfaces. J. Phys. Oceanogr., 17, 1950-1964. doi: 10.1175/1520-0485(1987)017<1950:NS>2.0.CO;2

Examples

>>> # Data from McDougall 1987
>>> import seawater as sw
>>> s, t, p = 40, 10, 4000
>>> sw.beta(s, t, p, pt=True)
0.0007208766174161893
seawater.eos80.cp(s, t, p)

Heat Capacity of Sea Water using UNESCO 1983 polynomial.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature [℃ (ITS-90)]

parray_like

pressure [db].

Returns:
cparray_like

specific heat capacity [J kg -1 C -1]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

Examples

>>> # Data from Pond and Pickard Intro. Dyn. Oceanography 2nd ed. 1986.
>>> import seawater as sw
>>> from seawater.library import T90conv
>>> t = T90conv(
...     [
...         [0, 0, 0, 0, 0, 0],
...         [10, 10, 10, 10, 10, 10],
...         [20, 20, 20, 20, 20, 20],
...         [30, 30, 30, 30, 30, 30],
...         [40, 40, 40, 40, 40, 40],
...     ]
... )
>>> s = [
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
... ]
>>> p = [0, 5000, 10000, 0, 5000, 10000]
>>> sw.cp(s, t, p)
array([[4048.4405375 , 3896.25585   , 3807.7330375 , 3986.53309476,
        3849.26094605, 3769.11791286],
       [4041.8276691 , 3919.5550066 , 3842.3111366 , 3986.34061786,
        3874.72665865, 3804.415624  ],
       [4044.8438591 , 3938.5978466 , 3866.7400391 , 3993.85441786,
        3894.99294519, 3828.29059113],
       [4049.0984351 , 3952.0375476 , 3882.9855526 , 4000.68382238,
        3909.24271128, 3844.32151784],
       [4051.2244911 , 3966.1132036 , 3905.9162711 , 4003.46192541,
        3923.89463092, 3868.28959814]])
seawater.eos80.dens(s, t, p)

Density of Sea Water using UNESCO 1983 (EOS 80) polynomial.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature [℃ (ITS-90)]

parray_like

pressure [db].

Returns:
densarray_like

density [kg m 3]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

[2]

Millero, F.J., Chen, C.T., Bradshaw, A., and Schleicher, K. A new high pressure equation of state for seawater. Deap-Sea Research., 1980, Vol27A, pp255-264. doi:10.1016/0198-0149(80)90016-3

Examples

>>> # Data from Unesco Tech. Paper in Marine Sci. No. 44, p22.
>>> import seawater as sw
>>> from seawater.library import T90conv
>>> s = [0, 0, 0, 0, 35, 35, 35, 35]
>>> t = T90conv([0, 0, 30, 30, 0, 0, 30, 30])
>>> p = [0, 10000, 0, 10000, 0, 10000, 0, 10000]
>>> sw.dens(s, t, p)
array([ 999.842594  , 1045.33710972,  995.65113374, 1036.03148891,
       1028.10633141, 1070.95838408, 1021.72863949, 1060.55058771])
seawater.eos80.dens0(s, t)

Density of Sea Water at atmospheric pressure.

Parameters:
s(p=0)array_like

salinity [psu (PSS-78)]

t(p=0)array_like

temperature [℃ (ITS-90)]

Returns:
dens0(s, t)array_like

density [kg m 3] of salt water with properties (s, t, p=0) 0 db gauge pressure

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

[2]

Millero, F.J. and Poisson, A. International one-atmosphere equation of state of seawater. Deep-Sea Res. 1981. Vol28A(6) pp625-629. doi:10.1016/0198-0149(81)90122-9

Examples

>>> # Data from UNESCO Tech. Paper in Marine Sci. No. 44, p22
>>> import seawater as sw
>>> from seawater.library import T90conv
>>> s = [0, 0, 0, 0, 35, 35, 35, 35]
>>> t = T90conv([0, 0, 30, 30, 0, 0, 30, 30])
>>> sw.dens0(s, t)
array([ 999.842594  ,  999.842594  ,  995.65113374,  995.65113374,
       1028.10633141, 1028.10633141, 1021.72863949, 1021.72863949])
seawater.eos80.dpth(p, lat)

Calculates depth in meters from pressure in dbars.

Parameters:
parray_like

pressure [db].

latnumber or array_like

latitude in decimal degrees north [-90..+90].

Returns:
zarray_like

depth [meters]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

Examples

>>> # UNESCO 1983 data p30.
>>> import seawater as sw
>>> lat = [0, 30, 45, 90]
>>> p = [
...     [500, 500, 500, 500],
...     [5000, 5000, 5000, 5000],
...     [10000, 10000, 10000, 10000],
... ]
>>> sw.dpth(p, lat)
array([[ 496.65299239,  495.99772917,  495.3427354 ,  494.03357499],
       [4915.04099112, 4908.55954332, 4902.08075214, 4889.13132561],
       [9725.47087508, 9712.6530721 , 9699.84050403, 9674.23144056]])
seawater.eos80.fp(s, p)

Freezing point of Sea Water using UNESCO 1983 polynomial.

Parameters:
sarray_like

salinity [psu (PSS-78)]

parray_like

pressure [db]

Returns:
fparray_like

freezing point temperature [℃ (ITS-90)]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

Examples

>>> # UNESCO DATA p.30.
>>> import seawater as sw
>>> s = [[5, 10, 15, 20, 25, 30, 35, 40], [5, 10, 15, 20, 25, 30, 35, 40]]
>>> p = [[0, 0, 0, 0, 0, 0, 0, 0], [500, 500, 500, 500, 500, 500, 500, 500]]
>>> sw.fp(s, p)
array([[-0.27369757, -0.54232831, -0.81142026, -1.0829461 , -1.35804594,
        -1.63748903, -1.9218401 , -2.2115367 ],
       [-0.65010724, -0.91873798, -1.18782992, -1.45935577, -1.73445561,
        -2.01389869, -2.29824976, -2.58794636]])
seawater.eos80.g(lat, z=0)

Calculates acceleration due to gravity as function of latitude.

Parameters:
latarray_like

latitude in decimal degrees north [-90..+90].

znumber or array_like. Default z = 0

height in meters (+ve above sea surface, -ve below).

Returns:
garray_like

gravity [m s -2]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

[2]

A.E. Gill 1982. p.54 Eqn. 3.7.15 “Atmosphere-Ocean Dynamics” Academic Press: New York. ISBN: 0-12-283522-0

Examples

>>> import seawater as sw
>>> sw.g(45, z=0)
9.8061898752054
seawater.eos80.pden(s, t, p, pr=0)

Calculates potential density of water mass relative to the specified reference pressure by pden = dens(S, ptmp, PR).

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature [℃ (ITS-90)]

parray_like

pressure [db].

prnumber

reference pressure [db], default = 0

Returns:
pdenarray_like

potential density relative to the ref. pressure [kg m :sup:3]

References

[1]

A.E. Gill 1982. p.54 Eqn. 3.7.15 “Atmosphere-Ocean Dynamics” Academic Press: New York. ISBN: 0-12-283522-0

Examples

>>> # Data from Unesco Tech. Paper in Marine Sci. No. 44, p22.
>>> import seawater as sw
>>> from seawater.library import T90conv
>>> s = [0, 0, 0, 0, 35, 35, 35, 35]
>>> t = T90conv([0, 0, 30, 30, 0, 0, 30, 30])
>>> p = [0, 10000, 0, 10000, 0, 10000, 0, 10000]
>>> sw.pden(s, t, p)
array([ 999.842594  ,  999.79523994,  995.65113374,  996.36115932,
       1028.10633141, 1028.15738545, 1021.72863949, 1022.59634627])

\(\\sigma_{4}\) (at 4000 db)

>>> sw.pden(s, t, p, 4000) - 1000
array([19.2895493 , 19.33422519, 12.43271053, 13.27563816, 46.30976432,
       46.48818851, 37.76150878, 38.74500757])
seawater.eos80.pres(depth, lat)

Calculates pressure in dbars from depth in meters.

Parameters:
deptharray_like

depth [meters]

latarray_like

latitude in decimal degrees north [-90..+90]

Returns:
parray_like

pressure [db]

References

[1]

Saunders, Peter M., 1981: Practical Conversion of Pressure to Depth. J. Phys. Oceanogr., 11, 573-574. doi: 10.1175/1520-0485(1981)011<0573:PCOPTD>2.0.CO;2

Examples

>>> import seawater as sw
>>> depth, lat = 7321.45, 30
>>> sw.pres(depth, lat)
7500.006513011802
seawater.eos80.ptmp(s, t, p, pr=0)

Calculates potential temperature as per UNESCO 1983 report.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature [℃ (ITS-90)]

parray_like

pressure [db].

prarray_like

reference pressure [db], default = 0

Returns:
ptarray_like

potential temperature relative to PR [℃ (ITS-90)]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

[2]

Bryden, H. 1973. New Polynomials for thermal expansion, adiabatic temperature gradient and potential temperature of sea water. Deep-Sea Res. Vol20,401-408. doi:10.1016/0011-7471(73)90063-6

Examples

>>> import seawater as sw
>>> from seawater.library import T90conv, T68conv
>>> t = T90conv(
...     [
...         [0, 0, 0, 0, 0, 0],
...         [10, 10, 10, 10, 10, 10],
...         [20, 20, 20, 20, 20, 20],
...         [30, 30, 30, 30, 30, 30],
...         [40, 40, 40, 40, 40, 40],
...     ]
... )
>>> s = [
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
... ]
>>> p = [0, 5000, 10000, 0, 5000, 10000]
>>> T68conv(sw.ptmp(s, t, p, pr=0))
array([[ 0.        , -0.30614418, -0.96669485,  0.        , -0.3855565 ,
        -1.09741136],
       [10.        ,  9.35306331,  8.46840949, 10.        ,  9.29063461,
         8.36425752],
       [20.        , 19.04376281, 17.94265   , 20.        , 18.99845171,
        17.86536441],
       [30.        , 28.75124632, 27.43529911, 30.        , 28.72313484,
        27.38506197],
       [40.        , 38.46068173, 36.92544552, 40.        , 38.44979906,
        36.90231661]])
seawater.eos80.salt(r, t, p)

Calculates Salinity from conductivity ratio. UNESCO 1983 polynomial.

Parameters:
rarray_like

conductivity ratio \(R = \\frac{C(S,T,P)}{C(35,15(IPTS-68),0)}\)

tarray_like

temperature [℃ (ITS-90)]

parray_like

pressure [db]

Returns:
sarray_like

salinity [psu (PSS-78)]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

Examples

Data from UNESCO 1983 p9. >>> import seawater as sw >>> from seawater.library import T90conv >>> r = [1, 1.2, 0.65] >>> t = T90conv([15, 20, 5]) >>> p = [0, 2000, 1500] >>> sw.salt(r, t, p) array([34.99999992, 37.24562765, 27.99534693])

seawater.eos80.svel(s, t, p)

Sound Velocity in sea water using UNESCO 1983 polynomial.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

t(p)array_like

temperature [℃ (ITS-90)]

parray_like

pressure [db].

Returns:
svelarray_like

sound velocity [m/s]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

Examples

Data from Pond and Pickard Intro. Dynamical Oceanography 2nd ed. 1986

>>> import seawater as sw
>>> from seawater.library import T90conv
>>> t = T90conv(
...     [
...         [0, 0, 0, 0, 0, 0],
...         [10, 10, 10, 10, 10, 10],
...         [20, 20, 20, 20, 20, 20],
...         [30, 30, 30, 30, 30, 30],
...         [40, 40, 40, 40, 40, 40],
...     ]
... )
>>> s = [
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
...     [25, 25, 25, 35, 35, 35],
... ]
>>> p = [0, 5000, 10000, 0, 5000, 10000]
>>> sw.svel(s, t, p)
array([[1435.789875  , 1520.358725  , 1610.4074    , 1449.13882813,
        1533.96863705, 1623.15007097],
       [1477.68316464, 1561.30635914, 1647.39267114, 1489.82233602,
        1573.40946928, 1658.99115504],
       [1510.31388348, 1593.59671798, 1676.80967748, 1521.4619731 ,
        1604.4762822 , 1687.18305631],
       [1535.21434752, 1618.95631952, 1700.60547902, 1545.59485539,
        1628.97322783, 1710.06294277],
       [1553.44506636, 1638.02522336, 1719.15088536, 1563.20925247,
        1647.29949576, 1727.83176404]])
seawater.eos80.temp(s, pt, p, pr=0)

Calculates temperature from potential temperature at the reference pressure PR and in situ pressure P.

Parameters:
s(p)array_like

salinity [psu (PSS-78)]

pt(p)array_like

potential temperature [℃ (ITS-90)]

parray_like

pressure [db].

prarray_like

reference pressure [db]

Returns:
temparray_like

temperature [℃ (ITS-90)]

References

[1]

Fofonoff, P. and Millard, R.C. Jr UNESCO 1983. Algorithms for computation of fundamental properties of seawater. UNESCO Tech. Pap. in Mar. Sci., No. 44, 53 pp. Eqn.(31) p.39. https://unesdoc.unesco.org/ark:/48223/pf0000059832_eng

[2]

Bryden, H. 1973. New Polynomials for thermal expansion, adiabatic temperature gradient and potential temperature of sea water. Deep-Sea Res. Vol20,401-408. doi:10.1016/0011-7471(73)90063-6

Examples

>>> import seawater as sw
>>> s, t, p = 35, 15, 100
>>> sw.temp(s, sw.ptmp(s, t, p), p)
15.0