bdy_3d.F90 6.54 KB
Newer Older
kbk's avatar
kbk committed
1
!$Id: bdy_3d.F90,v 1.10 2006-03-01 15:54:08 kbk Exp $
gotm's avatar
gotm committed
2 3 4 5
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
hb's avatar
hb committed
6
! !MODULE:  bdy_3d - 3D boundary conditions \label{bdy-3d}
gotm's avatar
gotm committed
7 8 9 10 11
!
! !INTERFACE:
   module bdy_3d
!
! !DESCRIPTION:
12 13 14
!  
! Here, the three-dimensional boundary
! conditions for temperature and salinity are handled.
gotm's avatar
gotm committed
15 16
!
! !USES:
kbk's avatar
kbk committed
17
   use halo_zones, only : H_TAG,U_TAG,V_TAG
gotm's avatar
gotm committed
18 19
   use domain, only: imin,jmin,imax,jmax,H,az,au,av
   use domain, only: iimin,jjmin,iimax,jjmax,kmax
kbk's avatar
kbk committed
20
   use domain, only: nsbv,NWB,NNB,NEB,NSB,bdy_index
gotm's avatar
gotm committed
21 22 23 24 25 26 27 28
   use domain, only: wi,wfj,wlj,nj,nfi,nli,ei,efj,elj,sj,sfi,sli
   use variables_3d
   IMPLICIT NONE
!
   private
!
! !PUBLIC DATA MEMBERS:
   public init_bdy_3d, do_bdy_3d
kbk's avatar
kbk committed
29
   REALTYPE, public, allocatable       :: S_bdy(:,:),T_bdy(:,:)
gotm's avatar
gotm committed
30 31 32
!
! !PRIVATE DATA MEMBERS:
!
kbk's avatar
kbk committed
33 34 35
! !REVISION HISTORY:
!  Original author(s): Karsten Bolding & Hans Burchard
!
gotm's avatar
gotm committed
36 37 38 39 40 41 42 43 44
! !LOCAL VARIABLES:
!EOP
!-----------------------------------------------------------------------

   contains

!-----------------------------------------------------------------------
!BOP
!
hb's avatar
hb committed
45 46
! !IROUTINE: init_bdy_3d - initialising 3D boundary conditions
! \label{sec-init-bdy-3d}
gotm's avatar
gotm committed
47 48 49 50 51
!
! !INTERFACE:
   subroutine init_bdy_3d()
!
! !DESCRIPTION:
52 53 54
!
! Here, the necessary fields {\tt S\_bdy} and {\tt T\_bdy} for
! salinity and temperature, respectively, are allocated.
gotm's avatar
gotm committed
55 56 57 58 59 60 61 62 63 64 65
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
!
! !INPUT/OUTPUT PARAMETERS:
!
! !OUTPUT PARAMETERS:
!
! !LOCAL VARIABLES:
kbk's avatar
kbk committed
66
   integer                   :: rc,i,j,k,n
gotm's avatar
gotm committed
67 68 69 70 71 72 73 74 75 76
!EOP
!-------------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'init_bdy_3d() # ',Ncall
#endif

   LEVEL2 'init_bdy_3d()'
77
   allocate(S_bdy(0:kmax,nsbv),stat=rc)
gotm's avatar
gotm committed
78 79
   if (rc /= 0) stop 'init_init_bdy_3d: Error allocating memory (S_bdy)'

80
   allocate(T_bdy(0:kmax,nsbv),stat=rc)
gotm's avatar
gotm committed
81 82 83 84 85 86 87 88 89 90 91 92 93
   if (rc /= 0) stop 'init_init_bdy_3d: Error allocating memory (T_bdy)'

#ifdef DEBUG
   write(debug,*) 'Leaving init_bdy_3d()'
   write(debug,*)
#endif
   return
   end subroutine init_bdy_3d
!EOC

!-----------------------------------------------------------------------
!BOP
!
94
! !IROUTINE:  do_bdy_3d  - updating 3D boundary conditions
hb's avatar
hb committed
95
! \label{sec-do-bdy-3d}
gotm's avatar
gotm committed
96 97 98 99 100
!
! !INTERFACE:
   subroutine do_bdy_3d(tag,field)
