メイン(ドライバールーチン)
sigma_theta.f90 ダウンロード
下請けルーチン
theta.f (ポテンシャル温度) ダウンロード
atg.f (adiabatic temperature gradient) ダウンロード
eqst2.f (海水の状態方程式) ダウンロード
theta.f90 (ポテンシャル温度)
function theta(s,t0,p0,pr)
!
! ***********************************
! to compute local potential temperature at pr
! using bryden 1973 polynomial for adiabatic lapse rate
! and runge-kutta 4-th order integration algorithm.
! ref: bryden,h.,1973,deep-sea res.,20,401-408
! fofonoff,n.,1977,deep-sea res.,24,489-491
! units:
! pressure p0 decibars
! temperature t0 deg celsius (ipts-68)
! salinity s (ipss-78)
! reference prs pr decibars
! potential tmp. theta deg celsius
! checkvalue: theta= 36.89073 c,s=40 (ipss-78),t0=40 deg c,
! p0=10000 decibars,pr=0 decibars
!
! set-up intermediate temperature and pressure variables
!
! http://www.marine.csiro.au/~jackett/NeutralDensity/
real p,p0,t,t0,h,pr,xk,s,q,theta,atg
p=p0
t=t0
h = pr - p
xk = h*atg(s,t,p)
t = t + 0.5*xk
q = xk
p = p + 0.5*h
xk = h*atg(s,t,p)
t = t + 0.29289322*(xk-q)
q = 0.58578644*xk + 0.121320344*q
xk = h*atg(s,t,p)
t = t + 1.707106781*(xk-q)
q = 3.414213562*xk - 4.121320344*q
p = p + 0.5*h
xk = h*atg(s,t,p)
theta = t + (xk-2.0*q)/6.0
return
end function theta
atg.f90
function atg(s,t,p)
!
! ****************************
! adiabatic temperature gradient deg c per decibar
! ref: bryden,h.,1973,deep-sea res.,20,401-408
! units:
! pressure p decibars
! temperature t deg celsius (ipts-68)
! salinity s (ipss-78)
! adiabatic atg deg. c/decibar
! checkvalue: atg=3.255976e-4 c/dbar for s=40 (ipss-78),
! t=40 deg c,p0=10000 decibars
!
! http://www.marine.csiro.au/~jackett/NeutralDensity/
ds = s - 35.0
atg = (((-2.1687d-16*t+1.8676d-14)*t-4.6206d-13)*p &
+ ((2.7759d-12*t-1.1351d-10)*ds+((-5.4481d-14*t &
+ 8.733d-12)*t-6.7795d-10)*t+1.8741d-8))*p &
+ (-4.2393d-8*t+1.8932d-6)*ds &
+ ((6.6228d-10*t-6.836d-8)*t+8.5258d-6)*t+3.5803d-5
return
end function atg