vertical_coordinates.F90 7.99 KB
Newer Older
gotm's avatar
gotm committed
1
#include "cppdefs.h"
2 3 4 5 6 7 8 9 10 11 12 13
!-----------------------------------------------------------------------
!BOP
!
! !MODULE: vertical_coordinates
!
! !INTERFACE:
   module vertical_coordinates
!
! !DESCRIPTION:
!
! !USES:
#ifdef SLICE_MODEL
14
   use variables_3d, only: kvmin
15
#endif
16 17
   use domain, only: imin,imax,jmin,jmax,kmax
   use domain, only: au,av
Knut's avatar
Knut committed
18
   use domain, only: H,vert_cord,maxdepth
Knut's avatar
Knut committed
19
   use halo_zones, only: U_TAG,V_TAG
20
   use variables_3d, only: ho,hn,huo,hun,hvo,hvn,hvel
Knut's avatar
Knut committed
21
   use variables_3d, only: zwn,zcn
22
   use variables_3d, only: Dun,Dvn
23 24 25 26 27
   use exceptions
   IMPLICIT NONE
!
! !PUBLIC DATA MEMBERS:
   public coordinates
Knut's avatar
Knut committed
28 29
   logical,public  :: restart_with_ho=.false.
   logical,public  :: restart_with_hn=.false.
30 31 32 33 34 35 36 37 38 39 40 41
   REALTYPE,public :: cord_relax=_ZERO_
!
! !PRIVATE DATA MEMBERS:
!
! !REVISION HISTORY:
!  Original author(s): Richard Hofmeister & Knut Klingbeil
!
!EOP
!-----------------------------------------------------------------------

   contains

gotm's avatar
gotm committed
42 43 44
!-----------------------------------------------------------------------
!BOP
!
hb's avatar
hb committed
45 46
! !ROUTINE:  coordinates - defines the vertical coordinate
! \label{sec-coordinates}
gotm's avatar
gotm committed
47 48
!
! !INTERFACE:
49
   subroutine coordinates(hotstart)
gotm's avatar
gotm committed
50 51 52
!
! !DESCRIPTION:
!
hb's avatar
hb committed
53 54 55 56
! Here, the vertical layer distribution in T-, U- and V-points is updated
! during every macro time step. This is done for the old and the new
! layer thicknesses at every point. Calculation of the layer distribution
! in the U- and V-points is done indepently from the calculation in the
57
! T-points, since different methods for the calculation of the
hb's avatar
hb committed
58 59 60
! bathymetry values in the U- and V-points are possible, see routine
! {\tt uv\_depths} described on page \pageref{sec-uv-depth}.
!
61 62 63 64 65 66 67 68 69
! The different methods for the vertical layer distribution
! are initialised and called to be chosen by the namelist paramter {\tt vert\_cord}:\\
! \\
! {\tt vert\_cord=1}: sigma coordinates (section~\ref{sec-sigma-coordinates}) \\
! {\tt vert\_cord=2}: z-level (not coded yet) \\
! {\tt vert\_cord=3}: general vertical coordinates (gvc, section~\ref{sec-general-coordinates})
! \\
! {\tt vert\_cord=5}: adaptive vertical coordinates (section~\ref{sec-adaptive-coordinates}) \\
! \\
hb's avatar
hb committed
70 71
!
!
gotm's avatar
gotm committed
72
! !USES:
73
   use variables_3d, only: zc
bjb's avatar
bjb committed
74
   use getm_timers, only: tic, toc,TIM_COORDS
gotm's avatar
gotm committed
75 76 77
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
78 79 80 81
!   integer, intent(in)                 :: cord_type
!   REALTYPE, intent(in)                :: cord_relax
!   REALTYPE, intent(in)                :: maxdepth
   logical, intent(in)                 :: hotstart
gotm's avatar
gotm committed
82
!
kbk's avatar
kbk committed
83 84 85
! !REVISION HISTORY:
!  Original author(s): Hans Burchard & Karsten Bolding
!
gotm's avatar
gotm committed
86
! !LOCAL VARIABLES:
87
   logical, save   :: first=.true.
88 89
   integer         :: ii
