- Timestamp:
- 2011-11-15T21:55:40+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_NEMO_MERGE_2011/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn.F90
r2528 r3116 15 15 !! 'key_bdy' : Unstructured Open Boundary Condition 16 16 !!---------------------------------------------------------------------- 17 !! bdy_dyn _frs : relaxation of velocities on unstructured open boundary18 !! bdy_dyn _fla : Flather condition for barotropic solution17 !! bdy_dyn3d : apply open boundary conditions to baroclinic velocities 18 !! bdy_dyn3d_frs : apply Flow Relaxation Scheme 19 19 !!---------------------------------------------------------------------- 20 20 USE oce ! ocean dynamics and tracers 21 21 USE dom_oce ! ocean space and time domain 22 USE dynspg_oce 22 23 USE bdy_oce ! ocean open boundary conditions 23 USE dynspg_oce ! for barotropic variables24 USE phycst ! physical constants24 USE bdydyn2d ! open boundary conditions for barotropic solution 25 USE bdydyn3d ! open boundary conditions for baroclinic velocities 25 26 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 26 USE bdytides ! for tidal harmonic forcing at boundary27 27 USE in_out_manager ! 28 28 … … 30 30 PRIVATE 31 31 32 PUBLIC bdy_dyn_frs ! routine called in dynspg_flt (free surface case ONLY) 33 # if defined key_dynspg_exp || defined key_dynspg_ts 34 PUBLIC bdy_dyn_fla ! routine called in dynspg_flt (free surface case ONLY) 35 # endif 32 PUBLIC bdy_dyn ! routine called in dynspg_flt (if lk_dynspg_flt) or 33 ! dyn_nxt (if lk_dynspg_ts or lk_dynspg_exp) 36 34 35 # include "domzgr_substitute.h90" 37 36 !!---------------------------------------------------------------------- 38 37 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 42 41 CONTAINS 43 42 44 SUBROUTINE bdy_dyn _frs( kt)43 SUBROUTINE bdy_dyn( kt, dyn3d_only ) 45 44 !!---------------------------------------------------------------------- 46 !! *** SUBROUTINE bdy_dyn _frs***45 !! *** SUBROUTINE bdy_dyn *** 47 46 !! 48 !! ** Purpose : - Apply the Flow Relaxation Scheme for dynamic in the 49 !! case of unstructured open boundaries. 47 !! ** Purpose : - Wrapper routine for bdy_dyn2d and bdy_dyn3d. 50 48 !! 51 !! References :- Engedahl H., 1995: Use of the flow relaxation scheme in52 !! a three-dimensional baroclinic ocean model with realistic53 !! topography. Tellus, 365-382.54 49 !!---------------------------------------------------------------------- 55 INTEGER, INTENT( in ) :: kt ! Main time step counter 50 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 51 USE wrk_nemo, ONLY: wrk_2d_7, wrk_2d_8 ! 2D workspace 56 52 !! 57 INTEGER :: jb, jk ! dummy loop indices 58 INTEGER :: ii, ij, igrd ! local integers 59 REAL(wp) :: zwgt ! boundary weight 60 !!---------------------------------------------------------------------- 61 ! 62 IF(ln_dyn_frs) THEN ! If this is false, then this routine does nothing. 63 ! 64 IF( kt == nit000 ) THEN 65 IF(lwp) WRITE(numout,*) 66 IF(lwp) WRITE(numout,*) 'bdy_dyn_frs : Flow Relaxation Scheme on momentum' 67 IF(lwp) WRITE(numout,*) '~~~~~~~' 68 ENDIF 69 ! 70 igrd = 2 ! Relaxation of zonal velocity 71 DO jb = 1, nblen(igrd) 72 DO jk = 1, jpkm1 73 ii = nbi(jb,igrd) 74 ij = nbj(jb,igrd) 75 zwgt = nbw(jb,igrd) 76 ua(ii,ij,jk) = ( ua(ii,ij,jk) * ( 1.- zwgt ) + ubdy(jb,jk) * zwgt ) * umask(ii,ij,jk) 77 END DO 78 END DO 79 ! 80 igrd = 3 ! Relaxation of meridional velocity 81 DO jb = 1, nblen(igrd) 82 DO jk = 1, jpkm1 83 ii = nbi(jb,igrd) 84 ij = nbj(jb,igrd) 85 zwgt = nbw(jb,igrd) 86 va(ii,ij,jk) = ( va(ii,ij,jk) * ( 1.- zwgt ) + vbdy(jb,jk) * zwgt ) * vmask(ii,ij,jk) 87 END DO 88 END DO 89 CALL lbc_lnk( ua, 'U', -1. ) ; CALL lbc_lnk( va, 'V', -1. ) ! Boundary points should be updated 90 ! 91 ENDIF ! ln_dyn_frs 92 ! 93 END SUBROUTINE bdy_dyn_frs 53 INTEGER, INTENT( in ) :: kt ! Main time step counter 54 LOGICAL, INTENT( in ), OPTIONAL :: dyn3d_only ! T => only update baroclinic velocities 55 !! 56 INTEGER :: jk,ii,ij,ib,igrd ! Loop counter 57 LOGICAL :: ll_dyn2d, ll_dyn3d 58 !! 94 59 60 IF(wrk_in_use(2, 7,8) ) THEN 61 CALL ctl_stop('bdy_dyn: ERROR: requested workspace arrays are unavailable.') ; RETURN 62 END IF 95 63 96 # if defined key_dynspg_exp || defined key_dynspg_ts 97 !!---------------------------------------------------------------------- 98 !! 'key_dynspg_exp' OR explicit sea surface height 99 !! 'key_dynspg_ts ' split-explicit sea surface height 100 !!---------------------------------------------------------------------- 101 102 !! Option to use Flather with dynspg_flt not coded yet... 64 ll_dyn2d = .true. 65 ll_dyn3d = .true. 103 66 104 SUBROUTINE bdy_dyn_fla( pssh ) 105 !!---------------------------------------------------------------------- 106 !! *** SUBROUTINE bdy_dyn_fla *** 107 !! 108 !! - Apply Flather boundary conditions on normal barotropic velocities 109 !! (ln_dyn_fla=.true. or ln_tides=.true.) 110 !! 111 !! ** WARNINGS about FLATHER implementation: 112 !!1. According to Palma and Matano, 1998 "after ssh" is used. 113 !! In ROMS and POM implementations, it is "now ssh". In the current 114 !! implementation (tested only in the EEL-R5 conf.), both cases were unstable. 115 !! So I use "before ssh" in the following. 116 !! 117 !!2. We assume that the normal ssh gradient at the bdy is zero. As a matter of 118 !! fact, the model ssh just inside the dynamical boundary is used (the outside 119 !! ssh in the code is not updated). 120 !! 121 !! References: Flather, R. A., 1976: A tidal model of the northwest European 122 !! continental shelf. Mem. Soc. R. Sci. Liege, Ser. 6,10, 141-164. 123 !!---------------------------------------------------------------------- 124 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pssh 67 IF( PRESENT(dyn3d_only) ) THEN 68 IF( dyn3d_only ) ll_dyn2d = .false. 69 ENDIF 125 70 126 INTEGER :: jb, igrd ! dummy loop indices 127 INTEGER :: ii, ij, iim1, iip1, ijm1, ijp1 ! 2D addresses 128 REAL(wp) :: zcorr ! Flather correction 129 REAL(wp) :: zforc ! temporary scalar 130 !!---------------------------------------------------------------------- 71 !------------------------------------------------------- 72 ! Set pointers 73 !------------------------------------------------------- 131 74 132 ! ---------------------------------!133 ! Flather boundary conditions :!134 ! ---------------------------------!135 136 IF(ln_dyn_fla .OR. ln_tides) THEN ! If these are both false, then this routine does nothing.75 pssh => sshn 76 phur => hur 77 phvr => hvr 78 pu2d => wrk_2d_7 79 pv2d => wrk_2d_8 137 80 138 ! Fill temporary array with ssh data (here spgu): 139 igrd = 4 140 spgu(:,:) = 0.0 141 DO jb = 1, nblenrim(igrd) 142 ii = nbi(jb,igrd) 143 ij = nbj(jb,igrd) 144 IF( ln_dyn_fla ) spgu(ii, ij) = sshbdy(jb) 145 IF( ln_tides ) spgu(ii, ij) = spgu(ii, ij) + sshtide(jb) 146 END DO 147 ! 148 igrd = 5 ! Flather bc on u-velocity; 149 ! ! remember that flagu=-1 if normal velocity direction is outward 150 ! ! I think we should rather use after ssh ? 151 DO jb = 1, nblenrim(igrd) 152 ii = nbi(jb,igrd) 153 ij = nbj(jb,igrd) 154 iim1 = ii + MAX( 0, INT( flagu(jb) ) ) ! T pts i-indice inside the boundary 155 iip1 = ii - MIN( 0, INT( flagu(jb) ) ) ! T pts i-indice outside the boundary 156 ! 157 zcorr = - flagu(jb) * SQRT( grav * hur_e(ii, ij) ) * ( pssh(iim1, ij) - spgu(iip1,ij) ) 158 zforc = ubtbdy(jb) + utide(jb) 159 ua_e(ii,ij) = zforc + zcorr * umask(ii,ij,1) 160 END DO 161 ! 162 igrd = 6 ! Flather bc on v-velocity 163 ! ! remember that flagv=-1 if normal velocity direction is outward 164 DO jb = 1, nblenrim(igrd) 165 ii = nbi(jb,igrd) 166 ij = nbj(jb,igrd) 167 ijm1 = ij + MAX( 0, INT( flagv(jb) ) ) ! T pts j-indice inside the boundary 168 ijp1 = ij - MIN( 0, INT( flagv(jb) ) ) ! T pts j-indice outside the boundary 169 ! 170 zcorr = - flagv(jb) * SQRT( grav * hvr_e(ii, ij) ) * ( pssh(ii, ijm1) - spgu(ii,ijp1) ) 171 zforc = vbtbdy(jb) + vtide(jb) 172 va_e(ii,ij) = zforc + zcorr * vmask(ii,ij,1) 173 END DO 174 CALL lbc_lnk( ua_e, 'U', -1. ) ! Boundary points should be updated 175 CALL lbc_lnk( va_e, 'V', -1. ) ! 176 ! 177 ENDIF ! ln_dyn_fla .or. ln_tides 178 ! 179 END SUBROUTINE bdy_dyn_fla 180 #endif 81 !------------------------------------------------------- 82 ! Split velocities into barotropic and baroclinic parts 83 !------------------------------------------------------- 84 85 pu2d(:,:) = 0.e0 86 pv2d(:,:) = 0.e0 87 DO jk = 1, jpkm1 !! Vertically integrated momentum trends 88 pu2d(:,:) = pu2d(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ua(:,:,jk) 89 pv2d(:,:) = pv2d(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * va(:,:,jk) 90 END DO 91 pu2d(:,:) = pu2d(:,:) * phur(:,:) 92 pv2d(:,:) = pv2d(:,:) * phvr(:,:) 93 DO jk = 1 , jpkm1 94 ua(:,:,jk) = ua(:,:,jk) - pu2d(:,:) 95 va(:,:,jk) = va(:,:,jk) - pv2d(:,:) 96 END DO 97 98 !------------------------------------------------------- 99 ! Apply boundary conditions to barotropic and baroclinic 100 ! parts separately 101 !------------------------------------------------------- 102 103 IF( ll_dyn2d ) CALL bdy_dyn2d( kt ) 104 105 IF( ll_dyn3d ) CALL bdy_dyn3d( kt ) 106 107 !------------------------------------------------------- 108 ! Recombine velocities 109 !------------------------------------------------------- 110 111 DO jk = 1 , jpkm1 112 ua(:,:,jk) = ( ua(:,:,jk) + pu2d(:,:) ) * umask(:,:,jk) 113 va(:,:,jk) = ( va(:,:,jk) + pv2d(:,:) ) * vmask(:,:,jk) 114 END DO 115 116 IF(wrk_not_released(2, 7,8) ) CALL ctl_stop('bdy_dyn: ERROR: failed to release workspace arrays.') 117 118 END SUBROUTINE bdy_dyn 181 119 182 120 #else … … 185 123 !!---------------------------------------------------------------------- 186 124 CONTAINS 187 SUBROUTINE bdy_dyn_frs( kt ) ! Empty routine 188 WRITE(*,*) 'bdy_dyn_frs: You should not have seen this print! error?', kt 189 END SUBROUTINE bdy_dyn_frs 190 SUBROUTINE bdy_dyn_fla( pssh ) ! Empty routine 191 REAL :: pssh(:,:) 192 WRITE(*,*) 'bdy_dyn_fla: You should not have seen this print! error?', pssh(1,1) 193 END SUBROUTINE bdy_dyn_fla 125 SUBROUTINE bdy_dyn( kt ) ! Empty routine 126 WRITE(*,*) 'bdy_dyn: You should not have seen this print! error?', kt 127 END SUBROUTINE bdy_dyn 194 128 #endif 195 129
Note: See TracChangeset
for help on using the changeset viewer.