密度は海洋の動きに密接している | 海洋観測

2023-11-192023-11-23

海洋観測において塩分の測定が最も基本で、重要な要素となります。それはなぜか。塩分が分かれば密度を求めることができるからです。海洋には熱延循環という、重い（密度が高い）水が深層まで沈み込む大きな流れがあります。また海水は同じ密度面に沿って流れる性質もあるため、海水の密度を測定することは海洋の流れを追ううえでとても重要です。

密度の計算式

本計算式の定義は古く、現在の温度の定義ITS-90とは異なるIPTS-68の時代のものです。このため本計算を行うためには、ITS-90とIPTS-68の温度換算が必要です。本ページ内において「$T$」はITS-90の温度を、「$t$」はIPTS-68の温度を表します。
1MPa = 100dbar ≒ 100m

あわせて読みたい

水温の定義 | 海洋観測 ITS-90とIPTS-68 海水の温度は圧力、塩分とともに海水の状態を調べるための最も基本的な要素です。温度の定義は実は何度か変わっています。現在温度というとITS-90（In…

水温を$T$ [℃]、圧力を$P$ [dbar]、実用塩分を$S$とすると、海水の密度$\rho$ [kg/m³]は下記式¹により求められます。

$$ \begin{aligned}
\rho(S,t,p) &= \frac{ \rho(S, t, 0) }{ 1-\displaystyle\frac{ p }{ K(S, t, p) } } \\
t &= 1.00024 \cdot T \\
p &= P / 10 \\
\rho(S,t,0) &= \rho_w + \left( b_{0} + b_{1} \cdot t + b_{2} \cdot t^{2} + b_{3} \cdot t^{3} + b_{4} \cdot t^{4} \right) \cdot S + \left( c_{0} + c_{1} \cdot t + c_{2} \cdot t^{2} \right) \cdot S^{3/2} + d_{0} \cdot S^{2} \\
\rho_{w} &= a_{0} + a_{1} \cdot t + a_{2} \cdot t^{2} + a_{3} \cdot t^{3} + a_{4} \cdot t^{4} + a_{5}
\cdot t^{5} \\
K(S,t,p) &= K(S, t, 0) + A \cdot p + B \cdot p^{2} \\
K(S,t,0) &= K_{w} + \left( f_{0} + f_{1} \cdot t + f_{2} \cdot t^{2} + f_{3} \cdot t^{3} \right) \cdot S + \left( g_{0} + g_{1} \cdot t + g_{2} \cdot t^{2} \right) \cdot S^{3/2} \\
A &= A_{w} + \left( i_{0} + i_{1} \cdot t + i_{2} \cdot t^{2} \right) \cdot S + j_{0} \cdot S^{3/2} \\
B &= B_{w} + \left( m_{0} + m_{1} \cdot t + m_{2} \cdot t^{2} \right) \cdot S \\
K_{w} &= e_{0} + e_{1} \cdot t + e_{2} \cdot t^{2} + e_{3} \cdot t^{3} + e_{4} \cdot t^{4} \\
A_{w} &= h_{0} + h_{1} \cdot t + h_{2} \cdot t^{2} + h_{3} \cdot t^{3} \\
B_{w} &= k_{0} + k_{1} \cdot t + k_{2} \cdot t^{2}
\end{aligned} $$

上記計算式で使用する係数は下記のとおりです。