!   integer         :: preadapt=0
kbk's avatar
kbk committed
90
   integer          :: i,j,k
gotm's avatar
gotm committed
91 92 93 94 95 96 97 98
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'coordinates() # ',Ncall
#endif
bjb's avatar
bjb committed
99
   call tic(TIM_COORDS)
gotm's avatar
gotm committed
100 101

   if (first) then
102 103 104 105 106 107 108 109 110 111 112 113
      if (hotstart) then
         if ( .not.restart_with_ho .or. .not.restart_with_hn ) then
            STDERR LINE
            LEVEL3 "ho and hn missing in restart file!!!"
            LEVEL3 "This might be ok for some specific settings, but in"
            LEVEL3 "general you should do a zero-length simulation with"
            LEVEL3 "your previous coordinate settings to create a valid"
            LEVEL3 "restart file."
            STDERR LINE
         end if
      end if

114 115
      select case (vert_cord)
         case (_SIGMA_COORDS_) ! sigma coordinates
Knut's avatar
Knut committed
116
            LEVEL2 'using ',kmax,' sigma layers'
117
            call sigma_coordinates(.true.,hotstart)
118
         case (_Z_COORDS_) ! z-level
119
            call getm_error("coordinates()","z-levels not implemented yet")
120
         case (_GENERAL_COORDS_) ! general vertical coordinates
Knut's avatar
Knut committed
121
            LEVEL2 'using ',kmax,' gvc layers'
122
            call general_coordinates(.true.,hotstart,cord_relax,maxdepth)
123
         case (_HYBRID_COORDS_) ! hybrid vertical coordinates
Knut's avatar
Knut committed
124
            LEVEL2 'using ',kmax,' hybrid layers'
125 126 127
            call hybrid_coordinates(.true.)
STDERR 'coordinates(): hybrid_coordinates not coded yet'
stop
128
         case (_ADAPTIVE_COORDS_) ! adaptive vertical coordinates
Knut's avatar
Knut committed
129
            LEVEL2 'using ',kmax,' adaptive layers'
130
            call adaptive_coordinates(.true.,hotstart)
131 132
         case default
      end select
Knut's avatar
Knut committed
133 134 135
      if (.not. hotstart) then
         ho = hn
      end if
gotm's avatar
gotm committed
136
      first = .false.
137
   else
Knut's avatar
Knut committed
138
      ho  = hn  ! ho before advection (already including rivers and fwf)
139 140
      huo = hun
      hvo = hvn
141 142
      select case (vert_cord)
         case (_SIGMA_COORDS_) ! sigma coordinates
143
            call sigma_coordinates(.false.,hotstart)
144 145
         case (_Z_COORDS_) ! z-level
         case (_GENERAL_COORDS_) ! general vertical coordinates
146
            call general_coordinates(.false.,hotstart,cord_relax,maxdepth)
147
         case (_HYBRID_COORDS_) ! hybrid vertical coordinates
148
            call hybrid_coordinates(.false.)
149 150
         case (_ADAPTIVE_COORDS_) ! adaptive vertical coordinates
            call adaptive_coordinates(.false.,hotstart)
gotm's avatar
gotm committed
151 152 153 154
         case default
      end select
   end if ! first

Knut's avatar
Knut committed
155 156 157 158 159 160
   zwn(:,:,0) = -H
   do k=1,kmax
      zwn(:,:,k) = zwn(:,:,k-1) + hn(:,:,k)
   end do
   zcn(:,:,1:kmax) = _HALF_ * ( zwn(:,:,0:kmax-1) + zwn(:,:,1:kmax) )

Knut's avatar
Knut committed
161 162
   hvel = _HALF_ * ( ho + hn )

163 164 165 166 167 168 169 170 171 172 173 174
   if (first .and. hotstart) then
   hun(_IRANGE_HALO_-1,:,1:kmax) = &
      _HALF_ * ( hvel(_IRANGE_HALO_-1,:,1:kmax) + hvel(1+_IRANGE_HALO_,:,1:kmax) )
   hvn(:,_JRANGE_HALO_-1,1:kmax) = &
      _HALF_ * ( hvel(:,_JRANGE_HALO_-1,1:kmax) + hvel(:,1+_JRANGE_HALO_,1:kmax) )