!
! !DESCRIPTION:
101 102 103 104
!
! Here, the boundary conditions for salinity and temperature are
! copied to the boundary points and relaxed to the near boundary points
! by means of the flow relaxation scheme by \cite{MARTINSENea87}.
gotm's avatar
gotm committed
105 106 107 108 109
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
kbk's avatar
kbk committed
110
   integer, intent(in)                 :: tag
gotm's avatar
gotm committed
111 112
!
! !INPUT/OUTPUT PARAMETERS:
kbk's avatar
kbk committed
113
   REALTYPE, intent(inout)             :: field(I3DFIELD)
gotm's avatar
gotm committed
114 115 116 117
!
! !OUTPUT PARAMETERS:
!
! !LOCAL VARIABLES:
kbk's avatar
kbk committed
118
   integer                   :: i,j,k,l,n,ii,jj
kbk's avatar
kbk committed
119
   REALTYPE                  :: sp(1:4),rat
gotm's avatar
gotm committed
120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'do_bdy_3d() # ',Ncall
#endif

#if 0
   select case (tag)
      case (1)
      ! Lateral zero-gradient boundary condition (north & south)
         do k=1,kmax
            do i=iimin,iimax
               if (au(i,jjmin) .eq. 3) field(i,jjmin,k)=field(i,jjmin+1,k)
               if (au(i,jjmax) .eq. 3) field(i,jjmax,k)=field(i,jjmax-1,k)
            end do
         end do
      case (2)
      ! Lateral zero-gradient boundary conditions (west & east)
         do k=1,kmax
            do j=jjmin,jjmax
               if (av(iimin,j) .eq. 3) field(iimin,j,k)=field(iimin+1,j,k)
               if (av(iimax,j) .eq. 3) field(iimax,j,k)=field(iimax-1,j,k)
            end do
         end do
      case default
         FATAL 'Non valid tag'
kbk's avatar
kbk committed
149
         stop 'do_bdy_3d'
gotm's avatar
gotm committed
150 151 152
   end select
#endif

kbk's avatar
kbk committed
153
#ifndef NO_BAROCLINIC
gotm's avatar
gotm committed
154 155 156 157 158 159
! Sponge layer factors according to Martinsen and Engedahl, 1987.
   sp(1)=1.0
   sp(2)=0.5625
   sp(3)=0.25
   sp(4)=0.0625

kbk's avatar
kbk committed
160
   l = 0
gotm's avatar
gotm committed
161 162
   k = 0
   do n=1,NWB
kbk's avatar
kbk committed
163 164
      l = l+1
      k = bdy_index(l)
gotm's avatar
gotm committed
165
      i = wi(n)
gotm's avatar
gotm committed
166
      do j=wfj(n),wlj(n)
gotm's avatar
gotm committed
167
         do ii=1,4
kbk's avatar
kbk committed
168
            if (az(i-1+ii,j).gt.0) then
169 170
               S(i-1+ii,j,:) = sp(ii)*S_bdy(:,k)+(1.-sp(ii))*S(i-1+ii,j,:)
               T(i-1+ii,j,:) = sp(ii)*T_bdy(:,k)+(1.-sp(ii))*T(i-1+ii,j,:)
kbk's avatar
kbk committed
171
            end if
172 173 174 175 176
#ifdef GETM_BIO
            if ( allocated(cc3d) ) then
               cc3d(:,i,j,:) = cc3d(:,i+1,j,:)
            end if
#endif
gotm's avatar
gotm committed
177
         end do
kbk's avatar
kbk committed
178
         k = k+1
gotm's avatar
gotm committed
179 180
      end do
   end do
kbk's avatar
kbk committed
181

gotm's avatar
gotm committed
182
   do n = 1,NNB
kbk's avatar
kbk committed
183 184
      l = l+1
      k = bdy_index(l)
gotm's avatar
gotm committed
185 186 187
      j = nj(n)
      do i = nfi(n),nli(n)
         do jj=1,4
kbk's avatar
kbk committed
188
            if (az(i,j+1-jj).gt.0) then