$a_{0} = 999.842594$、$a_{1} = 6.793952 \times 10^{-2}$、$a_{2} = -9.095290 \times 10^{-3}$、$a_{3} = 1.001685 \times 10^{-4}$、$a_{4} = -1.120083 \times 10^{-6}$、$a_{5} = 6.536332 \times 10^{-9}$
$b_{0} = 8.24493 \times 10^{-1}$、$b_{1} = -4.0899 \times 10^{-3}$、$b_{2} = 7.6438 \times 10^{-5}$、$b_{3} = -8.2467 \times 10^{-7}$、$b_{4} = 5.3875 \times 10^{-9}$
$c_{0} = -5.72466 \times 10^{-3}$、$c_{1} = 1.0227 \times 10^{-4}$、$c_{2} = -1.6546 \times 10^{-6}$
$d_{0} = 4.8314 \times 10^{-4}$
$e_{0} = 19652.21$、$e_{1} = 148.4206$、$e_{2} = -2.327105$、$e_{3} = 1.360477 \times 10^{-2}$、$e_{4} = -5.155288 \times 10^{-5}$
$f_{0} = 54.6746$、$f_{1} = -0.603459$、$f_{2} = 1.09987 \times 10^{-2}$、$f_{3} = -6.1670 \times 10^{-5}$
$g_{0} = 7.944 \times 10^{-2}$、$g_{1} = 1.6483 \times 10^{-2}$、$g_{2} = -5.3009 \times 10^{-4}$
$h_{0} = 3.239908$、$h_{1} = 1.43713 \times 10^{-3}$、$h_{2} = 1.16092 \times 10^{-4}$、$h_{3} = -5.77905 \times 10^{-7}$
$i_{0} = 2.2838 \times 10^{-3}$、$i_{1} = -1.0981 \times 10^{-5}$、$i_{2} = -1.6078 \times 10^{-6}$
$j_{0} = 1.91075 \times 10^{-4}$
$k_{0} = 8.50935 \times 10^{-5}$、$k_{1} = -6.12293 \times 10^{-6}$、$k_{2} = 5.2787 \times 10^{-8}$
$m_{0} = -9.9348 \times 10^{-7}$、$m_{1} = 2.0816 \times 10^{-8}$、$m_{2} = 9.1697 \times 10^{-10}$

なお海水の密度は、ほとんどが1,000～1,030kg/m³の間です。したがって重要なのは下2桁の変化となります。このため海洋では、$\sigma_{t}$を密度として表すこともあります。

$\sigma_{t} = \rho(S, t, 0)-1000$

なお上記式は実用塩分$S$が0～42、水温$t$が-2～40℃、圧力$P$が0～10,000dbarの範囲で有効です。

計算例

スクロールできます

$S$	$t$ [℃]	$P$ [dbar]	$\rho$ [kg/m³]
0	5	0	999.96675
0	5	10000	1044.12802
35	25	0	1023.34306
35	25	10000	1062.53817

サンプルコード

[density.c]

/**
 * @file density.c
 * @brief 密度計算
 * @brief 1MPa = 100dbar
 */

//------------------------------------------------------------------------------
// include
//------------------------------------------------------------------------------
#include <float.h>
#include <math.h>

#include "density.h"

/** ----------------------------------------------------------------------------
 * @brief 密度の計算
 * @details UNESCO (1983) "Specific Volume Anomaly and Density Anomaly of Seawater"
 *
 * @param temp : 水温[℃] (ITS-90) [-2 - 40 (IPTS-68)]
 * @param sal : 塩分 [0 - 42]
 * @param press : 圧力[dbar] [0 - 10000]
 * @param [out] *sigma_t : σt
 * @return 密度[kg/m3]
 */