!  KK-TODO: as long as hvel is based on "new" ho (including rivers)
!           hun and hvn do not coincide with Dun and Dvn
   call hcheck(hun,Dun,au)
   call hcheck(hvn,Dvn,av)
   end if

Knut's avatar
Knut committed
175 176 177 178 179 180 181
#ifdef _MIRROR_BDY_EXTRA_
!  Note (KK): required for calculation of SS
!             with non-zero velocity behind open bdy
   call mirror_bdy_3d(hun,U_TAG)
   call mirror_bdy_3d(hvn,V_TAG)
#endif

182
!  KK-TODO: remove because we already have zcn
Knut's avatar
Knut committed
183 184 185
   ! calculate the z-coordinate of the cell centers
   ! references to mean sea level
   zc(:,:,0)=-H(:,:)
Knut's avatar
Knut committed
186
   zc(:,:,1)=-H(:,:) + _HALF_*hn(:,:,1)
Knut's avatar
Knut committed
187
   do k=2,kmax
Knut's avatar
Knut committed
188
      zc(:,:,k)=zc(:,:,k-1)+_HALF_*(hn(:,:,k-1)+hn(:,:,k))
Knut's avatar
Knut committed
189 190
   end do

191
#ifdef SLICE_MODEL
Knut's avatar
Knut committed
192
   j = jmax/2
193
   do i=imin,imax
Knut's avatar
Knut committed
194 195 196
      do k=kvmin(i,j),kmax
         hvn(i,j-1,k)=hvn(i,j,k)
         hvn(i,j+1,k)=hvn(i,j,k)
197 198 199 200
      end do
   end do
#endif

bjb's avatar
bjb committed
201
   call toc(TIM_COORDS)
gotm's avatar
gotm committed
202
#ifdef DEBUG
203
   write(debug,*) 'Leaving coordinates()'
gotm's avatar
gotm committed
204 205 206 207 208
   write(debug,*)
#endif
   return
   end subroutine coordinates
!EOC
209 210 211 212 213 214
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE:  hcheck -
!
! !INTERFACE:
Knut's avatar
Knut committed
215
   subroutine hcheck(hn,Dn,mask)
216 217 218 219 220 221 222 223
!
! !DESCRIPTION:
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
   REALTYPE,intent(in)        :: Dn(I2DFIELD)
Knut's avatar
Knut committed
224
   integer,intent(in)         :: mask(E2DFIELD)
225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245
!
! !INPUT/OUTPUT PARAMETERS:
   REALTYPE,intent(inout)     :: hn(I3DFIELD)
!
! !REVISION HISTORY:
!  Original author(s): Richard Hofmeister
!
! !LOCAL VARIABLES:
   REALTYPE        :: HH,depthmin
   integer         :: i,j,k
!
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'hcheck() # ',Ncall
#endif

! Final check of layer thicnkess thoug not necessary if zpos treated correctly
Knut's avatar
Knut committed
246 247
   do j=jmin-HALO,jmax+HALO
      do i=imin-HALO,imax+HALO
Knut's avatar
Knut committed
248
         if (mask(i,j) .ne. 0) then
Knut's avatar
Knut committed
249 250 251 252 253 254 255
            HH=0.
            do k=1,kmax
               HH=HH+hn(i,j,k)
            end do
            do k=1,kmax
               hn(i,j,k)=hn(i,j,k)* Dn(i,j)/HH
            end do
Knut's avatar
Knut committed
256
         end if
257 258 259 260 261 262 263 264 265 266
      end do
   end do

#ifdef DEBUG
   write(debug,*) 'Leaving hcheck()'
   write(debug,*)
#endif
   return
   end subroutine hcheck
!EOC
267 268 269
!-----------------------------------------------------------------------

   end module vertical_coordinates
gotm's avatar
gotm committed
270 271

!-----------------------------------------------------------------------
272
! Copyright (C) 2012 - Hans Burchard and Karsten Bolding               !
gotm's avatar
gotm committed
273
!-----------------------------------------------------------------------