Go to the documentation of this file.
13 implicit none ;
private
36 real,
parameter ::
r00 = 999.842594,
r10 = 6.793952e-2,
r20 = -9.095290e-3, &
37 r30 = 1.001685e-4,
r40 = -1.120083e-6,
r50 = 6.536332e-9,
r01 = 0.824493, &
38 r11 = -4.0899e-3,
r21 = 7.6438e-5,
r31 = -8.2467e-7,
r41 = 5.3875e-9, &
39 r032 = -5.72466e-3,
r132 = 1.0227e-4,
r232 = -1.6546e-6,
r02 = 4.8314e-4
45 real,
parameter ::
s00 = 1.965933e4,
s10 = 1.444304e2,
s20 = -1.706103, &
46 s30 = 9.648704e-3,
s40 = -4.190253e-5,
s01 = 52.84855,
s11 = -3.101089e-1, &
47 s21 = 6.283263e-3,
s31 = -5.084188e-5,
s032 = 3.886640e-1,
s132 = 9.085835e-3, &
48 s232 = -4.619924e-4,
sp00 = 3.186519,
sp10 = 2.212276e-2,
sp20 = -2.984642e-4, &
49 sp30 = 1.956415e-6,
sp01 = 6.704388e-3,
sp11 = -1.847318e-4,
sp21 = 2.059331e-7, &
62 real,
intent(in) :: pressure
63 real,
intent(out) :: rho
64 real,
optional,
intent(in) :: rho_ref
67 real,
dimension(1) :: t0, s0, pressure0
68 real,
dimension(1) :: rho0
72 pressure0(1) = pressure
83 real,
dimension(:),
intent(in) :: t
84 real,
dimension(:),
intent(in) :: s
85 real,
dimension(:),
intent(in) :: pressure
86 real,
dimension(:),
intent(out) :: rho
87 integer,
intent(in) :: start
88 integer,
intent(in) :: npts
89 real,
optional,
intent(in) :: rho_ref
92 real :: t_local, t2, t3, t4, t5
93 real :: s_local, s32, s2
100 do j=start,start+npts-1
101 if (s(j) < -1.0e-10)
then
106 p1 = pressure(j)*1.0e-5; p2 = p1*p1
107 t_local = t(j); t2 = t_local*t_local; t3 = t_local*t2; t4 = t2*t2; t5 = t3*t2
108 s_local = s(j); s2 = s_local*s_local; s32 = s_local*sqrt(s_local)
125 if (
present(rho_ref))
then
126 rho(j) = ((
r00 - rho_ref)*ks + (sig0*ks + p1*rho_ref)) / (ks - p1)
128 rho(j) = rho0*ks / (ks - p1)
138 real,
intent(in) :: T
140 real,
intent(in) :: S
141 real,
intent(in) :: pressure
142 real,
intent(out) :: specvol
143 real,
optional,
intent(in) :: spv_ref
146 real,
dimension(1) :: T0, S0, pressure0, spv0
148 t0(1) = t ; s0(1) = s ; pressure0(1) = pressure
159 real,
dimension(:),
intent(in) :: T
161 real,
dimension(:),
intent(in) :: S
162 real,
dimension(:),
intent(in) :: pressure
163 real,
dimension(:),
intent(out) :: specvol
164 integer,
intent(in) :: start
165 integer,
intent(in) :: npts
166 real,
optional,
intent(in) :: spv_ref
169 real :: t_local, t2, t3, t4, t5
170 real :: s_local, s32, s2
176 do j=start,start+npts-1
177 if (s(j) < -1.0e-10)
then
179 if (
present(spv_ref)) specvol(j) = 0.001 - spv_ref
183 p1 = pressure(j)*1.0e-5; p2 = p1*p1
184 t_local = t(j); t2 = t_local*t_local; t3 = t_local*t2; t4 = t2*t2; t5 = t3*t2
185 s_local = s(j); s2 = s_local*s_local; s32 = s_local*sqrt(s_local)
201 if (
present(spv_ref))
then
202 specvol(j) = (ks*(1.0 - (rho0*spv_ref)) - p1) / (rho0*ks)
204 specvol(j) = (ks - p1) / (rho0*ks)
213 real,
intent(in),
dimension(:) :: t
215 real,
intent(in),
dimension(:) :: s
216 real,
intent(in),
dimension(:) :: pressure
217 real,
intent(out),
dimension(:) :: drho_dt
219 real,
intent(out),
dimension(:) :: drho_ds
221 integer,
intent(in) :: start
222 integer,
intent(in) :: npts
225 real :: t_local, t2, t3, t4, t5
226 real :: s12, s_local, s32, s2
237 do j=start,start+npts-1
238 if (s(j) < -1.0e-10)
then
239 drho_dt(j) = 0.0 ; drho_ds(j) = 0.0
243 p1 = pressure(j)*1.0e-5; p2 = p1*p1
244 t_local = t(j); t2 = t_local*t_local; t3 = t_local*t2; t4 = t2*t2; t5 = t3*t2
245 s_local = s(j); s2 = s_local*s_local; s12 = sqrt(s_local); s32 = s_local*s12
252 drho0_dt =
r10 + 2.0*
r20*t_local + 3.0*
r30*t2 + 4.0*
r40*t3 + 5.0*
r50*t4 + &
253 s_local*(
r11 + 2.0*
r21*t_local + 3.0*
r31*t2 + 4.0*
r41*t3) + &
265 dks_dt =
s10 + 2.0*
s20*t_local + 3.0*
s30*t2 + 4.0*
s40*t3 + &
273 denom = 1.0 / (ks - p1)
274 drho_dt(j) = denom*(ks*drho0_dt - rho0*p1*denom*dks_dt)
275 drho_ds(j) = denom*(ks*drho0_ds - rho0*p1*denom*dks_ds)
284 real,
intent(in),
dimension(:) :: t
286 real,
intent(in),
dimension(:) :: s
287 real,
intent(in),
dimension(:) :: pressure
288 real,
intent(out),
dimension(:) :: rho
289 real,
intent(out),
dimension(:) :: drho_dp
292 integer,
intent(in) :: start
293 integer,
intent(in) :: npts
296 real :: t_local, t2, t3, t4, t5
297 real :: s_local, s32, s2
301 real :: ks_0, ks_1, ks_2
306 do j=start,start+npts-1
307 if (s(j) < -1.0e-10)
then
308 rho(j) = 1000.0 ; drho_dp(j) = 0.0
312 p1 = pressure(j)*1.0e-5; p2 = p1*p1
313 t_local = t(j); t2 = t_local*t_local; t3 = t_local*t2; t4 = t2*t2; t5 = t3*t2
314 s_local = s(j); s2 = s_local*s_local; s32 = s_local*sqrt(s_local)
329 ks = ks_0 + p1*ks_1 + p2*ks_2
330 dks_dp = ks_1 + 2.0*p1*ks_2
332 rho(j) = rho0*ks / (ks - p1)
334 drho_dp(j) = 1.0e-5 * (rho(j) / (ks - p1)) * (1.0 - dks_dp*p1/ks)
Compute the in situ specific volume of sea water (in [m3 kg-1]), or an anomaly with respect to a refe...
real, parameter sp011
Parameters in the UNESCO equation of state.
real, parameter sp000
Parameters in the UNESCO equation of state.
The equation of state using the Jackett and McDougall fits to the UNESCO EOS.
real, parameter s01
Parameters in the UNESCO equation of state.
real, parameter s31
Parameters in the UNESCO equation of state.
real, parameter r40
Parameters in the UNESCO equation of state.
real, parameter r50
Parameters in the UNESCO equation of state.
real, parameter s30
Parameters in the UNESCO equation of state.
real, parameter r032
Parameters in the UNESCO equation of state.
real, parameter s032
Parameters in the UNESCO equation of state.
real, parameter r20
Parameters in the UNESCO equation of state.
real, parameter r41
Parameters in the UNESCO equation of state.
real, parameter sp10
Parameters in the UNESCO equation of state.
real, parameter r21
Parameters in the UNESCO equation of state.
subroutine calculate_spec_vol_array_unesco(T, S, pressure, specvol, start, npts, spv_ref)
This subroutine computes the in situ specific volume of sea water (specvol in [m3 kg-1]) from salinit...
real, parameter r31
Parameters in the UNESCO equation of state.
real, parameter s00
Parameters in the UNESCO equation of state.
real, parameter sp21
Parameters in the UNESCO equation of state.
real, parameter r11
Parameters in the UNESCO equation of state.
real, parameter s20
Parameters in the UNESCO equation of state.
real, parameter r10
Parameters in the UNESCO equation of state.
real, parameter r00
Parameters in the UNESCO equation of state.
real, parameter sp11
Parameters in the UNESCO equation of state.
real, parameter r232
Parameters in the UNESCO equation of state.
real, parameter r02
Parameters in the UNESCO equation of state.
real, parameter r30
Parameters in the UNESCO equation of state.
real, parameter s132
Parameters in the UNESCO equation of state.
real, parameter sp021
Parameters in the UNESCO equation of state.
real, parameter sp001
Parameters in the UNESCO equation of state.
real, parameter sp010
Parameters in the UNESCO equation of state.
real, parameter sp020
Parameters in the UNESCO equation of state.
subroutine, public calculate_compress_unesco(T, S, pressure, rho, drho_dp, start, npts)
This subroutine computes the in situ density of sea water (rho) and the compressibility (drho/dp == C...
real, parameter s40
Parameters in the UNESCO equation of state.
real, parameter sp032
Parameters in the UNESCO equation of state.
real, parameter sp01
Parameters in the UNESCO equation of state.
real, parameter r132
Parameters in the UNESCO equation of state.
real, parameter sp00
Parameters in the UNESCO equation of state.
real, parameter r01
Parameters in the UNESCO equation of state.
subroutine, public calculate_density_scalar_unesco(T, S, pressure, rho, rho_ref)
This subroutine computes the in situ density of sea water (rho in [kg m-3]) from salinity (S [PSU]),...
real, parameter s11
Parameters in the UNESCO equation of state.
real, parameter s232
Parameters in the UNESCO equation of state.
real, parameter sp30
Parameters in the UNESCO equation of state.
real, parameter sp20
Parameters in the UNESCO equation of state.
Compute the in situ density of sea water (in [kg m-3]), or its anomaly with respect to a reference de...
subroutine calculate_spec_vol_scalar_unesco(T, S, pressure, specvol, spv_ref)
This subroutine computes the in situ specific volume of sea water (specvol in [m3 kg-1]) from salinit...
real, parameter s10
Parameters in the UNESCO equation of state.
subroutine, public calculate_density_derivs_unesco(T, S, pressure, drho_dT, drho_dS, start, npts)
This subroutine calculates the partial derivatives of density with potential temperature and salinity...
real, parameter s21
Parameters in the UNESCO equation of state.
subroutine, public calculate_density_array_unesco(T, S, pressure, rho, start, npts, rho_ref)
This subroutine computes the in situ density of sea water (rho in [kg m-3]) from salinity (S [PSU]),...