coordinates.F90 4.39 KB
Newer Older
gotm's avatar
gotm committed
1 2 3 4
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
hb's avatar
hb committed
5 6
! !ROUTINE:  coordinates - defines the vertical coordinate
! \label{sec-coordinates}
gotm's avatar
gotm committed
7 8
!
! !INTERFACE:
9
   subroutine coordinates(hotstart)
gotm's avatar
gotm committed
10 11 12
!
! !DESCRIPTION:
!
hb's avatar
hb committed
13 14 15 16
! 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
17
! T-points, since different methods for the calculation of the
hb's avatar
hb committed
18 19 20
! bathymetry values in the U- and V-points are possible, see routine
! {\tt uv\_depths} described on page \pageref{sec-uv-depth}.
!
21 22 23 24 25 26 27 28 29
! 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
30 31
!
!
gotm's avatar
gotm committed
32
! !USES:
Knut's avatar
Knut committed
33
   use domain, only: imin,imax,jmin,jmax,kmax,H
Knut's avatar
Knut committed
34
   use variables_3d, only: ho,hn,hvel
kbk's avatar
kbk committed
35 36 37
#ifdef SLICE_MODEL
   use variables_3d, only: kvmin,hvo,hvn
#endif
bjb's avatar
bjb committed
38
   use getm_timers, only: tic, toc,TIM_COORDS
39 40
   use m3d
   use domain, only: vert_cord
gotm's avatar
gotm committed
41 42 43
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
44 45 46 47
!   integer, intent(in)                 :: cord_type
!   REALTYPE, intent(in)                :: cord_relax
!   REALTYPE, intent(in)                :: maxdepth
   logical, intent(in)                 :: hotstart
gotm's avatar
gotm committed
48
!
kbk's avatar
kbk committed
49 50 51
! !REVISION HISTORY:
!  Original author(s): Hans Burchard & Karsten Bolding
!
gotm's avatar
gotm committed
52
! !LOCAL VARIABLES:
53
   logical, save   :: first=.true.
54 55
   integer         :: ii
!   integer         :: preadapt=0
kbk's avatar
kbk committed
56
   integer          :: i,j,k
gotm's avatar
gotm committed
57 58 59 60 61 62 63 64
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'coordinates() # ',Ncall
#endif
bjb's avatar
bjb committed
65
   call tic(TIM_COORDS)
gotm's avatar
gotm committed
66 67

   if (first) then
68 69
      select case (vert_cord)
         case (_SIGMA_COORDS_) ! sigma coordinates
70 71
            LEVEL2 'using sigma vertical coordinates'
            call sigma_coordinates(.true.)
72 73
         case (_Z_COORDS_) ! z-level
         case (_GENERAL_COORDS_) ! general vertical coordinates
74 75
            LEVEL2 'using general vertical coordinates'
            call general_coordinates(.true.,cord_relax,maxdepth)
76
         case (_HYBRID_COORDS_) ! hybrid vertical coordinates
77 78 79 80
            LEVEL2 'using hybrid vertical coordinates'
            call hybrid_coordinates(.true.)
STDERR 'coordinates(): hybrid_coordinates not coded yet'
stop
81
         case (_ADAPTIVE_COORDS_) ! adaptive vertical coordinates
82
            LEVEL2 'using adaptive vertical coordinates'
83
            call adaptive_coordinates(.true.,hotstart)
84 85
         case default
      end select
gotm's avatar
gotm committed
86
      first = .false.
87
   else
88 89
      select case (vert_cord)
         case (_SIGMA_COORDS_) ! sigma coordinates
90
            call sigma_coordinates(.false.)
91 92
         case (_Z_COORDS_) ! z-level
         case (_GENERAL_COORDS_) ! general vertical coordinates
93
            call general_coordinates(.false.,cord_relax,maxdepth)
94
         case (_HYBRID_COORDS_) ! hybrid vertical coordinates
95
            call hybrid_coordinates(.false.)
96 97
         case (_ADAPTIVE_COORDS_) ! adaptive vertical coordinates
            call adaptive_coordinates(.false.,hotstart)
gotm's avatar
gotm committed
98 99 100 101
         case default
      end select
   end if ! first

Knut's avatar
Knut committed
102 103
   hvel = _HALF_ * ( ho + hn )

Knut's avatar
Knut committed
104 105 106
   ! calculate the z-coordinate of the cell centers
   ! references to mean sea level
   zc(:,:,0)=-H(:,:)
Knut's avatar
Knut committed
107
   zc(:,:,1)=-H(:,:) + _HALF_*hn(:,:,1)
Knut's avatar
Knut committed
108
   do k=2,kmax
Knut's avatar
Knut committed
109
      zc(:,:,k)=zc(:,:,k-1)+_HALF_*(hn(:,:,k-1)+hn(:,:,k))
Knut's avatar
Knut committed
110 111
   end do

112
#ifdef SLICE_MODEL
113
   do i=imin,imax
114 115 116 117 118 119 120 121 122
      do k=kvmin(i,2),kmax
         hvo(i,1,k)=hvo(i,2,k)
         hvo(i,3,k)=hvo(i,2,k)
         hvn(i,1,k)=hvn(i,2,k)
         hvn(i,3,k)=hvn(i,2,k)
      end do
   end do
#endif

bjb's avatar
bjb committed
123
   call toc(TIM_COORDS)
gotm's avatar
gotm committed
124
#ifdef DEBUG
125
   write(debug,*) 'Leaving coordinates()'
gotm's avatar
gotm committed
126 127 128 129 130 131 132
   write(debug,*)
#endif
   return
   end subroutine coordinates
!EOC

!-----------------------------------------------------------------------
133
! Copyright (C) 2007 - Hans Burchard and Karsten Bolding               !
gotm's avatar
gotm committed
134
!-----------------------------------------------------------------------