Browse Source

Merge branch 'time_microseconds' into iow_master

master
Knut 1 year ago
parent
commit
0e2be6621b
  1. 10
      doc/getm_pub.bib
  2. 9
      schemas/getm-2.5.defaults
  3. 3
      schemas/getm-2.5.schema
  4. 2
      src/2d/adv_split_u.F90
  5. 2
      src/2d/adv_split_v.F90
  6. 10
      src/2d/advection.F90
  7. 15
      src/2d/variables_2d.F90
  8. 41
      src/3d/advection_3d.F90
  9. 2
      src/3d/getm_fabm.F90
  10. 2
      src/3d/m3d.F90
  11. 76
      src/3d/rivers.F90
  12. 334
      src/3d/salinity.F90
  13. 5
      src/3d/tracer_diffusion.F90
  14. 257
      src/3d/variables_3d.F90
  15. 1
      src/CMakeLists.txt
  16. 3
      src/domain/domain.F90
  17. 1
      src/futils/Makefile
  18. 2
      src/futils/getm_timers.F90
  19. 71
      src/futils/rebin.F90
  20. 3
      src/getm/compilation_options.F90
  21. 11
      src/getm/initialise.F90
  22. 16
      src/getm/integration.F90
  23. 10
      src/getm/register_all_variables.F90
  24. 4
      src/ncdf/ncdf_rivers.F90
  25. 306
      src/output/output_processing.F90

10
doc/getm_pub.bib