189 190
               S(i,j+1-jj,:) = sp(jj)*S_bdy(:,k)+(1.-sp(jj))*S(i,j+1-jj,:)
               T(i,j+1-jj,:) = sp(jj)*T_bdy(:,k)+(1.-sp(jj))*T(i,j+1-jj,:)
kbk's avatar
kbk committed
191
            end if
192 193 194 195 196
#ifdef GETM_BIO
            if ( allocated(cc3d) ) then
               cc3d(:,i,j,:) = cc3d(:,i,j-1,:)
            end if
#endif
gotm's avatar
gotm committed
197
         end do
kbk's avatar
kbk committed
198
         k = k+1
gotm's avatar
gotm committed
199 200
      end do
   end do
kbk's avatar
kbk committed
201

gotm's avatar
gotm committed
202
   do n=1,NEB
kbk's avatar
kbk committed
203 204
      l = l+1
      k = bdy_index(l)
gotm's avatar
gotm committed
205
      i = ei(n)
gotm's avatar
gotm committed
206
      do j=efj(n),elj(n)
gotm's avatar
gotm committed
207
         do ii=1,4
kbk's avatar
kbk committed
208
            if (az(i+1-ii,j).gt.0) then
209 210
               S(i+1-ii,j,:) = sp(ii)*S_bdy(:,k)+(1.-sp(ii))*S(i+1-ii,j,:)
               T(i+1-ii,j,:) = sp(ii)*T_bdy(:,k)+(1.-sp(ii))*T(i+1-ii,j,:)
kbk's avatar
kbk committed
211
            end if
212 213 214 215 216
#ifdef GETM_BIO
            if ( allocated(cc3d) ) then
               cc3d(:,i,j,:) = cc3d(:,i-1,j,:)
            end if
#endif
gotm's avatar
gotm committed
217
         end do
kbk's avatar
kbk committed
218
         k = k+1
gotm's avatar
gotm committed
219 220
      end do
   end do
kbk's avatar
kbk committed
221

gotm's avatar
gotm committed
222
   do n = 1,NSB
kbk's avatar
kbk committed
223 224
      l = l+1
      k = bdy_index(l)
gotm's avatar
gotm committed
225 226 227
      j = sj(n)
      do i = sfi(n),sli(n)
         do jj=1,4
kbk's avatar
kbk committed
228
            if (az(i,j-1+jj).gt.0) then
229 230
               S(i,j-1+jj,:) = sp(jj)*S_bdy(:,k)+(1.-sp(jj))*S(i,j-1+jj,:)
               T(i,j-1+jj,:) = sp(jj)*T_bdy(:,k)+(1.-sp(jj))*T(i,j-1+jj,:)
kbk's avatar
kbk committed
231
            end if
232 233 234 235 236
#ifdef GETM_BIO
            if ( allocated(cc3d) ) then
               cc3d(:,i,j,:) = cc3d(:,i,j+1,:)
            end if
#endif
gotm's avatar
gotm committed
237
         end do
kbk's avatar
kbk committed
238
         k = k+1
gotm's avatar
gotm committed
239 240
      end do
   end do
kbk's avatar
kbk committed
241 242 243

   call mirror_bdy_3d(T,H_TAG)
   call mirror_bdy_3d(S,H_TAG)
gotm's avatar
gotm committed
244

245 246 247 248 249 250 251 252 253
#ifdef GETM_BIO
   if ( allocated(cc3d) ) then
      do n=1,size(cc3d,1)
         call mirror_bdy_3d(cc3d(n,:,:,:),H_TAG)
      end do
   end if
#endif


kbk's avatar
kbk committed
254 255
#endif

gotm's avatar
gotm committed
256 257 258 259 260 261 262 263 264 265 266 267 268 269 270
#ifdef DEBUG
   write(debug,*) 'leaving do_bdy_3d()'
   write(debug,*)
#endif
   return
   end subroutine do_bdy_3d
!EOC

!-----------------------------------------------------------------------

   end module bdy_3d

!-----------------------------------------------------------------------
! Copyright (C) 2001 - Karsten Bolding and Hans Burchard               !
!-----------------------------------------------------------------------