double Calc_Density(double temp, double sal, double press, double *sigma_t)
{
  static const double a[6] = {999.842594, 6.793952e-2, -9.095290e-3, 1.001685e-4, -1.120083e-6, 6.536332e-9};
  static const double b[5] = {8.24493e-1, -4.0899e-3, 7.6438e-5, -8.2467e-7, 5.3875e-9};
  static const double c[3] = {-5.72466e-3, 1.0227e-4, -1.6546e-6};
  static const double d0   = 4.8314e-4;
  static const double e[5] = {19652.21, 148.4206, -2.327105, 1.360477e-2, -5.155288e-5};
  static const double f[4] = {54.6746, -0.603459, 1.09987e-2, -6.1670e-5};
  static const double g[3] = {7.944e-2, 1.6483e-2, -5.3009e-4};
  static const double h[4] = {3.239908, 1.43713e-3, 1.16092e-4, -5.77905e-7};
  static const double i[3] = {2.2838e-3, -1.0981e-5, -1.6078e-6};
  static const double j0   = 1.91075e-4;
  static const double k[3] = {8.50935e-5, -6.12293e-6, 5.2787e-8};
  static const double m[3] = {-9.9348e-7, 2.0816e-8, 9.1697e-10};

  double rho, rho_0, rho_w, K, K_0, K_w, A, A_w, B, B_w;
  double t, t2, t3, t4, t5;
  double density, sal15;
  double p, p2;
  double aa, bb, cc, ee, ff, gg, hh, ii, kk, mm;

  *sigma_t = 0.0;

  t  = temp * 1.00024;  // ITS-90 ⇒ IPTS-68
  t2 = t * t;
  t3 = t2 * t;
  t4 = t3 * t;
  t5 = t4 * t;

  if (sal < DBL_EPSILON) {
    sal = 0.0;
  }
  sal15 = sqrt(sal) * sal;

  p = press;
  if (p < DBL_EPSILON) {
    p = 0.0;
  }
  p /= 10.0;  // dbar ⇒ bar
  p2 = p * p;

  aa = a[0] + a[1] * t + a[2] * t2 + a[3] * t3 + a[4] * t4 + a[5] * t5;
  bb = b[0] + b[1] * t + b[2] * t2 + b[3] * t3 + b[4] * t4;
  cc = c[0] + c[1] * t + c[2] * t2;
  ee = e[0] + e[1] * t + e[2] * t2 + e[3] * t3 + e[4] * t4;
  ff = f[0] + f[1] * t + f[2] * t2 + f[3] * t3;
  gg = g[0] + g[1] * t + g[2] * t2;
  hh = h[0] + h[1] * t + h[2] * t2 + h[3] * t3;
  kk = k[0] + k[1] * t + k[2] * t2;
  ii = i[0] + i[1] * t + i[2] * t2;
  mm = m[0] + m[1] * t + m[2] * t2;

  K_w = ee;
  A_w = hh;
  B_w = kk;
  K_0 = K_w + ff * sal + gg * sal15;

  A = A_w + ii * sal + j0 * sal15;
  B = B_w + mm * sal;
  K = K_0 + A * p + B * p2;

  rho_w = aa;
  rho_0 = rho_w + bb * sal + cc * sal15 + d0 * sal * sal;

  if (fabs(K) < DBL_EPSILON) {
    rho = 1.0;
  } else {
    rho = 1.0 - p / K;
  }
  if (fabs(rho) < DBL_EPSILON) {
    return 0.0;
  }

  density = rho_0 / rho;
  if (density < DBL_EPSILON) {
    density = 0.0;
  }
  if (rho_0 > (1000.0 + DBL_EPSILON)) {
    *sigma_t = rho_0 - 1000.0;
  } else {
    *sigma_t = 0.0;
  }

  return density;
}

[density.h]

#ifndef DENSITY_H
#define DENSITY_H

#ifdef __cplusplus
extern "C" {
#endif

//------------------------------------------------------------------------------
// function
//------------------------------------------------------------------------------
double Calc_Density(double temp, double sal, double press, double *sigma_t);

#ifdef __cplusplus
}
#endif

#endif

サンプルコードの使用例

使用例代わりに、CppUTestによるテストコードを掲載しておきます。

#include "CppUTest/TestHarness.h"

#include "density.h"

// clang-format off
TEST_GROUP(density){
  TEST_SETUP(){
  }

  TEST_TEARDOWN(){
  }
};
// clang-format on

TEST(density, Test_density)
{
  double den, sigma_t;

  den = Calc_Density(5.0 / 1.00024, 0.0, 0.0, &sigma_t);
  DOUBLES_EQUAL(den, 999.96675, 5e-6);

  den = Calc_Density(5.0 / 1.00024, -10.0, -10.0, &sigma_t);
  DOUBLES_EQUAL(den, 999.96675, 5e-6);

  den = Calc_Density(5.0 / 1.00024, 0.0, 10000.0, &sigma_t);
  DOUBLES_EQUAL(den, 1044.12802, 5e-6);

  den = Calc_Density(25.0 / 1.00024, 35.0, 0.0, &sigma_t);
  DOUBLES_EQUAL(den, 1023.34306, 5e-6);

  den = Calc_Density(25.0 / 1.00024, 35.0, 10000.0, &sigma_t);
  DOUBLES_EQUAL(den, 1062.53817, 5e-6);
}

脚注

UNESCO (1983) “Specific Volume Anomaly and Density Anomaly of Seawater” Unesco technical papers in marine science No.44 ↩︎