@ -42,6 +42,16 @@
doi = {10.1029/2020JC016411},
}
@Article{fofonova2021pst,
author = {Fofonova, Vera and K\"arn\"a, Tuomas and Klingbeil, Knut and Androsov, Alexey and Kuznetsov, Ivan and Sidorenko, Dmitry and Danilov, Sergey and Burchard, Hans and Wiltshire, Karen Helen},
journal = {Geoscientific Model Development},
title = {{Plume spreading test case for coastal ocean models}},
year = {2021},
pages = {6945--6975},
volume = {14},
doi = {10.5194/gmd-14-6945-2021},
}
@article{Friedland2021,
author = {Friedland, Ren{\'{e}} and Macias, Diego and Cossarini, Gianpiero and Daewel, Ute and Estournel, Claude and Garcia-Gorriz, Elisa and Gr{\'{e}}goire, Marilaure and Gustafson, Bo and Grizzetti, Bruna and Kalaroni, Sofia and Kerimoglu, Onur and Lazzari, Paolo and Lenhart, Hermann and Lessin, Gennadi and Maljutenko, Ilja and Miladinova, Svetla and M{\"{u}}ller-Karulis, B{\"{a}}rbel and Neumann, Thomas and Parn, Ove and P{\"{a}}tsch, Johannes and Piroddi, Chiara and Raudsepp, Urmas and Schrum, Corinna and Stegert, Christoph and Stips, Adolf and Tsiaras, Kostas and Ulses, Caroline and Vandenbulcke, Luc},
doi = {10.3389/fmars.2021.596126},

9
schemas/getm-2.5.defaults

@ -620,6 +620,15 @@
<sss_nudging_time>
-1.0
</sss_nudging_time>
<saltbins>
0
</saltbins>
<saltbinmin>
0.0
</saltbinmin>
<saltbinmax>
40.0
</saltbinmax>
<salt_check>
0
</salt_check>

3
schemas/getm-2.5.schema

@ -990,6 +990,9 @@
<condition type="eq" variable="./salt_AH_method" value="3"/>
</element>
<element name="sss_nudging_time" type="float" minExclusive="0.0" label="nudging time for sss read from files given in ssts_files.dat" unit="s"/>
<element name="saltbins" type="int" label="number of salinity bins for TEF"/>
<element name="saltbinmin" type="float" label="minimum salinity for TEF"/>
<element name="saltbinmax" type="float" label="maximum salinity for TEF"/>
<element name="salt_check" type="int" label="interval between checks for out-of-bounds salinities" description="Interval between checks for out-of-bounds salinities. 0: disabled, &gt;0: check every salt_check step, abort if check fails, &lt;0: check every abs(salt_check) step, warn if check fails"/>
<element name="min_salt" type="float" label="minimum valid salinity" unit="psu">
<condition type="ne" variable="./salt_check" value="0"/>

2
src/2d/adv_split_u.F90

@ -151,6 +151,7 @@
! !USES:
use domain, only: imin,imax,jmin,jmax
use util, only: adv_reconstruct
use advection, only: ffluxu2
!$ use omp_lib
IMPLICIT NONE
!
@ -273,6 +274,7 @@
if (associated(ffluxu)) then
ffluxu(_U2DFIELD_) = ffluxu(_U2DFIELD_) + splitfac*uflux(_U2DFIELD_)
end if
ffluxu2(_U2DFIELD_) = ffluxu2(_U2DFIELD_) + splitfac*uflux2(_U2DFIELD_)
fio = fi
Dio = Di
!$OMP END WORKSHARE

2
src/2d/adv_split_v.F90

@ -21,6 +21,7 @@
! !USES:
use domain, only: imin,imax,jmin,jmax
use util, only: adv_reconstruct
use advection, only: ffluxv2
!$ use omp_lib
IMPLICIT NONE
!
@ -127,6 +128,7 @@
if (associated(ffluxv)) then
ffluxv(_V2DFIELD_) = ffluxv(_V2DFIELD_) + splitfac*vflux(_V2DFIELD_)
end if
ffluxv2(_V2DFIELD_) = ffluxv2(_V2DFIELD_) + splitfac*vflux2(_V2DFIELD_)
fio = fi
Dio = Di
!$OMP END WORKSHARE

10
src/2d/advection.F90

@ -82,10 +82,12 @@
logical,dimension(E2DFIELD),target :: mask_updateH
logical,dimension(E2DFIELD),target :: mask_uflux,mask_vflux,mask_xflux
logical,dimension(E2DFIELD),target :: mask_uupdateU,mask_vupdateV
REALTYPE,dimension(E2DFIELD),public :: ffluxu2,ffluxv2
#else
logical,dimension(:,:),allocatable,target :: mask_updateH
logical,dimension(:,:),allocatable,target :: mask_uflux,mask_vflux,mask_xflux
logical,dimension(:,:),allocatable,target :: mask_uupdateU,mask_vupdateV
REALTYPE,dimension(:,:),allocatable,public :: ffluxu2,ffluxv2
#endif
!
! !REVISION HISTORY:
@ -234,6 +236,12 @@
allocate(mask_vupdateV(E2DFIELD),stat=rc) ! work array
if (rc /= 0) stop 'init_advection: Error allocating memory (mask_vupdateV)'
allocate(ffluxu2(E2DFIELD),stat=rc) ! work array
if (rc /= 0) stop 'init_advection: Error allocating memory (ffluxu2)'
allocate(ffluxv2(E2DFIELD),stat=rc) ! work array
if (rc /= 0) stop 'init_advection: Error allocating memory (ffluxv2)'
#endif
! Note (KK): In this module pointers are used extensively.
@ -468,6 +476,8 @@
end if
end if
ffluxu2 = _ZERO_ ; ffluxv2 = _ZERO_
if (present(Dires)) then
p_Di => Dires
else

15
src/2d/variables_2d.F90

@ -34,6 +34,7 @@
logical :: calc_taubmax=.false.
logical :: save_numdis_2d=.false.
logical :: save_phydis_2d=.false.
logical :: save_rvol=.false.
#ifdef STATIC
#include "static_2d.h"
@ -61,6 +62,8 @@
REALTYPE,dimension(:,:),contiguous,pointer :: numdis_2d=>null()
REALTYPE,dimension(:,:),contiguous,pointer :: phydis_2d=>null()
REALTYPE,dimension(:,:),allocatable :: rvol
integer :: size2d_field
integer :: mem2d
!
@ -317,6 +320,12 @@
taubmax = _ZERO_
end if
if (save_rvol) then
allocate(rvol(E2DFIELD),stat=rc)
if (rc /= 0) stop 'postinit_2d: Error allocating memory (rvol)'
rvol = _ZERO_
end if
end if
#ifdef DEBUG
@ -340,6 +349,7 @@
! !USES:
use field_manager
use meteo, only: fwf_method
use domain, only: nriverl
IMPLICIT NONE
!
! !INPUT PARAMETERS:
@ -412,6 +422,10 @@
call fm%register('fwf', 'm', 'surface freshwater fluxes', standard_name='', data2d=fwf(_2D_FO_), category='2d', add_to_mean_before_save=.true., output_level=output_level_debug)
end if
if (nriverl .gt. 0) then
call fm%register('rvol', 'm**3', 'river discharge', category='2d', fill_value=_ZERO_, output_level=output_level_debug, used=save_rvol)
end if
return
end subroutine register_2d_variables
!EOC
@ -445,6 +459,7 @@
if (associated(taubmax)) call fm%send_data('taubmax', taubmax(_2D_FO_))
if (associated(numdis_2d)) call fm%send_data('numdis_2d', numdis_2d(_2D_FO_))
if (associated(phydis_2d)) call fm%send_data('phydis_2d', phydis_2d(_2D_FO_))
if (allocated(rvol)) call fm%send_data('rvol', rvol(_2D_FO_))
return
end subroutine finalize_register_2d_variables

41
src/3d/advection_3d.F90

@ -45,6 +45,12 @@
"full step splitting: u + v + w ", &
"half step splitting: u/2 + v/2 + w + v/2 + u/2", &
"hor/ver splitting: uv + w "/)
#ifdef STATIC
REALTYPE,dimension(I3DFIELD),public :: ffluxu2_3d, ffluxv2_3d
#else
REALTYPE,dimension(:,:,:),allocatable,public :: ffluxu2_3d, ffluxv2_3d
#endif
!
! !REVISION HISTORY:
! Original author(s): Karsten Bolding & Hans Burchard
@ -86,6 +92,9 @@
! !USES:
IMPLICIT NONE
!
! !LOCAL VARIABLES:
integer :: rc
!
!EOP
!-------------------------------------------------------------------------
!BOC
@ -106,6 +115,14 @@
end if
#endif
#ifndef STATIC
allocate(ffluxu2_3d(I3DFIELD),stat=rc) ! work array
if (rc /= 0) stop 'init_advection_3d: Error allocating memory (ffluxu2_3d)'
allocate(ffluxv2_3d(I3DFIELD),stat=rc) ! work array
if (rc /= 0) stop 'init_advection_3d: Error allocating memory (ffluxv2_3d)'
#endif
#ifdef DEBUG
write(debug,*) 'Leaving init_advection_3d()'
write(debug,*)
@ -376,6 +393,8 @@
end if
end if
ffluxu2_3d = _ZERO_ ; ffluxv2_3d = _ZERO_
if (present(hires)) then
p_hi => hires
else
@ -405,12 +424,15 @@
if (calc_ffluxu) pa_ffluxu2d(k)%p2d = _ZERO_
if (calc_ffluxv) pa_ffluxv2d(k)%p2d = _ZERO_
#endif
ffluxu2 = _ZERO_ ; ffluxv2 = _ZERO_
call do_advection(dt,f(:,:,k),uu(:,:,k),vv(:,:,k), &
hu(:,:,k),hv(:,:,k),ho(:,:,k),hn(:,:,k), &
split,hscheme,AH,tag2d, &
Dires=p_hi(:,:,k),advres=p_adv(:,:,k), &
ffluxu=pa_ffluxu2d(k)%p2d,ffluxv=pa_ffluxv2d(k)%p2d, &
nvd=pa_nvd2d(k)%p2d)
ffluxu2_3d(:,:,k) = ffluxu2
ffluxv2_3d(:,:,k) = ffluxv2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxu) ffluxu(:,:,k) = pa_ffluxu2d(k)%p2d
@ -435,19 +457,23 @@
if (calc_ffluxu) pa_ffluxu2d(k)%p2d = _ZERO_
if (calc_ffluxv) pa_ffluxv2d(k)%p2d = _ZERO_
#endif
ffluxu2 = _ZERO_
call adv_split_u(dt,f(:,:,k),fi(:,:,k),p_hi(:,:,k), &
p_adv(:,:,k),uu(:,:,k),hu(:,:,k), &
adv_grid%dxu,adv_grid%dyu,adv_grid%arcd1, &
_ONE_,hscheme,AH, &
adv_grid%mask_uflux,adv_grid%mask_uupdate, &
pa_ffluxu2d(k)%p2d,pa_nvd2d(k)%p2d)
ffluxu2_3d(:,:,k) = ffluxu2
#ifndef SLICE_MODEL
ffluxv2 = _ZERO_
call adv_split_v(dt,f(:,:,k),fi(:,:,k),p_hi(:,:,k), &
p_adv(:,:,k),vv(:,:,k),hv(:,:,k), &
adv_grid%dxv,adv_grid%dyv,adv_grid%arcd1, &
_ONE_,hscheme,AH, &
adv_grid%mask_vflux,adv_grid%mask_vupdate, &
pa_ffluxv2d(k)%p2d,pa_nvd2d(k)%p2d)
ffluxv2_3d(:,:,k) = ffluxv2
#endif
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
@ -532,12 +558,15 @@
if (calc_ffluxu) pa_ffluxu2d(k)%p2d = _ZERO_
if (calc_ffluxv) pa_ffluxv2d(k)%p2d = _ZERO_
#endif
ffluxu2 = _ZERO_ ; ffluxv2 = _ZERO_
call do_advection(dt,f(:,:,k),uu(:,:,k),vv(:,:,k), &
hu(:,:,k),hv(:,:,k),ho(:,:,k),hn(:,:,k), &
FULLSPLIT,hscheme,AH,tag2d, &
Dires=p_hi(:,:,k),advres=p_adv(:,:,k), &
ffluxu=pa_ffluxu2d(k)%p2d,ffluxv=pa_ffluxv2d(k)%p2d, &
nvd=pa_nvd2d(k)%p2d)
ffluxu2_3d(:,:,k) = ffluxu2
ffluxv2_3d(:,:,k) = ffluxv2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxu) ffluxu(:,:,k) = pa_ffluxu2d(k)%p2d
@ -560,12 +589,14 @@
if (calc_nvd) pa_nvd2d (k)%p2d = _ZERO_
if (calc_ffluxu) pa_ffluxu2d(k)%p2d = _ZERO_
#endif
ffluxu2 = _ZERO_
call adv_split_u(dt,f(:,:,k),f(:,:,k),p_hi(:,:,k), &
p_adv(:,:,k),uu(:,:,k),hu(:,:,k), &
adv_grid%dxu,adv_grid%dyu,adv_grid%arcd1, &
_HALF_,hscheme,AH, &
adv_grid%mask_uflux,adv_grid%mask_uupdate, &
pa_ffluxu2d(k)%p2d,pa_nvd2d(k)%p2d)
ffluxu2_3d(:,:,k) = ffluxu2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxu) ffluxu(:,:,k) = pa_ffluxu2d(k)%p2d
@ -585,12 +616,14 @@
if (calc_nvd) pa_nvd2d (k)%p2d = nvd(:,:,k)
if (calc_ffluxv) pa_ffluxv2d(k)%p2d = _ZERO_
#endif
ffluxv2 = _ZERO_
call adv_split_v(dt,f(:,:,k),f(:,:,k),p_hi(:,:,k), &
p_adv(:,:,k),vv(:,:,k),hv(:,:,k), &
adv_grid%dxv,adv_grid%dyv,adv_grid%arcd1, &
_HALF_,hscheme,AH, &
adv_grid%mask_vflux,adv_grid%mask_vupdate, &
pa_ffluxv2d(k)%p2d,pa_nvd2d(k)%p2d)
ffluxv2_3d(:,:,k) = ffluxv2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxv) ffluxv(:,:,k) = pa_ffluxv2d(k)%p2d
@ -616,12 +649,14 @@
if (calc_nvd) pa_nvd2d (k)%p2d = nvd (:,:,k)
if (calc_ffluxv) pa_ffluxv2d(k)%p2d = ffluxv(:,:,k)
#endif
ffluxv2 = ffluxv2_3d(:,:,k)
call adv_split_v(dt,f(:,:,k),f(:,:,k),p_hi(:,:,k), &
p_adv(:,:,k),vv(:,:,k),hv(:,:,k), &
adv_grid%dxv,adv_grid%dyv,adv_grid%arcd1, &
_HALF_,hscheme,AH, &
adv_grid%mask_vflux,adv_grid%mask_vupdate, &
pa_ffluxv2d(k)%p2d,pa_nvd2d(k)%p2d)
ffluxv2_3d(:,:,k) = ffluxv2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxv) ffluxv(:,:,k) = pa_ffluxv2d(k)%p2d
@ -641,12 +676,14 @@
if (calc_nvd) pa_nvd2d (k)%p2d = nvd (:,:,k)
if (calc_ffluxu) pa_ffluxu2d(k)%p2d = ffluxu(:,:,k)
#endif
ffluxu2 = ffluxu2_3d(:,:,k)
call adv_split_u(dt,f(:,:,k),f(:,:,k),p_hi(:,:,k), &
p_adv(:,:,k),uu(:,:,k),hu(:,:,k), &
adv_grid%dxu,adv_grid%dyu,adv_grid%arcd1, &
_HALF_,hscheme,AH, &
adv_grid%mask_uflux,adv_grid%mask_uupdate, &
pa_ffluxu2d(k)%p2d,pa_nvd2d(k)%p2d)
ffluxu2_3d(:,:,k) = ffluxu2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxu) ffluxu(:,:,k) = pa_ffluxu2d(k)%p2d
@ -671,13 +708,15 @@
if (calc_ffluxu) pa_ffluxu2d(k)%p2d = _ZERO_
if (calc_ffluxv) pa_ffluxv2d(k)%p2d = _ZERO_
#endif
ffluxu2 = _ZERO_ ; ffluxv2 = _ZERO_
call do_advection(dt,f(:,:,k),uu(:,:,k),vv(:,:,k), &
hu(:,:,k),hv(:,:,k),ho(:,:,k),hn(:,:,k), &
NOSPLIT,hscheme,AH,tag2d, &
Dires=p_hi(:,:,k),advres=p_adv(:,:,k), &
ffluxu=pa_ffluxu2d(k)%p2d,ffluxv=pa_ffluxv2d(k)%p2d, &
nvd=pa_nvd2d(k)%p2d)
ffluxu2_3d(:,:,k) = ffluxu2
ffluxv2_3d(:,:,k) = ffluxv2
#ifndef _POINTER_REMAP_
if (calc_nvd) nvd (:,:,k) = pa_nvd2d (k)%p2d
if (calc_ffluxu) ffluxu(:,:,k) = pa_ffluxu2d(k)%p2d

2
src/3d/getm_fabm.F90

@ -70,7 +70,7 @@
type t_pa3d
REALTYPE,dimension(:,:,:),pointer,contiguous :: p3d
end type t_pa3d
type(t_pa3d),dimension(:),allocatable :: pa_nummix_fabm_pel,pa_phymix_fabm_pel
type(t_pa3d),dimension(:),allocatable,public :: pa_nummix_fabm_pel,pa_phymix_fabm_pel
integer :: fabm_adv_split=0
integer :: fabm_adv_hor=1
integer :: fabm_adv_ver=1

2
src/3d/m3d.F90

@ -329,7 +329,7 @@
else
T = _ZERO_ ; S = _ZERO_ ; rho = _ZERO_
if(calc_temp) call init_temperature(hotstart)
if(calc_salt) call init_salinity(hotstart)
if(calc_salt) call init_salinity(runtype,hotstart)
if (runtype .eq. 4) then
LEVEL2 'AH_min_method = ', AH_min_method
select case (AH_min_method)

76
src/3d/rivers.F90

@ -28,11 +28,13 @@
! !USES:
use domain, only: imin,jmin,imax,jmax,ioff,joff
use domain, only: H,az,kmax,arcd1
use domain, only: nriverl,ri,rj
use time, only: write_time_string,timestr
use variables_2d, only: dtm,z
#ifndef NO_BAROCLINIC
use m3d, only: update_salt,update_temp
use m3d, only: update_salt,update_temp,dt
use variables_3d, only: hn,ssen,T,S
use variables_3d, only: phymix_S
#endif
#ifdef GETM_BIO
use bio, only: bio_calc
@ -76,6 +78,7 @@
#ifdef _FABM_
REALTYPE, public, allocatable :: river_fabm(:,:)
#endif
integer, public, allocatable :: rnl(:)
!
! !PRIVATE DATA MEMBERS:
integer :: river_format=2
@ -222,9 +225,11 @@
else if (.not. outside) then
xxx = ' inside'
ok(n) = 1
nriverl = nriverl + 1
else
xxx = ' in halo'
ok(n) = 2
nriverl = nriverl + 1
end if
bathy = H(i,j)
if (rzu(n) .gt. rzl(n)) then
@ -249,6 +254,23 @@
LEVEL3 trim(line)
end do
allocate(ri(nriverl),stat=rc)
if (rc /= 0) stop 'rivers: Error allocating memory (ri)'
allocate(rj(nriverl),stat=rc)
if (rc /= 0) stop 'rivers: Error allocating memory (rj)'
allocate(rnl(nriverl),stat=rc)
if (rc /= 0) stop 'rivers: Error allocating memory (rnl)'
nn = 0
do n=1,nriver
if (ok(n).eq.1 .or. ok(n).eq.2) then
nn = nn + 1
ri(nn) = ir(n) - ioff
rj(nn) = jr(n) - joff
rnl(nn) = n
end if
end do
! Calculate the number of used gridboxes.
! This particular section is prepared for multi-cell rivers not-in-sequence.
LEVEL3 'Number of unique rivers: ',rriver
@ -511,6 +533,9 @@
! the layer heights are increased proportionally.
!
! !USES:
use variables_2d, only: rvol_2d=>rvol
use variables_3d, only: saltbins,Si_s,rvol_s
use getm_timers , only: tic,toc,TIM_TEF
IMPLICIT NONE
!
! !INPUT PARAMETERS:
@ -518,10 +543,12 @@
logical, intent(in) :: do_3d
!
! !LOCAL VARIABLES:
integer :: i,j,k,m,n
integer :: i,j,k,l,m,n,nn
integer :: kl,kh
REALTYPE :: rvol,height
REALTYPE :: river_depth,x
REALTYPE :: S_old(0:kmax),S_riv(0:kmax)
integer :: binidx(0:kmax)
!EOP
!-----------------------------------------------------------------------
!BOC
@ -548,16 +575,18 @@
select case (river_method)
case(0)
case(1,2)
nn = 0
do n=1,nriver
if(ok(n) .gt. 0) then
nn = nn + 1
i = ir(n)-ioff; j = jr(n)-joff
rvol = ramp * dtm * river_flow(n) * flow_fraction(n)
if (allocated(rvol_2d)) rvol_2d(i,j) = rvol
if (ok(n) .eq. 1) then
irr(n) = irr(n) + rvol
end if
height = rvol * ARCD1
z(i,j) = z(i,j) + height
#ifndef NO_BAROCLINIC
macro_height(n)=macro_height(n)+height
! on macrotime step adjust 3d fields
if (do_3d) then
@ -582,16 +611,40 @@
kh = m
end if
river_depth = sum(hn(i,j,kl:kh))
! Changes of total and layer height due to river inflow:
ssen(i,j) = ssen(i,j)+macro_height(n)
x = macro_height(n) / river_depth
hn(i,j,kl:kh) = hn(i,j,kl:kh) * ( _ONE_ + x )
#ifndef NO_BAROCLINIC
if (macro_height(n).gt._ZERO_ .or. .not.river_outflow_properties_follow_source_cell) then
if (update_salt ) then
S_old(:) = S(i,j,:)
if ( river_salt(n) .ne. salt_missing) then
S(i,j,kl:kh) = (S(i,j,kl:kh)*river_depth &
+ river_salt(n)*macro_height(n)) &
/ (river_depth+macro_height(n))
S(i,j,kl:kh) = ( S(i,j,kl:kh) + x*river_salt(n) ) / ( _ONE_ + x )
S_riv(kl:kh) = river_salt(n)
else
S(i,j,kl:kh) = S(i,j,kl:kh)*river_depth &
/ (river_depth+macro_height(n))
S(i,j,kl:kh) = S(i,j,kl:kh) / ( _ONE_ + x )
S_riv(kl:kh) = _ZERO_
end if
! only need to calculate mixing here (no mixing for zero-gradient)
if (associated(phymix_S)) then
phymix_S(i,j,kl:kh) = &
hn(i,j,kl:kh)*( ( S_old(kl:kh)**2 + x*S_riv(kl:kh)**2 ) / ( _ONE_ + x ) - S(i,j,kl:kh)**2 ) / dt
end if
if (allocated(rvol_s)) then
call tic(TIM_TEF)
rvol_s(i,j,:) = _ZERO_
call rebin(kmax,S_old,saltbins,Si_s,binidx)
do k=kl,kh
l = binidx(k)
if (l .lt. 0) cycle
rvol_s(i,j,l) = rvol_s(i,j,l) + x*hn(i,j,k)/((_ONE_+x)*ARCD1)
end do
call toc(TIM_TEF)
end if
end if
if (update_temp .and. river_temp(n) .ne. temp_missing) then
T(i,j,kl:kh) = (T(i,j,kl:kh)*river_depth &
@ -623,14 +676,11 @@
end if
#endif
end if
#endif
! Changes of total and layer height due to river inflow:
hn(i,j,kl:kh) = hn(i,j,kl:kh)/river_depth &
*(river_depth+macro_height(n))
ssen(i,j) = ssen(i,j)+macro_height(n)
macro_height(n) = _ZERO_
end if
#endif
end if
end do
case default

334
src/3d/salinity.F90

@ -20,6 +20,7 @@
!KB use get_field, only: get_3d_field
use variables_2d, only: fwf_int
use variables_3d, only: rk,S,hn,kmin,sss
use variables_3d, only: saltbins,saltbinmin,saltbinmax,Si_s
use halo_zones, only: update_3d_halo,wait_halo,D_TAG,H_TAG
IMPLICIT NONE
@ -61,6 +62,7 @@
integer :: salt_check=0
REALTYPE :: min_salt=_ZERO_,max_salt=40*_ONE_
integer,public :: nonnegsalt_method=0
integer, allocatable :: binidx_kij(:,:,:)
!
! !REVISION HISTORY:
! Original author(s): Karsten Bolding & Hans Burchard
@ -77,7 +79,7 @@
! \label{sec-init-salinity}
!
! !INTERFACE:
subroutine init_salinity(hotstart)
subroutine init_salinity(runtype,hotstart)
!
! !DESCRIPTION:
!
@ -104,10 +106,11 @@
IMPLICIT NONE
!
! !INPUT PARAMETERS:
integer,intent(in) :: runtype
logical,intent(in) :: hotstart
!
! !LOCAL VARIABLES:
integer :: i,j,k,n,rc
integer :: i,j,k,l,n,rc
integer :: status
NAMELIST /salt/ &
salt_method,salt_const,salt_file, &
@ -115,6 +118,7 @@
salt_adv_split,salt_adv_hor,salt_adv_ver, &
avmols,salt_AH_method,salt_AH_const,salt_AH_Prt, &
salt_AH_stirr_const,sss_nudging_time, &
saltbins,saltbinmin,saltbinmax, &
salt_check,min_salt,max_salt, &
nonnegsalt_method
!EOP
@ -136,10 +140,6 @@
LEVEL3 'avmols = ',real(avmols)
end if
if (.not. hotstart) then
call init_salinity_field()
end if
! Sanity checks for advection specifications
LEVEL3 'Advection of salinity'
if (salt_adv_hor .eq. J7) stop 'init_salinity: J7 not implemented yet'
@ -188,6 +188,22 @@
allocate(sss(E2DFIELD),stat=rc)
end if
if (runtype.eq.4 .and. saltbins.gt.0) then
LEVEL3 'salinity space for TEF:',real(saltbinmin),real(saltbinmax),'(',saltbins,'bins)'
allocate(Si_s(0:saltbins),stat=rc)
if (rc /= 0) stop 'postinit_3d: Error allocating memory (Si_s)'
do l=0,saltbins
Si_s(l) = saltbinmin + (saltbinmax-saltbinmin)/saltbins * l
end do
allocate(binidx_kij(_KRANGE_,I2DFIELD),stat=rc)
if (rc /= 0) stop 'postinit_3d: Error allocating memory (binidx_kij)'
binidx_kij = -1
end if
if (.not. hotstart) then
call init_salinity_field()
end if
LEVEL3 'salt_check=',salt_check
if (salt_check .ne. 0) then
LEVEL4 'doing sanity check on salinity'
@ -298,6 +314,14 @@
call wait_halo(D_TAG)
call mirror_bdy_3d(S,D_TAG)
if (allocated(binidx_kij)) then
do j = jmin-HALO,jmax+HALO
do i = imin-HALO,imax+HALO
if (az(i,j).ne.0) call rebin(kmax,S(i,j,:),saltbins,Si_s,binidx_kij(:,i,j))
end do
end do
end if
#ifdef DEBUG
write(debug,*) 'Leaving init_salinity_field()'
write(debug,*)
@ -312,7 +336,7 @@
! !IROUTINE: do_salinity - salinity equation \label{sec-do-salinity}
!
! !INTERFACE:
subroutine do_salinity(n)
subroutine do_salinity(loop)
!
! !DESCRIPTION:
!
@ -333,24 +357,40 @@
!
! !USES:
use advection_3d, only: do_advection_3d
use advection_3d, only: ffluxu2_3d,ffluxv2_3d
use variables_3d, only: dt,cnpar,hn,ho,nuh,uu,vv,ww,hun,hvn
use domain, only: imin,imax,jmin,jmax,kmax,az
use getm_timers, only: tic,toc,TIM_SALT,TIM_SALTH,TIM_MIXANALYSIS
use domain, only: imin,imax,jmin,jmax,kmax,az,au,av
use domain, only: dyu,dxv
use domain, only: nriverl,ri,rj
use getm_timers, only: tic,toc,TIM_SALT,TIM_SALTH,TIM_MIXANALYSIS,TIM_TEF
use variables_3d, only: Sfluxu,Sfluxv,Sfluxw
use variables_3d, only: nummix_S, phymix_S
use variables_3d, only: Sfluxu2,Sfluxv2
use variables_3d, only: counts_s,flags_s,h_s,hS_s,hS2_s,hpmS_s,hnmS_s
use variables_3d, only: hu_s,uu_s,Sfluxu_s,S2fluxu_s
use variables_3d, only: hv_s,vv_s,Sfluxv_s,S2fluxv_s
use variables_3d, only: fwf_s,fwfS2_s, phymix_S_fwf_s, phymix_S_riv_s
use variables_3d, only: nummix_S, phymix_S, phymix_S_fwf
!$ use omp_lib
IMPLICIT NONE
!
! !INPUT PARAMETERS:
integer, intent(in) :: n
integer, intent(in) :: loop
!
! !LOCAL VARIABLES:
integer :: i,j,k,rc
integer :: i,j,k,l,n,rc
REALTYPE, POINTER :: Res(:)
REALTYPE, POINTER :: auxn(:),auxo(:)
REALTYPE, POINTER :: a1(:),a2(:),a3(:),a4(:)
REALTYPE, pointer :: fluxw(:)
REALTYPE :: sss_old(I2DFIELD)
REALTYPE, dimension(I3DFIELD) :: Sold,Sfluxu_tmp,Sfluxv_tmp
REALTYPE, dimension(I2DFIELD) :: wrk2d
REALTYPE :: rivmix(0:kmax,1:nriverl)
REALTYPE :: deltaS,hu_k,uu_k
REALTYPE :: S_k(0:kmax)
integer :: binidx(0:kmax)
integer :: status
logical :: calc_phymix_fwf
!EOP
!-----------------------------------------------------------------------
!BOC
@ -361,6 +401,15 @@
#endif
call tic(TIM_SALT)
if (allocated(phymix_S_riv_s)) then
do n=1,nriverl
rivmix(:,n) = phymix_S(ri(n),rj(n),:)
end do
end if
calc_phymix_fwf = ( associated(phymix_S) .or. allocated(phymix_S_fwf) )
if (calc_phymix_fwf) sss_old = S(:,:,kmax)
do j=jmin-HALO,jmax+HALO
do i=imin-HALO,imax+HALO
if (az(i,j) .eq. 1) then
@ -369,28 +418,202 @@
end do
end do
if (calc_phymix_fwf) then
! Note (KK): division by hn because we already refer to new volume
! where diffusion will be added
where (az .eq. 1) wrk2d = sss_old * S(:,:,kmax) * fwf_int / dt
if (associated(phymix_S)) then
where (az .eq. 1) phymix_S(:,:,kmax) = phymix_S(:,:,kmax) + wrk2d
end if
if (allocated(phymix_S_fwf)) then
where (az .eq. 1) phymix_S_fwf = wrk2d
end if
end if
if (allocated(fwf_s).or.allocated(fwfS2_s)) then
call tic(TIM_TEF)
if (allocated(fwf_s )) fwf_s = _ZERO_
if (allocated(fwfS2_s)) fwfS2_s = _ZERO_
do j = jmin,jmax
do i = imin,imax
if (az(i,j) .eq. 1) then
l = binidx_kij(kmax,i,j) ! binning based on old salinity
if (l .ge. 0) then
if (allocated(fwf_s )) fwf_s (i,j,l) = fwf_int(i,j)
if (allocated(fwfS2_s)) fwfS2_s(i,j,l) = fwf_int(i,j)*S(i,j,kmax)*S(i,j,kmax)
end if
end if
end do
end do
call toc(TIM_TEF)
end if
call do_advection_3d(dt,S,uu,vv,ww,hun,hvn,ho,hn, &
salt_adv_split,salt_adv_hor,salt_adv_ver,_ZERO_,H_TAG, &
ffluxu=Sfluxu,ffluxv=Sfluxv,ffluxw=Sfluxw,nvd=nummix_S)
if (allocated(Sfluxu2)) Sfluxu2 = ffluxu2_3d
if (allocated(Sfluxv2)) Sfluxv2 = ffluxv2_3d
! KK-TODO: do we have to treat salt fluxes from different
! directional-split steps separately?
if (allocated(hu_s).or.allocated(uu_s).or.allocated(Sfluxu_s).or.allocated(S2fluxu_s)) then
call tic(TIM_TEF)
if (allocated(hu_s )) hu_s = _ZERO_
if (allocated(uu_s )) uu_s = _ZERO_
if (allocated(Sfluxu_s )) Sfluxu_s = _ZERO_
if (allocated(S2fluxu_s)) S2fluxu_s = _ZERO_
do j = jmin,jmax
do i = imin,imax
if (au(i,j).eq.1 .or. au(i,j).eq.2) then
where (uu(i,j,1:kmax) .ne. _ZERO_)
S_k(1:kmax) = Sfluxu(i,j,1:kmax) / uu(i,j,1:kmax) / DYU
else where
S_k(1:kmax) = saltbinmax + 9999.0_rk
end where
call rebin(kmax,S_k,saltbins,Si_s,binidx)
do k=1,kmax
l = binidx(k)
if (l .lt. 0) cycle
if (allocated(hu_s )) hu_s (i,j,l) = hu_s (i,j,l) + hun(i,j,k)
if (allocated(uu_s )) uu_s (i,j,l) = uu_s (i,j,l) + uu(i,j,k)
if (allocated(Sfluxu_s )) Sfluxu_s (i,j,l) = Sfluxu_s (i,j,l) + Sfluxu(i,j,k)
if (allocated(S2fluxu_s)) S2fluxu_s(i,j,l) = S2fluxu_s(i,j,l) + Sfluxu2(i,j,k)
end do
end if
end do
end do
call toc(TIM_TEF)
end if
if (allocated(hv_s).or.allocated(vv_s).or.allocated(Sfluxv_s).or.allocated(S2fluxv_s)) then
call tic(TIM_TEF)
if (allocated(hv_s )) hv_s = _ZERO_
if (allocated(vv_s )) vv_s = _ZERO_
if (allocated(Sfluxv_s )) Sfluxv_s = _ZERO_
if (allocated(S2fluxv_s)) S2fluxv_s = _ZERO_
do j = jmin,jmax
do i = imin,imax
if (av(i,j).eq.1 .or. av(i,j).eq.2) then
where (vv(i,j,1:kmax) .ne. _ZERO_)
S_k(1:kmax) = Sfluxv(i,j,1:kmax) / vv(i,j,1:kmax) / DXV
else where
S_k(1:kmax) = saltbinmax + 9999.0_rk
end where
call rebin(kmax,S_k,saltbins,Si_s,binidx)
do k=1,kmax
l = binidx(k)
if (l .lt. 0) cycle
if (allocated(hv_s )) hv_s (i,j,l) = hv_s (i,j,l) + hvn(i,j,k)
if (allocated(vv_s )) vv_s (i,j,l) = vv_s (i,j,l) + vv(i,j,k)
if (allocated(Sfluxv_s )) Sfluxv_s (i,j,l) = Sfluxv_s (i,j,l) + Sfluxv(i,j,k)
if (allocated(S2fluxv_s)) S2fluxv_s(i,j,l) = S2fluxv_s(i,j,l) + Sfluxv2(i,j,k)
end do
end if
end do
end do
call toc(TIM_TEF)
end if
if (salt_AH_method .gt. 0) then
! S is not halo updated after advection
call tic(TIM_SALTH)
call update_3d_halo(S,S,az,imin,jmin,imax,jmax,kmax,D_TAG)
call wait_halo(D_TAG)
call toc(TIM_SALTH)
if ( allocated(hu_s).or.allocated(uu_s).or.allocated(Sfluxu_s).or.allocated(S2fluxu_s) &
.or.allocated(hv_s).or.allocated(vv_s).or.allocated(Sfluxv_s).or.allocated(S2fluxv_s)) then
call tic(TIM_TEF)
Sold = S
do j = jmin,jmax+1
do i = imin,imax+1
if (az(i,j).ne.0) call rebin(kmax,S(i,j,:),saltbins,Si_s,binidx_kij(:,i,j))
end do
end do
call toc(TIM_TEF)
end if
if (allocated(hu_s).or.allocated(uu_s).or.allocated(Sfluxu_s)) Sfluxu_tmp = Sfluxu
if (allocated(hv_s).or.allocated(vv_s).or.allocated(Sfluxv_s)) Sfluxv_tmp = Sfluxv
call tracer_diffusion(S,hn,salt_AH_method,salt_AH_const,salt_AH_Prt,salt_AH_stirr_const, &
ffluxu=Sfluxu,ffluxv=Sfluxv, &
phymix=phymix_S)
end if
if (allocated(Sfluxu2)) Sfluxu2 = Sfluxu2 + ffluxu2_3d
if (allocated(Sfluxv2)) Sfluxv2 = Sfluxv2 + ffluxv2_3d
if (allocated(hu_s).or.allocated(uu_s).or.allocated(Sfluxu_s).or.allocated(S2fluxu_s)) then
call tic(TIM_TEF)
Sfluxu_tmp = Sfluxu - Sfluxu_tmp
do j = jmin,jmax
do i = imin,imax
if (au(i,j).eq.1 .or. au(i,j).eq.2) then
do k=1,kmax
if (Sfluxu_tmp(i,j,k) .ne. _ZERO_) then
deltaS = Sold(i+1,j,k) - Sold(i,j,k)
!if (deltaS .ne. _ZERO_) then
hu_k = _HALF_ * ( hn(i,j,k) + hn(i+1,j,k) )
uu_k = Sfluxu_tmp(i,j,k) / deltaS / DYU
! fluxes with S(i)
l = binidx_kij(k,i,j)
if (l .ge. 0) then
if (allocated(hu_s )) hu_s (i,j,l) = hu_s (i,j,l) + hu_k
if (allocated(uu_s )) uu_s (i,j,l) = uu_s (i,j,l) - uu_k
if (allocated(Sfluxu_s )) Sfluxu_s (i,j,l) = Sfluxu_s (i,j,l) - uu_k*DYU*Sold(i,j,k)
if (allocated(S2fluxu_s)) S2fluxu_s(i,j,l) = S2fluxu_s(i,j,l) - uu_k*DYU*Sold(i,j,k)*Sold(i,j,k)
end if
! fluxes with S(i+1)
l = binidx_kij(k,i+1,j)
if (l .ge. 0) then
if (allocated(hu_s )) hu_s (i,j,l) = hu_s (i,j,l) + hu_k
if (allocated(uu_s )) uu_s (i,j,l) = uu_s (i,j,l) + uu_k
if (allocated(Sfluxu_s )) Sfluxu_s (i,j,l) = Sfluxu_s (i,j,l) + uu_k*DYU*Sold(i+1,j,k)
if (allocated(S2fluxu_s)) S2fluxu_s(i,j,l) = S2fluxu_s(i,j,l) + uu_k*DYU*Sold(i+1,j,k)*Sold(i+1,j,k)
end if
end if
end do
end if
end do
end do
call toc(TIM_TEF)
end if
if (allocated(hv_s).or.allocated(vv_s).or.allocated(Sfluxv_s).or.allocated(S2fluxv_s)) then
call tic(TIM_TEF)
Sfluxv_tmp = Sfluxv - Sfluxv_tmp
do j = jmin,jmax
do i = imin,imax
if (av(i,j).eq.1 .or. av(i,j).eq.2) then
do k=1,kmax
if (Sfluxv_tmp(i,j,k) .ne. _ZERO_) then
deltaS = Sold(i,j+1,k) - Sold(i,j,k)
!if (deltaS .ne. _ZERO_) then
hu_k = _HALF_ * ( hn(i,j,k) + hn(i,j+1,k) )
uu_k = Sfluxv_tmp(i,j,k) / deltaS / DXV
! fluxes with S(j)
l = binidx_kij(k,i,j)
if (l .ge. 0) then
if (allocated(hv_s )) hv_s (i,j,l) = hv_s (i,j,l) + hu_k
if (allocated(vv_s )) vv_s (i,j,l) = vv_s (i,j,l) - uu_k
if (allocated(Sfluxv_s )) Sfluxv_s (i,j,l) = Sfluxv_s (i,j,l) - uu_k*DXV*Sold(i,j,k)
if (allocated(S2fluxv_s)) S2fluxv_s(i,j,l) = S2fluxv_s(i,j,l) - uu_k*DXV*Sold(i,j,k)*Sold(i,j,k)
end if
! fluxes with S(j+1)
l = binidx_kij(k,i,j+1)
if (l .ge. 0) then
if (allocated(hv_s )) hv_s (i,j,l) = hv_s (i,j,l) + hu_k
if (allocated(vv_s )) vv_s (i,j,l) = vv_s (i,j,l) + uu_k
if (allocated(Sfluxv_s )) Sfluxv_s (i,j,l) = Sfluxv_s (i,j,l) + uu_k*DXV*Sold(i,j+1,k)
if (allocated(S2fluxv_s)) S2fluxv_s(i,j,l) = S2fluxv_s(i,j,l) + uu_k*DXV*Sold(i,j+1,k)*Sold(i,j+1,k)
end if
end if
end do
end if
end do
end do
call toc(TIM_TEF)
end if
if (associated(phymix_S)) then
call toc(TIM_SALT)
call tic(TIM_MIXANALYSIS)
call physical_mixing(S,avmols,phymix_S,salt_AH_method)
call toc(TIM_MIXANALYSIS)
call tic(TIM_SALT)
end if
@ -417,7 +640,7 @@
if (rc /= 0) stop 'do_salinity: Error allocating memory (auxo)'
fluxw => null()
if (associated(Sfluxw)) then
if (associated(Sfluxw) .or. associated(phymix_S)) then
allocate(fluxw(0:kmax),stat=rc) ! work array
if (rc /= 0) stop 'do_salinity: Error allocating memory (fluxw)'
end if
@ -486,6 +709,28 @@
if (associated(Sfluxw)) then
Sfluxw(i,j,1:kmax-1) = Sfluxw(i,j,1:kmax-1) + fluxw(1:kmax-1)
end if
if (associated(phymix_S)) then
! KK-TODO: only consider changes of S2 due to diffusion!
! (no sss_nudging or radiation for temperature)
do k=1,kmax-1
#if 0
fluxw(k) = fluxw(k) * &
( cnpar *(Res (k)+Res (k+1)) &
+ (_ONE_-cnpar)*(S(i,j,k)+S(i,j,k+1)))
#else
fluxw(k) = - ( &
auxo(k) * ( S(i,j,k+1) - S(i,j,k) ) * (S(i,j,k)+S(i,j,k+1)) &
+ auxn(k) * ( Res (k+1) - Res (k) ) * (Res (k)+Res (k+1)) &
) / dt
#endif
end do
do k=1,kmax
phymix_S(i,j,k) = phymix_S(i,j,k) &
- ( hn(i,j,k)*(Res(k)*Res(k) - S(i,j,k)*S(i,j,k))/dt &
+ (fluxw(k) - fluxw(k-1) ) )
end do
end if
do k=1,kmax
S(i,j,k)=Res(k)
@ -520,8 +765,55 @@
!$OMP END PARALLEL
if (allocated(counts_s).or.allocated(flags_s).or.allocated(h_s).or.allocated(hS_s).or.allocated(hS2_s).or.allocated(hpmS_s).or.allocated(hnmS_s).or.allocated(phymix_S_fwf_s).or.allocated(phymix_S_riv_s)) then
call tic(TIM_TEF)
if (allocated(counts_s)) counts_s = 0
if (allocated(flags_s )) flags_s = 0
if (allocated(h_s )) h_s = _ZERO_
if (allocated(hS_s )) hS_s = _ZERO_
if (allocated(hS2_s )) hS2_s = _ZERO_
if (allocated(hpmS_s )) hpmS_s = _ZERO_
if (allocated(hnmS_s )) hnmS_s = _ZERO_
if (allocated(phymix_S_fwf_s)) phymix_S_fwf_s = _ZERO_
if (allocated(phymix_S_riv_s)) phymix_S_riv_s = _ZERO_
do j = jmin,jmax
do i = imin,imax
if (az(i,j) .eq. 1) then
call rebin(kmax,S(i,j,:),saltbins,Si_s,binidx)
binidx_kij(:,i,j) = binidx(:)
do k=1,kmax
l = binidx(k)
if (l .lt. 0) cycle
if (allocated(counts_s)) counts_s(i,j,l) = counts_s(i,j,l) + 1
if (allocated(flags_s )) flags_s (i,j,l) = 1
if (allocated(h_s )) h_s (i,j,l) = h_s (i,j,l) + hn(i,j,k)
if (allocated(hS_s )) hS_s (i,j,l) = hS_s (i,j,l) + hn(i,j,k)*S(i,j,k)
if (allocated(hS2_s )) hS2_s (i,j,l) = hS2_s (i,j,l) + hn(i,j,k)*S(i,j,k)*S(i,j,k)
if (allocated(hpmS_s)) hpmS_s(i,j,l) = hpmS_s(i,j,l) + phymix_S(i,j,k)
if (allocated(hnmS_s)) hnmS_s(i,j,l) = hnmS_s(i,j,l) + nummix_S(i,j,k)
end do
l = binidx(kmax)
if (l .ge. 0) then
if (allocated(phymix_S_fwf_s)) phymix_S_fwf_s(i,j,l) = phymix_S_fwf(i,j)
end if
end if
end do
end do
if (allocated(phymix_S_riv_s)) then
do n=1,nriverl
j = rj(n)
i = ri(n)
do k=1,kmax
l = binidx_kij(k,i,j)
if (l .lt. 0) cycle
phymix_S_riv_s(i,j,l) = phymix_S_riv_s(i,j,l) + rivmix(k,n)
end do
end do
end if
call toc(TIM_TEF)
end if
if (salt_check .ne. 0 .and. mod(n,abs(salt_check)) .eq. 0) then
if (salt_check .ne. 0 .and. mod(loop,abs(salt_check)) .eq. 0) then
call check_3d_fields(imin,jmin,imax,jmax,kmin,kmax,az, &
S,min_salt,max_salt,status)
if (status .gt. 0) then

5
src/3d/tracer_diffusion.F90

@ -23,6 +23,7 @@
use variables_les, only: AmU_3d,AmV_3d
use m3d, only: AH_min_method
use variables_3d, only: AH_min
use advection_3d, only: ffluxu2_3d,ffluxv2_3d
use getm_timers, only: tic,toc,TIM_TRACEDIFF
!$ use omp_lib
IMPLICIT NONE
@ -130,6 +131,8 @@
calc_phymix = .false.
end if
ffluxu2_3d = _ZERO_ ; ffluxv2_3d = _ZERO_
! Note (KK): diffusion only within new layer heigts
!$OMP PARALLEL DEFAULT(SHARED) &
@ -420,6 +423,8 @@
end do
#ifndef SLICE_MODEL
end do
ffluxu2_3d(:,:,k) = ffluxu2
ffluxv2_3d(:,:,k) = ffluxv2
#endif
end if

257
src/3d/variables_3d.F90

@ -111,13 +111,13 @@
! {\tt init\_variables\_3d}) and cleanup (see {\tt clean\_variables\_3d}).
!
! !USES:
use parameters, only: rk
use domain, only: imin,imax,jmin,jmax,kmax,az,bottfric_method,rdrag
use waves , only: waveforcing_method,waves_method,NO_WAVES,WAVES_VF
use waves , only: waves_bbl_method,NO_WBBL
IMPLICIT NONE
!
! !PUBLIC DATA MEMBERS:
integer, parameter :: rk = kind(_ONE_)
REALTYPE :: dt,cnpar=0.9
REALTYPE :: avmback=_ZERO_,avhback=_ZERO_
logical :: deformC_3d=.false.
@ -129,6 +129,8 @@
logical :: save_Sfluxu=.false.
logical :: save_Sfluxv=.false.
logical :: save_Sfluxw=.false.
logical :: save_Sfluxu2=.false.
logical :: save_Sfluxv2=.false.
logical :: save_Tfluxu=.false.
logical :: save_Tfluxv=.false.
logical :: save_Tfluxw=.false.
@ -136,8 +138,12 @@
logical :: save_phydis_3d=.false.
logical :: save_nummix_S=.false.
logical :: save_phymix_S=.false.
logical :: save_phymix_S_fwf=.false.
logical :: save_nummix_T=.false.
logical :: save_phymix_T=.false.
integer :: saltbins = 0
REALTYPE :: saltbinmin=_ZERO_,saltbinmax=40.0_rk
integer :: id_dim_salt
!
#ifdef STATIC
#include "static_3d.h"
@ -160,6 +166,15 @@
REALTYPE,dimension(:,:,:),pointer,contiguous :: Tfluxu=>null()
REALTYPE,dimension(:,:,:),pointer,contiguous :: Tfluxv=>null()
REALTYPE,dimension(:,:,:),pointer,contiguous :: Tfluxw=>null()
REALTYPE,dimension(:,:,:),allocatable :: Sfluxu2
REALTYPE,dimension(:,:,:),allocatable :: Sfluxv2
REALTYPE,dimension(:) ,allocatable :: Si_s
REALTYPE,dimension(:,:,:),allocatable :: counts_s,flags_s,h_s,hS_s,hS2_s,hpmS_s,hnmS_s
REALTYPE,dimension(:,:,:),allocatable :: hu_s,uu_s,Sfluxu_s,S2fluxu_s
REALTYPE,dimension(:,:,:),allocatable :: hv_s,vv_s,Sfluxv_s,S2fluxv_s
REALTYPE,dimension(: