Changeset 2392 for branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/cla.F90
- Timestamp:
- 2010-11-15T22:20:05+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/cla.F90
r2287 r2392 1 1 MODULE cla 2 !!============================================================================== 3 !! *** MODULE cla *** 4 !! Cross Land Advection : parameterize ocean exchanges through straits by a 5 !! specified advection across land. 6 !!============================================================================== 2 !!====================================================================== 3 !! *** MODULE cla *** 4 !! Cross Land Advection : specific update of the horizontal divergence, 5 !! tracer trends and after velocity 6 !! 7 !! --- Specific to ORCA_R2 --- 8 !! 9 !!====================================================================== 10 !! History : 1.0 ! 2002-11 (A. Bozec) Original code 11 !! 3.2 ! 2009-07 (G. Madec) merge cla, cla_div, tra_cla, cla_dynspg 12 !! ! and correct a mpp bug reported by A.R. Porter 13 !!---------------------------------------------------------------------- 7 14 #if defined key_orca_r2 8 15 !!---------------------------------------------------------------------- 9 !! 'key_orca_r2' : ORCA R2 configuration16 !! 'key_orca_r2' global ocean model R2 10 17 !!---------------------------------------------------------------------- 11 !! tra_cla : update the tracer trend with the horizontal 12 !! and vertical advection trends at straits 13 !! tra_bab_el_mandeb : 14 !! tra_gibraltar : 15 !! tra_hormuz : 16 !! tra_cla_init : 18 !! cla_div : update of horizontal divergence at cla straits 19 !! tra_cla : update of tracers at cla straits 20 !! cla_dynspg : update of after horizontal velocities at cla straits 21 !! cla_init : initialisation - control check 22 !! cla_bab_el_mandeb : cross land advection for Bab-el-mandeb strait 23 !! cla_gibraltar : cross land advection for Gibraltar strait 24 !! cla_hormuz : cross land advection for Hormuz strait 17 25 !!---------------------------------------------------------------------- 18 !! * Modules used19 USE oce ! ocean dynamics and tracers variables20 USE dom_oce ! ocean space and time domain variables21 USE sbc_oce ! surface boundary condition: ocean22 USE in_out_manager 23 USE l bclnk ! ocean lateral boundary conditions (or mpp link)24 USE l ib_mpp ! distributed memory computing26 USE oce ! ocean dynamics and tracers 27 USE dom_oce ! ocean space and time domain 28 USE sbc_oce ! surface boundary condition: ocean 29 USE dynspg_oce ! ocean dynamics: surface pressure gradient variables 30 USE in_out_manager ! I/O manager 31 USE lib_mpp ! distributed memory computing library 32 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 25 33 26 34 IMPLICIT NONE 27 35 PRIVATE 28 29 !! * Routine accessibility 30 PUBLIC tra_cla ! routine called by step.F90 31 PUBLIC tra_cla_init ! routine called by opa.F90 32 33 !! * Modules variables 34 REAL(wp) :: zempmed, zempred 35 36 REAL(wp) :: zisw_rs, zurw_rs, zbrw_rs ! Imposed transport Red Sea 37 REAL(wp) :: zisw_ms, zmrw_ms, zurw_ms, zbrw_ms ! Imposed transport Med Sea 38 REAL(wp) :: zisw_pg, zbrw_pg ! Imposed transport Persic Gulf 39 40 REAL(wp), DIMENSION(jpk) :: & 41 zu1_rs_i, zu2_rs_i, zu3_rs_i, & ! Red Sea velocities 42 zu1_ms_i, zu2_ms_i, zu3_ms_i, & ! Mediterranean Sea velocities 43 zu_pg ! Persic Gulf velocities 44 REAL(wp), DIMENSION (jpk) :: zthor, zshor ! Temperature, salinity Hormuz 36 37 PUBLIC cla_init ! routine called by opa.F90 38 PUBLIC cla_div ! routine called by divcur.F90 39 PUBLIC cla_traadv ! routine called by traadv.F90 40 PUBLIC cla_dynspg ! routine called by dynspg_flt.F90 41 42 INTEGER :: nbab, ngib, nhor ! presence or not of required grid-point on local domain 43 ! ! for Bab-el-Mandeb, Gibraltar, and Hormuz straits 44 45 ! !!! profile of hdiv for some straits 46 REAL(wp), DIMENSION (jpk) :: hdiv_139_101, hdiv_139_101_kt ! Gibraltar strait, fixed & time evolving part (i,j)=(172,101) 47 REAL(wp), DIMENSION (jpk) :: hdiv_139_102 ! Gibraltar strait, fixed part only (i,j)=(139,102) 48 REAL(wp), DIMENSION (jpk) :: hdiv_141_102, hdiv_141_102_kt ! Gibraltar strait, fixed & time evolving part (i,j)=(141,102) 49 REAL(wp), DIMENSION (jpk) :: hdiv_161_88 , hdiv_161_88_kt ! Bab-el-Mandeb strait, fixed & time evolving part (i,j)=(161,88) 50 REAL(wp), DIMENSION (jpk) :: hdiv_161_87 ! Bab-el-Mandeb strait, fixed part only (i,j)=(161,87) 51 REAL(wp), DIMENSION (jpk) :: hdiv_160_89 , hdiv_160_89_kt ! Bab-el-Mandeb strait, fixed & time evolving part (i,j)=(160,89) 52 REAL(wp), DIMENSION (jpk) :: hdiv_172_94 ! Hormuz strait, fixed part only (i,j)=(172, 94) 53 54 REAL(wp), DIMENSION (jpk) :: t_171_94_hor, s_171_94_hor ! Temperature, salinity in the Hormuz strait 45 55 46 56 !! * Substitutions 47 57 # include "domzgr_substitute.h90" 48 # include "vectopt_loop_substitute.h90"49 58 !!---------------------------------------------------------------------- 50 59 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 51 60 !! $Id$ 52 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)61 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 53 62 !!---------------------------------------------------------------------- 54 55 63 CONTAINS 56 64 57 SUBROUTINE tra_cla( kt ) 65 SUBROUTINE cla_div( kt ) 66 !!---------------------------------------------------------------------- 67 !! *** ROUTINE div_cla *** 68 !! 69 !! ** Purpose : update the horizontal divergence of the velocity field 70 !! at some straits ( Gibraltar, Bab el Mandeb and Hormuz ). 71 !! 72 !! ** Method : - first time-step: initialisation of cla 73 !! - all time-step: using imposed transport at each strait, 74 !! the now horizontal divergence is updated 75 !! 76 !! ** Action : phdivn updted now horizontal divergence at cla straits 77 !!---------------------------------------------------------------------- 78 INTEGER, INTENT( in ) :: kt ! ocean time-step index 79 !!---------------------------------------------------------------------- 80 ! 81 IF( kt == nit000 ) THEN 82 ! 83 CALL cla_init ! control check 84 ! 85 IF(lwp) WRITE(numout,*) 86 IF(lwp) WRITE(numout,*) 'div_cla : cross land advection on hdiv ' 87 IF(lwp) WRITE(numout,*) '~~~~~~~~' 88 ! 89 IF( nbab == 1 ) CALL cla_bab_el_mandeb('ini') ! Bab el Mandeb ( Red Sea - Indian ocean ) 90 IF( ngib == 1 ) CALL cla_gibraltar ('ini') ! Gibraltar strait (Med Sea - Atlantic ocean) 91 IF( nhor == 1 ) CALL cla_hormuz ('ini') ! Hormuz Strait ( Persian Gulf - Indian ocean ) 92 ! 93 ENDIF 94 ! 95 IF( nbab == 1 ) CALL cla_bab_el_mandeb('div') ! Bab el Mandeb ( Red Sea - Indian ocean ) 96 IF( ngib == 1 ) CALL cla_gibraltar ('div') ! Gibraltar strait (Med Sea - Atlantic ocean) 97 IF( nhor == 1 ) CALL cla_hormuz ('div') ! Hormuz Strait ( Persian Gulf - Indian ocean ) 98 ! 99 !!gm lbc useless here, no? 100 !!gm CALL lbc_lnk( hdivn, 'T', 1. ) ! Lateral boundary conditions on hdivn 101 ! 102 END SUBROUTINE cla_div 103 104 105 SUBROUTINE cla_traadv( kt ) 58 106 !!---------------------------------------------------------------------- 59 107 !! *** ROUTINE tra_cla *** … … 63 111 !! at some straits ( Bab el Mandeb, Gibraltar, Hormuz ). 64 112 !! 65 !! ** Method : ... 66 !! Add this trend now to the general trend of tracer (ta,sa): 67 !! (ta,sa) = (ta,sa) + ( zta , zsa ) 68 !! 69 !! ** Action : update (ta,sa) with the now advective tracer trends 70 !! 71 !! History : 72 !! ! (A. Bozec) original code 73 !! 8.5 ! 02-11 (A. Bozec) F90: Free form and module 74 !!---------------------------------------------------------------------- 75 !! * Arguments 113 !! ** Method : using both imposed transport at each strait and T & S 114 !! budget, the now tracer trends is updated 115 !! 116 !! ** Action : (ta,sa) updated now tracer trends at cla straits 117 !!---------------------------------------------------------------------- 76 118 INTEGER, INTENT( in ) :: kt ! ocean time-step index 77 119 !!---------------------------------------------------------------------- 78 79 ! Bab el Mandeb strait horizontal advection 80 81 CALL tra_bab_el_mandeb 82 83 ! Gibraltar strait horizontal advection 84 85 CALL tra_gibraltar 86 87 ! Hormuz Strait ( persian Gulf) horizontal advection 88 89 CALL tra_hormuz 90 91 END SUBROUTINE tra_cla 92 93 94 SUBROUTINE tra_bab_el_mandeb 120 ! 121 IF( kt == nit000 ) THEN 122 IF(lwp) WRITE(numout,*) 123 IF(lwp) WRITE(numout,*) 'tra_cla : cross land advection on tracers ' 124 IF(lwp) WRITE(numout,*) '~~~~~~~~' 125 ENDIF 126 ! 127 IF( nbab == 1 ) CALL cla_bab_el_mandeb('tra') ! Bab el Mandeb strait 128 IF( ngib == 1 ) CALL cla_gibraltar ('tra') ! Gibraltar strait 129 IF( nhor == 1 ) CALL cla_hormuz ('tra') ! Hormuz Strait ( Persian Gulf) 130 ! 131 END SUBROUTINE cla_traadv 132 133 134 SUBROUTINE cla_dynspg( kt ) 135 !!---------------------------------------------------------------------- 136 !! *** ROUTINE cla_dynspg *** 137 !! 138 !! ** Purpose : Update the after velocity at some straits 139 !! (Bab el Mandeb, Gibraltar, Hormuz). 140 !! 141 !! ** Method : required to compute the filtered surface pressure gradient 142 !! 143 !! ** Action : (ua,va) after velocity at the cla straits 144 !!---------------------------------------------------------------------- 145 INTEGER, INTENT( in ) :: kt ! ocean time-step index 146 !!---------------------------------------------------------------------- 147 ! 148 IF( kt == nit000 ) THEN 149 IF(lwp) WRITE(numout,*) 150 IF(lwp) WRITE(numout,*) 'cla_dynspg : cross land advection on (ua,va) ' 151 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 152 ENDIF 153 ! 154 IF( nbab == 1 ) CALL cla_bab_el_mandeb('spg') ! Bab el Mandeb strait 155 IF( ngib == 1 ) CALL cla_gibraltar ('spg') ! Gibraltar strait 156 IF( nhor == 1 ) CALL cla_hormuz ('spg') ! Hormuz Strait ( Persian Gulf) 157 ! 158 !!gm lbc is needed here, not? 159 !!gm CALL lbc_lnk( hdivn, 'U', -1. ) ; CALL lbc_lnk( hdivn, 'V', -1. ) ! Lateral boundary conditions 160 ! 161 END SUBROUTINE cla_dynspg 162 163 164 SUBROUTINE cla_init 165 !! ------------------------------------------------------------------- 166 !! *** ROUTINE cla_init *** 167 !! 168 !! ** Purpose : control check for mpp computation 169 !! 170 !! ** Method : - All the strait grid-points must be inside one of the 171 !! local domain interior for the cla advection to work 172 !! properly in mpp (i.e. inside (2:jpim1,2:jpjm1) ). 173 !! Define the corresponding indicators (nbab, ngib, nhor) 174 !! - The profiles of cross-land fluxes are currently hard 175 !! coded for L31 levels. Stop if jpk/=31 176 !! 177 !! ** Action : nbab, ngib, nhor strait inside the local domain or not 95 178 !!--------------------------------------------------------------------- 96 !! *** ROUTINE tra_bab_el_mandeb *** 97 !! 98 !! ** Purpose : Update the horizontal advective trend of tracers 99 !! correction in Bab el Mandeb strait and 100 !! add it to the general trend of tracer equations. 179 REAL(wp) :: ztemp 180 !!--------------------------------------------------------------------- 181 ! 182 IF(lwp) WRITE(numout,*) 183 IF(lwp) WRITE(numout,*) 'cla_init : cross land advection initialisation ' 184 IF(lwp) WRITE(numout,*) '~~~~~~~~~' 185 ! 186 IF( .NOT.lk_dynspg_flt ) CALL ctl_stop( 'cla_init: Cross Land Advection works only with lk_dynspg_flt=T ' ) 187 ! 188 IF( lk_vvl ) CALL ctl_stop( 'cla_init: Cross Land Advection does not work with lk_vvl=T option' ) 189 ! 190 IF( jpk /= 31 ) CALL ctl_stop( 'cla_init: Cross Land Advection hard coded for ORCA_R2_L31' ) 191 ! 192 ! _|_______|_______|_ 193 ! 89 | |///////| 194 ! _|_______|_______|_ 195 ! ------------------------ ! 88 |///////| | 196 ! Bab el Mandeb strait ! _|_______|_______|_ 197 ! ------------------------ ! 87 |///////| | 198 ! _|_______|_______|_ 199 ! | 160 | 161 | 200 ! 201 ! The 6 Bab el Mandeb grid-points must be inside one of the interior of the 202 ! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1) 203 nbab = 0 204 IF( ( 1 <= mj0( 88) .AND. mj1( 89) <= jpj ) .AND. & !* (161,89), (161,88) and (161,88) on the local pocessor 205 & ( 1 <= mi0(160) .AND. mi1(161) <= jpi ) ) nbab = 1 206 ! 207 ! test if there is no local domain that includes all required grid-points 208 ztemp = REAL( nbab ) 209 IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value 210 IF( ztemp == 0 ) THEN ! Only 2 points in each direction, this should never be a problem 211 CALL ctl_stop( ' cross land advection at Bab-el_Mandeb does not work with your processor cutting: change it' ) 212 ENDIF 213 ! ___________________________ 214 ! ------------------------ ! 102 | |///////| | 215 ! Gibraltar strait ! _|_______|_______|_______|_ 216 ! ------------------------ ! 101 | |///////| | 217 ! _|_______|_______|_______|_ 218 ! | 139 | 140 | 141 | 219 ! 220 ! The 6 Gibraltar grid-points must be inside one of the interior of the 221 ! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1) 222 ngib = 0 223 IF( ( 2 <= mj0(101) .AND. mj1(102) <= jpjm1 ) .AND. & !* (139:141,101:102) on the local pocessor 224 & ( 2 <= mi0(139) .AND. mi1(141) <= jpim1 ) ) ngib = 1 225 ! 226 ! test if there is no local domain that includes all required grid-points 227 ztemp = REAL( ngib ) 228 IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value 229 IF( ztemp == 0 ) THEN ! 3 points in i-direction, this may be a problem with some cutting 230 CALL ctl_stop( ' cross land advection at Gibraltar does not work with your processor cutting: change it' ) 231 ENDIF 232 ! _______________ 233 ! ------------------------ ! 94 |/////| | 234 ! Hormuz strait ! _|_____|_____|_ 235 ! ------------------------ ! 171 172 236 ! 237 ! The 2 Hormuz grid-points must be inside one of the interior of the 238 ! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1) 239 nhor = 0 240 IF( 2 <= mj0( 94) .AND. mj1( 94) <= jpjm1 .AND. & 241 & 2 <= mi0(171) .AND. mi1(172) <= jpim1 ) nhor = 1 242 ! 243 ! test if there is no local domain that includes all required grid-points 244 ztemp = REAL( nhor ) 245 IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value 246 IF( ztemp == 0 ) THEN ! 3 points in i-direction, this may be a problem with some cutting 247 CALL ctl_stop( ' cross land advection at Hormuz does not work with your processor cutting: change it' ) 248 ENDIF 249 ! 250 END SUBROUTINE cla_init 251 252 253 SUBROUTINE cla_bab_el_mandeb( cd_td ) 254 !!---------------------------------------------------------------------- 255 !! *** ROUTINE cla_bab_el_mandeb *** 256 !! 257 !! ** Purpose : update the now horizontal divergence, the tracer tendancy 258 !! and the after velocity in vicinity of Bab el Mandeb ( Red Sea - Indian ocean). 259 !! 260 !! ** Method : compute the exchanges at each side of the strait : 261 !! 262 !! surf. zio_flow 263 !! (+ balance of emp) /\ |\\\\\\\\\\\| 264 !! || |\\\\\\\\\\\| 265 !! deep zio_flow || |\\\\\\\\\\\| 266 !! | || || |\\\\\\\\\\\| 267 !! 89 | || || |\\\\\\\\\\\| 268 !! |__\/_v_||__|____________ 269 !! !\\\\\\\\\\\| surf. zio_flow 270 !! |\\\\\\\\\\\|<=== (+ balance of emp) 271 !! |\\\\\\\\\\\u 272 !! 88 |\\\\\\\\\\\|<--- deep zrecirc (upper+deep at 2 different levels) 273 !! |___________|__________ 274 !! !\\\\\\\\\\\| 275 !! |\\\\\\\\\\\| ---\ deep zrecirc (upper+deep) 276 !! 87 !\\\\\\\\\\\u ===/ + deep zio_flow (all at the same level) 277 !! !\\\\\\\\\\\| 278 !! !___________|__________ 279 !! 160 161 280 !! 281 !!---------------------------------------------------------------------- 282 CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence 283 ! ! ='tra' update the tracers 284 ! ! ='spg' update after velocity 285 INTEGER :: ji, jj, jk ! dummy loop indices 286 REAL(wp) :: zemp_red ! temporary scalar 287 REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot 288 !!--------------------------------------------------------------------- 289 ! 290 SELECT CASE( cd_td ) 291 ! ! ---------------- ! 292 CASE( 'ini' ) ! initialisation ! 293 ! ! ---------------- ! 294 ! 295 zio_flow = 0.4e6 ! imposed in/out flow 296 zrecirc_upp = 0.2e6 ! imposed upper recirculation water 297 zrecirc_bot = 0.5e6 ! imposed bottom recirculation water 298 299 hdiv_161_88(:) = 0.e0 ! (161,88) Gulf of Aden side, north point 300 hdiv_161_87(:) = 0.e0 ! (161,87) Gulf of Aden side, south point 301 hdiv_160_89(:) = 0.e0 ! (160,89) Red sea side 302 303 DO jj = mj0(88), mj1(88) !** profile of hdiv at (161,88) (Gulf of Aden side, north point) 304 DO ji = mi0(161), mi1(161) !------------------------------ 305 DO jk = 1, 8 ! surface in/out flow (Ind -> Red) (div >0) 306 hdiv_161_88(jk) = + zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 307 END DO 308 ! ! recirculation water (Ind -> Red) (div >0) 309 hdiv_161_88(20) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,20) ) 310 hdiv_161_88(21) = + ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) ) 311 END DO 312 END DO 313 ! 314 DO jj = mj0(87), mj1(87) !** profile of hdiv at (161,88) (Gulf of Aden side, south point) 315 DO ji = mi0(161), mi1(161) !------------------------------ 316 ! ! deep out flow + recirculation (Red -> Ind) (div <0) 317 hdiv_161_87(21) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) ) 318 END DO 319 END DO 320 ! 321 DO jj = mj0(89), mj1(89) !** profile of hdiv at (161,88) (Red sea side) 322 DO ji = mi0(160), mi1(160) !------------------------------ 323 DO jk = 1, 8 ! surface inflow (Ind -> Red) (div <0) 324 hdiv_160_89(jk) = - zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 325 END DO 326 ! ! deep outflow (Red -> Ind) (div >0) 327 hdiv_160_89(16) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,16) ) 328 END DO 329 END DO 330 ! ! ---------------- ! 331 CASE( 'div' ) ! update hdivn ! (call by divcur module) 332 ! ! ---------=====-- ! 333 ! !** emp on the Red Sea (div >0) 334 zemp_red = 0.e0 !--------------------- 335 DO jj = mj0(87), mj1(96) ! sum over the Red sea 336 DO ji = mi0(148), mi1(160) 337 zemp_red = zemp_red + emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj) 338 END DO 339 END DO 340 IF( lk_mpp ) CALL mpp_sum( zemp_red ) ! sum with other processors value 341 zemp_red = zemp_red * 1.e-3 ! convert in m3 342 ! 343 ! !** Correct hdivn (including emp adjustment) 344 ! !------------------------------------------- 345 DO jj = mj0(88), mj1(88) !* profile of hdiv at (161,88) (Gulf of Aden side, north point) 346 DO ji = mi0(161), mi1(161) 347 hdiv_161_88_kt(:) = hdiv_161_88(:) 348 DO jk = 1, 8 ! increase the inflow from the Indian (div >0) 349 hdiv_161_88_kt(jk) = hdiv_161_88(jk) + zemp_red / (8. * e2u(ji,jj) * fse3u(ji,jj,jk) ) 350 END DO 351 hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_161_88_kt(:) 352 END DO 353 END DO 354 DO jj = mj0(87), mj1(87) !* profile of divergence at (161,87) (Gulf of Aden side, south point) 355 DO ji = mi0(161), mi1(161) 356 hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_161_87(:) 357 END DO 358 END DO 359 DO jj = mj0(89), mj1(89) !* profile of divergence at (160,89) (Red sea side) 360 DO ji = mi0(160), mi1(160) 361 hdiv_160_89_kt(:) = hdiv_160_89(:) 362 DO jk = 1, 18 ! increase the inflow from the Indian (div <0) 363 hdiv_160_89_kt(jk) = hdiv_160_89(jk) - zemp_red / (10. * e1v(ji,jj) * fse3v(ji,jj,jk) ) 364 END DO 365 hdivn(ji, jj,:) = hdivn(ji, jj,:) + hdiv_160_89_kt(:) 366 END DO 367 END DO 368 ! ! ---------------- ! 369 CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module) 370 ! ! --------=======- ! 371 ! 372 DO jj = mj0(88), mj1(88) !** (161,88) (Gulf of Aden side, north point) 373 DO ji = mi0(161), mi1(161) 374 DO jk = 1, jpkm1 ! surf inflow + reciculation (from Gulf of Aden) 375 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_161_88_kt(jk) * tn(ji,jj,jk) 376 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_161_88_kt(jk) * sn(ji,jj,jk) 377 END DO 378 END DO 379 END DO 380 DO jj = mj0(87), mj1(87) !** (161,87) (Gulf of Aden side, south point) 381 DO ji = mi0(161), mi1(161) 382 jk = 21 ! deep outflow + recirulation (combined flux) 383 ta(ji,jj,jk) = ta(ji,jj,jk) + hdiv_161_88(20) * tn(ji ,jj+1,20) & ! upper recirculation from Gulf of Aden 384 & + hdiv_161_88(21) * tn(ji ,jj+1,21) & ! deep recirculation from Gulf of Aden 385 & + hdiv_160_89(16) * tn(ji-1,jj+2,16) ! deep inflow from Red sea 386 sa(ji,jj,jk) = sa(ji,jj,jk) + hdiv_161_88(20) * sn(ji ,jj+1,20) & 387 & + hdiv_161_88(21) * sn(ji ,jj+1,21) & 388 & + hdiv_160_89(16) * sn(ji-1,jj+2,16) 389 END DO 390 END DO 391 DO jj = mj0(89), mj1(89) !** (161,88) (Red sea side) 392 DO ji = mi0(160), mi1(160) 393 DO jk = 1, 14 ! surface inflow (from Gulf of Aden) 394 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_160_89_kt(jk) * tn(ji+1,jj-1,jk) 395 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_160_89_kt(jk) * sn(ji+1,jj-1,jk) 396 END DO 397 ! ! deep outflow (from Red sea) 398 ta(ji,jj,16) = ta(ji,jj,16) - hdiv_160_89(jk) * tn(ji,jj,jk) 399 sa(ji,jj,16) = sa(ji,jj,16) - hdiv_160_89(jk) * sn(ji,jj,jk) 400 END DO 401 END DO 402 ! 403 ! ! ---------------- ! 404 CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module) 405 ! ! --------=======- ! 406 ! at this stage, (ua,va) are the after velocity, not the tendancy 407 ! compute the velocity from the divergence at T-point 408 ! 409 DO jj = mj0(88), mj1(88) !** (160,88) (Gulf of Aden side, north point) 410 DO ji = mi0(160), mi1(160) ! 160, not 161 as it is a U-point) 411 ua(ji,jj,:) = - hdiv_161_88_kt(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) & 412 & * e2u(ji,jj) * fse3u(ji,jj,:) 413 END DO 414 END DO 415 DO jj = mj0(87), mj1(87) !** (160,87) (Gulf of Aden side, south point) 416 DO ji = mi0(160), mi1(160) ! 160, not 161 as it is a U-point) 417 ua(ji,jj,:) = - hdiv_161_87(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) & 418 & * e2u(ji,jj) * fse3u(ji,jj,:) 419 END DO 420 END DO 421 DO jj = mj0(88), mj1(88) !** profile of divergence at (160,89) (Red sea side) 422 DO ji = mi0(160), mi1(160) ! 88, not 89 as it is a V-point) 423 va(ji,jj,:) = - hdiv_160_89_kt(:) / ( e1t(ji,jj+1) * e2t(ji,jj+1) * fse3t(ji,jj+1,:) ) & 424 & * e1v(ji,jj) * fse3v(ji,jj,:) 425 END DO 426 END DO 427 END SELECT 428 ! 429 END SUBROUTINE cla_bab_el_mandeb 430 431 432 SUBROUTINE cla_gibraltar( cd_td ) 433 !! ------------------------------------------------------------------- 434 !! *** ROUTINE cla_gibraltar *** 435 !! 436 !! ** Purpose : update the now horizontal divergence, the tracer 437 !! tendancyand the after velocity in vicinity of Gibraltar 438 !! strait ( Persian Gulf - Indian ocean ). 101 439 !! 102 440 !! ** Method : 103 !! We impose transport at Bab el Mandeb and knowing T and S in 104 !! surface and depth at each side of the strait, we deduce T and S 105 !! of the deep outflow of the Red Sea in the Indian ocean . 106 !! | 107 !! |/ \| N |\ /| 108 !! |_|_|______ | |___|______ 109 !! 88 | |<- W - - E 88 | |<- 110 !! 87 |___|______ | 87 |___|->____ 111 !! 160 161 S 160 161 112 !! horizontal view horizontal view 113 !! surface depth 114 !! 115 !! The horizontal advection is evaluated by a second order cen- 116 !! tered scheme using now fields (leap-frog scheme). In specific 117 !! areas (vicinity of major river mouths, some straits, or tn 118 !! approaching the freezing point) it is mixed with an upstream 119 !! scheme for stability reasons. 120 !! 121 !! C A U T I O N : the trend saved is the centered trend only. 122 !! It doesn't take into account the upstream part of the scheme. 123 !! 124 !! ** history : 125 !! ! 02-11 (A. Bozec) Original code 126 !! 8.5 ! 02-11 (A. Bozec) F90: Free form and module 441 !! _______________________ 442 !! deep zio_flow /====|///////|====> surf. zio_flow 443 !! + deep zrecirc \----|///////| (+balance of emp) 444 !! 102 u///////u 445 !! mid. recicul <--|///////|<==== deep zio_flow 446 !! _____|_______|_____ 447 !! surf. zio_flow ====>|///////| 448 !! (+balance of emp) |///////| 449 !! 101 u///////| 450 !! mid. recicul -->|///////| Caution: zrecirc split into 451 !! deep zrecirc ---->|///////| upper & bottom recirculation 452 !! _______|_______|_______ 453 !! 139 140 141 454 !! 127 455 !!--------------------------------------------------------------------- 128 !! * Local declarations 129 INTEGER :: ji, jj, jk ! dummy loop indices 130 REAL(wp) :: zsu, zvt 131 REAL(wp) :: zsumt, zsumt1, zsumt2, zsumt3, zsumt4 132 REAL(wp) :: zsums, zsums1, zsums2, zsums3, zsums4 133 REAL(wp) :: zt, zs 134 REAL(wp) :: zwei 135 REAL(wp), DIMENSION (jpk) :: zu1_rs, zu2_rs, zu3_rs 456 CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence 457 ! ! ='tra' update the tracers 458 ! ! ='spg' update after velocity 459 INTEGER :: ji, jj, jk ! dummy loop indices 460 REAL(wp) :: zemp_med ! temporary scalar 461 REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot 136 462 !!--------------------------------------------------------------------- 137 138 ! Initialization of vertical sum for T and S transport 139 ! ---------------------------------------------------- 140 141 zsumt = 0.e0 ! East Bab el Mandeb surface north point (T) 142 zsums = 0.e0 ! East Bab el Mandeb surface north point (S) 143 zsumt1 = 0.e0 ! East Bab el Mandeb depth south point (T) 144 zsums1 = 0.e0 ! East Bab el Mandeb depth south point (S) 145 zsumt2 = 0.e0 ! West Bab el Mandeb surface (T) 146 zsums2 = 0.e0 ! West Bab el Mandeb surface (S) 147 zsumt3 = 0.e0 ! West Bab el Mandeb depth (T) 148 zsums3 = 0.e0 ! West Bab el Mandeb depth (S) 149 zsumt4 = 0.e0 ! East Bab el Mandeb depth north point (T) 150 zsums4 = 0.e0 ! East Bab el Mandeb depth north point (S) 151 152 ! EMP of the Red Sea 153 ! ------------------ 154 155 zempred = 0.e0 156 zwei = 0.e0 157 DO jj = mj0(87), mj1(96) 158 DO ji = mi0(148), mi1(160) 159 zwei = tmask(ji,jj,1) * e1t(ji,jj) * e2t(ji,jj) 160 zempred = zempred + ( emp(ji,jj) - rnf(ji,jj) ) * zwei 161 END DO 162 END DO 163 IF( lk_mpp ) CALL mpp_sum( zempred ) ! sum with other processors value 164 165 ! convert in m3 166 zempred = zempred * 1.e-3 167 168 ! Velocity profile at each point 169 ! ------------------------------ 170 171 zu1_rs(:) = zu1_rs_i(:) 172 zu2_rs(:) = zu2_rs_i(:) 173 zu3_rs(:) = zu3_rs_i(:) 174 175 ! velocity profile at 161,88 East Bab el Mandeb North point 176 ! we imposed zisw_rs + EMP above the Red Sea 177 DO jk = 1, 8 178 DO jj = mj0(88), mj1(88) 179 DO ji = mi0(160), mi1(160) 180 zu1_rs(jk) = zu1_rs(jk) - ( zempred / 8. ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) 181 END DO 182 END DO 183 END DO 184 185 ! velocity profile at 161, 88 West Bab el Mandeb 186 ! we imposed zisw_rs + EMP above the Red Sea 187 DO jk = 1, 10 188 DO jj = mj0(88), mj1(88) 189 DO ji = mi0(160), mi1(160) 190 zu3_rs(jk) = zu3_rs(jk) + ( zempred / 10. ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) 191 END DO 192 END DO 193 END DO 194 195 ! Balance of temperature and salinity 196 ! ----------------------------------- 197 198 ! east Bab el Mandeb surface vertical sum of transport* S,T 199 DO jk = 1, 19 200 DO jj = mj0(88), mj1(88) 201 DO ji = mi0(161), mi1(161) 202 zsumt = zsumt + tn(ji,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * zu1_rs(jk) 203 zsums = zsums + sn(ji,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * zu1_rs(jk) 204 END DO 205 END DO 206 END DO 207 208 ! west Bab el Mandeb surface vertical sum of transport* S,T 209 DO jk = 1, 10 210 DO jj = mj0(88), mj1(88) 211 DO ji = mi0(161), mi1(161) 212 zsumt2 = zsumt2 + tn(ji,jj,jk) * e1v(ji-1,jj) * fse3v(ji-1,jj,jk) * zu3_rs(jk) 213 zsums2 = zsums2 + sn(ji,jj,jk) * e1v(ji-1,jj) * fse3v(ji-1,jj,jk) * zu3_rs(jk) 214 END DO 215 END DO 216 END DO 217 218 ! west Bab el Mandeb deeper 219 DO jj = mj0(89), mj1(89) 220 DO ji = mi0(160), mi1(160) 221 zsumt3 = tn(ji,jj,16) * e1v(ji,jj-1) * fse3v(ji,jj-1,16) * zu3_rs(16) 222 zsums3 = sn(ji,jj,16) * e1v(ji,jj-1) * fse3v(ji,jj-1,16) * zu3_rs(16) 223 END DO 224 END DO 225 226 ! east Bab el Mandeb deeper 227 DO jk = 20, 21 228 DO jj = mj0(88), mj1(88) 229 DO ji = mi0(161), mi1(161) 230 zsumt4 = zsumt4 + tn(ji,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * zu1_rs(jk) 231 zsums4 = zsums4 + sn(ji,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * zu1_rs(jk) 232 END DO 233 END DO 234 END DO 235 236 ! Total transport 237 zsumt1 = -( zsumt3 + zsumt2 + zsumt + zsumt4 ) 238 zsums1 = -( zsums3 + zsums2 + zsums + zsums4 ) 239 240 ! Temperature and Salinity at East Bab el Mandeb, Level 21 241 DO jj = mj0(88), mj1(88) 242 DO ji = mi0(160), mi1(160) 243 zt = zsumt1 / ( zu2_rs(21) * e2u(ji,jj-1) * fse3u(ji,jj-1,21) ) 244 zs = zsums1 / ( zu2_rs(21) * e2u(ji,jj-1) * fse3u(ji,jj-1,21) ) 245 END DO 246 END DO 247 248 ! New Temperature and Salinity at East Bab el Mandeb 249 ! -------------------------------------------------- 250 251 ! north point 252 DO jk = 1, jpk 253 DO jj = mj0(88), mj1(88) 254 DO ji = mi0(161), mi1(161) 255 zvt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 256 zsu = e2u(ji-1,jj) * fse3u(ji-1,jj,jk) 257 ta(ji,jj,jk) = ta(ji,jj,jk) + ( 1. / zvt ) * zsu * zu1_rs(jk) * tn(ji,jj,jk) 258 sa(ji,jj,jk) = sa(ji,jj,jk) + ( 1. / zvt ) * zsu * zu1_rs(jk) * sn(ji,jj,jk) 259 END DO 260 END DO 261 END DO 262 263 ! south point 264 jk = 21 265 DO jj = mj0(87), mj1(87) 266 DO ji = mi0(161), mi1(161) 267 zvt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 268 zsu = e2u(ji-1,jj) * fse3u(ji-1,jj,jk) 269 ta(ji,jj,jk) = ta(ji,jj,jk) + ( 1. / zvt ) * zsu * zu2_rs(jk) * zt 270 sa(ji,jj,jk) = sa(ji,jj,jk) + ( 1. / zvt ) * zsu * zu2_rs(jk) * zs 271 END DO 272 END DO 273 274 275 ! New Temperature and Salinity at West Bab el Mandeb 276 ! -------------------------------------------------- 277 278 ! surface 279 DO jk = 1, 10 280 DO jj = mj0(89), mj1(89) 281 DO ji = mi0(160), mi1(160) 282 zvt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 283 zsu = e1v(ji,jj-1) * fse3v(ji,jj-1,jk) 284 ta(ji,jj,jk) = ta(ji,jj,jk) + ( 1. / zvt ) * zsu * zu3_rs(jk) * tn(ji+1,jj-1,jk) 285 sa(ji,jj,jk) = sa(ji,jj,jk) + ( 1. / zvt ) * zsu * zu3_rs(jk) * sn(ji+1,jj-1,jk) 286 END DO 287 END DO 288 END DO 289 ! deeper 290 jk = 16 291 DO jj = mj0(89), mj1(89) 292 DO ji = mi0(160), mi1(160) 293 zvt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 294 zsu = e1v(ji,jj-1) * fse3v(ji,jj-1,jk) 295 ta(ji,jj,jk) = ta(ji,jj,jk) + ( 1. / zvt ) * zsu * zu3_rs(jk) * tn(ji,jj,jk) 296 sa(ji,jj,jk) = sa(ji,jj,jk) + ( 1. / zvt ) * zsu * zu3_rs(jk) * sn(ji,jj,jk) 297 END DO 298 END DO 299 300 END SUBROUTINE tra_bab_el_mandeb 301 302 303 SUBROUTINE tra_gibraltar 463 ! 464 SELECT CASE( cd_td ) 465 ! ! ---------------- ! 466 CASE( 'ini' ) ! initialisation ! 467 ! ! ---------------- ! 468 ! !** initialization of the velocity 469 hdiv_139_101(:) = 0.e0 ! 139,101 (Atlantic side, south point) 470 hdiv_139_102(:) = 0.e0 ! 139,102 (Atlantic side, north point) 471 hdiv_141_102(:) = 0.e0 ! 141,102 (Med sea side) 472 473 ! !** imposed transport 474 zio_flow = 0.8e6 ! inflow surface water 475 zrecirc_mid = 0.7e6 ! middle recirculation water 476 zrecirc_upp = 2.5e6 ! upper recirculation water 477 zrecirc_bot = 3.5e6 ! bottom recirculation water 478 ! 479 DO jj = mj0(101), mj1(101) !** profile of hdiv at 139,101 (Atlantic side, south point) 480 DO ji = mi0(139), mi1(139) !----------------------------- 481 DO jk = 1, 14 ! surface in/out flow (Atl -> Med) (div >0) 482 hdiv_139_101(jk) = + zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 483 END DO 484 DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div >0) 485 hdiv_139_101(jk) = + zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 486 END DO 487 ! ! upper reciculation (Atl 101 -> Atl 101) (div >0) 488 hdiv_139_101(21) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 489 ! 490 ! ! upper & bottom reciculation (Atl 101 -> Atl 101 & 102) (div >0) 491 hdiv_139_101(22) = ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 492 END DO 493 END DO 494 DO jj = mj0(102), mj1(102) !** profile of hdiv at 139,102 (Atlantic side, north point) 495 DO ji = mi0(139), mi1(139) !----------------------------- 496 DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0) 497 hdiv_139_102(jk) = - zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 498 END DO 499 ! ! outflow of Mediterranean sea + deep recirculation (div <0) 500 hdiv_139_102(22) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 501 END DO 502 END DO 503 DO jj = mj0(102), mj1(102) !** velocity profile at 141,102 (Med sea side) 504 DO ji = mi0(141), mi1(141) !------------------------------ 505 DO jk = 1, 14 ! surface inflow in the Med (div <0) 506 hdiv_141_102(jk) = - zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 507 END DO 508 ! ! deep outflow toward the Atlantic (div >0) 509 hdiv_141_102(21) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 510 END DO 511 END DO 512 ! ! ---------------- ! 513 CASE( 'div' ) ! update hdivn ! (call by divcur module) 514 ! ! ---------=====-- ! 515 ! !** emp on the Mediterranean Sea (div >0) 516 zemp_med = 0.e0 !------------------------------- 517 DO jj = mj0(96), mj1(110) ! sum over the Med sea 518 DO ji = mi0(141),mi1(181) 519 zemp_med = zemp_med + emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj) 520 END DO 521 END DO 522 DO jj = mj0(96), mj1(96) ! minus 2 points in Red Sea 523 DO ji = mi0(148),mi1(148) 524 zemp_med = zemp_med - emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj) 525 END DO 526 DO ji = mi0(149),mi1(149) 527 zemp_med = zemp_med - emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj) 528 END DO 529 END DO 530 IF( lk_mpp ) CALL mpp_sum( zemp_med ) ! sum with other processors value 531 zemp_med = zemp_med * 1.e-3 ! convert in m3 532 ! 533 ! !** Correct hdivn (including emp adjustment) 534 ! !------------------------------------------- 535 DO jj = mj0(101), mj1(101) !* 139,101 (Atlantic side, south point) 536 DO ji = mi0(139), mi1(139) 537 hdiv_139_101_kt(:) = hdiv_139_101(:) 538 DO jk = 1, 14 ! increase the inflow from the Atlantic (div >0) 539 hdiv_139_101_kt(jk) = hdiv_139_101(jk) + zemp_med / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 540 END DO 541 hdivn(ji, jj,:) = hdivn(ji, jj,:) + hdiv_139_101_kt(:) 542 END DO 543 END DO 544 DO jj = mj0(102), mj1(102) !* 139,102 (Atlantic side, north point) 545 DO ji = mi0(139), mi1(139) 546 hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_139_102(:) 547 END DO 548 END DO 549 DO jj = mj0(102), mj1(102) !* 141,102 (Med side) 550 DO ji = mi0(141), mi1(141) 551 hdiv_141_102(:) = hdiv_141_102(:) 552 DO jk = 1, 14 ! increase the inflow from the Atlantic (div <0) 553 hdiv_141_102_kt(jk) = hdiv_141_102(jk) - zemp_med / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 554 END DO 555 hdivn(ji, jj,:) = hdivn(ji, jj,:) + hdiv_141_102_kt(:) 556 END DO 557 END DO 558 ! ! ---------------- ! 559 CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module) 560 ! ! --------=======- ! 561 ! 562 DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point) (div >0) 563 DO ji = mi0(139), mi1(139) 564 DO jk = 1, jpkm1 ! surf inflow + mid. & bottom reciculation (from Atlantic) 565 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_139_101_kt(jk) * tn(ji,jj,jk) 566 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_139_101_kt(jk) * sn(ji,jj,jk) 567 END DO 568 END DO 569 END DO 570 ! 571 DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point) (div <0) 572 DO ji = mi0(139), mi1(139) 573 DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0) 574 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_139_102(jk) * tn(ji,jj-1,jk) ! middle Atlantic recirculation 575 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_139_102(jk) * sn(ji,jj-1,jk) 576 END DO 577 ! ! upper & bottom Atl. reciculation (Atl 101 -> Atl 102) - (div <0) 578 ! ! deep Med flow (Med 102 -> Atl 102) - (div <0) 579 ta(ji,jj,22) = ta(ji,jj,22) + hdiv_141_102(21) * tn(ji+2,jj ,21) & ! deep Med flow 580 & + hdiv_139_101(21) * tn(ji ,jj-1,21) & ! upper Atlantic recirculation 581 & + hdiv_139_101(22) * tn(ji ,jj-1,22) ! bottom Atlantic recirculation 582 sa(ji,jj,22) = sa(ji,jj,22) + hdiv_141_102(21) * sn(ji+2,jj ,21) & 583 & + hdiv_139_101(21) * sn(ji ,jj-1,21) & 584 & + hdiv_139_101(22) * sn(ji ,jj-1,22) 585 END DO 586 END DO 587 DO jj = mj0(102), mj1(102) !* 141,102 (Med side) (div <0) 588 DO ji = mi0(141), mi1(141) 589 DO jk = 1, 14 ! surface flow from Atlantic to Med sea 590 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_141_102_kt(jk) * tn(ji-2,jj-1,jk) 591 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_141_102_kt(jk) * sn(ji-2,jj-1,jk) 592 END DO 593 ! ! deeper flow from Med sea to Atlantic 594 ta(ji,jj,21) = ta(ji,jj,21) - hdiv_141_102(21) * tn(ji,jj,21) 595 sa(ji,jj,21) = sa(ji,jj,21) - hdiv_141_102(21) * sn(ji,jj,21) 596 END DO 597 END DO 598 ! ! ---------------- ! 599 CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module) 600 ! ! --------=======- ! 601 ! at this stage, (ua,va) are the after velocity, not the tendancy 602 ! compute the velocity from the divergence at T-point 603 ! 604 DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point) 605 DO ji = mi0(139), mi1(139) ! div >0 => ua >0, same sign 606 ua(ji,jj,:) = hdiv_139_101_kt(:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) & 607 & * e2u(ji,jj) * fse3u(ji,jj,:) 608 END DO 609 END DO 610 DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point) 611 DO ji = mi0(139), mi1(139) ! div <0 => ua <0, same sign 612 ua(ji,jj,:) = hdiv_139_102(:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) & 613 & * e2u(ji,jj) * fse3u(ji,jj,:) 614 END DO 615 END DO 616 DO jj = mj0(102), mj1(102) !** 140,102 (Med side) (140 not 141 as it is a U-point) 617 DO ji = mi0(140), mi1(140) ! div >0 => ua <0, opposite sign 618 ua(ji,jj,:) = - hdiv_141_102(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) & 619 & * e2u(ji,jj) * fse3u(ji,jj,:) 620 END DO 621 END DO 622 ! 623 END SELECT 624 ! 625 END SUBROUTINE cla_gibraltar 626 627 628 SUBROUTINE cla_hormuz( cd_td ) 629 !! ------------------------------------------------------------------- 630 !! *** ROUTINE div_hormuz *** 631 !! 632 !! ** Purpose : update the now horizontal divergence, the tracer 633 !! tendancyand the after velocity in vicinity of Hormuz 634 !! strait ( Persian Gulf - Indian ocean ). 635 !! 636 !! ** Method : Hormuz strait 637 !! ______________ 638 !! |/////|<== surface inflow 639 !! 94 |/////| 640 !! |/////|==> deep outflow 641 !! |_____|_______ 642 !! 171 172 304 643 !!--------------------------------------------------------------------- 305 !! *** ROUTINE tra_gibraltar *** 306 !! 307 !! ** Purpose : 308 !! Update the horizontal advective trend of tracers (t & s) 309 !! correction in Gibraltar and 310 !! add it to the general trend of tracer equations. 311 !! 312 !! ** Method : 313 !! We impose transport at Gibraltar and knowing T and S in 314 !! surface and deeper at each side of the strait, we deduce T and S 315 !! of the outflow of the Mediterranean Sea in the Atlantic ocean . 316 !! 317 !! ________________ N ________________ 318 !! 102 | |-> | <-| |<- 319 !! 101 ___->|____|_____ W - - E ___->|____|_____ 320 !! 139 140 141 | 139 140 141 321 !! horizontal view S horizontal view 322 !! surface depth 323 !! C A U T I O N : the trend saved is the centered trend only. 324 !! It doesn't take into account the upstream part of the scheme. 325 !! 326 !! ** history : 327 !! ! 02-06 (A. Bozec) Original code 328 !! 8.5 ! 02-11 (A. Bozec) F90: Free form and module 644 CHARACTER(len=1), INTENT(in) :: cd_td ! ='ini' initialisation 645 !! ! ='div' update the divergence 646 !! ! ='tra' update the tracers 647 !! ! ='spg' update after velocity 648 !! 649 INTEGER :: ji, jj, jk ! dummy loop indices 650 REAL(wp) :: zio_flow ! temporary scalar 329 651 !!--------------------------------------------------------------------- 330 !! * Local declarations 331 INTEGER :: ji, jj, jk ! dummy loop indices 332 REAL(wp) :: zsu, zvt 333 REAL(wp) :: zsumt, zsumt1, zsumt2, zsumt3, zsumt4 334 REAL(wp) :: zsums, zsums1, zsums2, zsums3, zsums4 335 REAL(wp) :: zt, zs 336 REAL(wp) :: zwei 337 REAL(wp), DIMENSION (jpk) :: zu1_ms, zu2_ms, zu3_ms 338 !!--------------------------------------------------------------------- 339 340 ! Initialization of vertical sum for T and S transport 341 ! ---------------------------------------------------- 342 343 zsumt = 0.e0 ! West Gib. surface south point ( T ) 344 zsums = 0.e0 ! West Gib. surface south point ( S ) 345 zsumt1 = 0.e0 ! East Gib. surface north point ( T ) 346 zsums1 = 0.e0 ! East Gib. surface north point ( S ) 347 zsumt2 = 0.e0 ! East Gib. depth north point ( T ) 348 zsums2 = 0.e0 ! East Gib. depth north point ( S ) 349 zsumt3 = 0.e0 ! West Gib. depth south point ( T ) 350 zsums3 = 0.e0 ! West Gib. depth south point ( S ) 351 zsumt4 = 0.e0 ! West Gib. depth north point ( T ) 352 zsums4 = 0.e0 ! West Gib. depth north point ( S ) 353 354 ! EMP of Mediterranean Sea 355 ! ------------------------ 356 357 zempmed = 0.e0 358 zwei = 0.e0 359 DO jj = mj0(96),mj1(110) 360 DO ji = mi0(141),mi1(181) 361 zwei = tmask(ji,jj,1) * e1t(ji,jj) * e2t(ji,jj) 362 zempmed = zempmed + ( emp(ji,jj) - rnf(ji,jj) ) * zwei 363 END DO 364 END DO 365 IF( lk_mpp ) CALL mpp_sum( zempmed ) ! sum with other processors value 366 367 368 ! minus 2 points in Red Sea and 3 in Atlantic ocean 369 DO jj = mj0(96),mj1(96) 370 DO ji = mi0(148),mi1(148) 371 zempmed = zempmed - ( emp(ji ,jj)-rnf(ji ,jj) ) * tmask(ji ,jj,1) * e1t(ji ,jj) * e2t(ji ,jj) & 372 - ( emp(ji+1,jj)-rnf(ji+1,jj) ) * tmask(ji+1,jj,1) * e1t(ji+1,jj) * e2t(ji+1,jj) 373 END DO 374 END DO 375 376 ! convert in m3 377 zempmed = zempmed * 1.e-3 378 379 ! Velocity profile at each point 380 ! ------------------------------ 381 382 zu1_ms(:) = zu1_ms_i(:) 383 zu2_ms(:) = zu2_ms_i(:) 384 zu3_ms(:) = zu3_ms_i(:) 385 386 ! velocity profile at 139,101 South point + (emp-rnf) on surface 387 DO jk = 1, 14 388 DO jj = mj0(102), mj1(102) 389 DO ji = mi0(140), mi1(140) 390 zu1_ms(jk) = zu1_ms(jk) + ( zempmed / 14. ) / ( e2u(ji-1, jj-1) * fse3u(ji-1, jj-1,jk) ) 391 END DO 392 END DO 393 END DO 394 395 ! profile at East Gibraltar 396 ! velocity profile at 141,102 + (emp-rnf) on surface 397 DO jk = 1, 14 398 DO jj = mj0(102), mj1(102) 399 DO ji = mi0(140), mi1(140) 400 zu3_ms(jk) = zu3_ms(jk) + ( zempmed / 14. ) / ( e2u(ji, jj) * fse3u(ji, jj,jk) ) 401 END DO 402 END DO 403 END DO 404 405 ! Balance of temperature and salinity 406 ! ----------------------------------- 407 408 ! west gibraltar surface vertical sum of transport* S,T 409 DO jk = 1, 14 410 DO jj = mj0(101), mj1(101) 411 DO ji = mi0(139), mi1(139) 412 zsumt = zsumt + tn(ji, jj,jk) * e2u(ji, jj) * fse3u(ji, jj,jk) * zu1_ms(jk) 413 zsums = zsums + sn(ji, jj,jk) * e2u(ji, jj) * fse3u(ji, jj,jk) * zu1_ms(jk) 414 END DO 415 END DO 416 END DO 417 418 ! east Gibraltar surface vertical sum of transport* S,T 419 DO jk = 1, 14 420 DO jj = mj0(101), mj1(101) 421 DO ji = mi0(139), mi1(139) 422 zsumt1 = zsumt1 + tn(ji, jj,jk) * e2u(ji+1, jj+1) * fse3u(ji+1, jj+1,jk) * zu3_ms(jk) 423 zsums1 = zsums1 + sn(ji, jj,jk) * e2u(ji+1, jj+1) * fse3u(ji+1, jj+1,jk) * zu3_ms(jk) 424 END DO 425 END DO 426 END DO 427 428 ! east Gibraltar deeper vertical sum of transport* S,T 429 DO jj = mj0(102), mj1(102) 430 DO ji = mi0(141), mi1(141) 431 zsumt2 = tn(ji, jj,21) * e2u(ji-1, jj) * fse3u(ji-1, jj,21) * zu3_ms(21) 432 zsums2 = sn(ji, jj,21) * e2u(ji-1, jj) * fse3u(ji-1, jj,21) * zu3_ms(21) 433 END DO 434 END DO 435 436 ! west Gibraltar deeper vertical sum of transport* S,T 437 DO jk = 21, 22 438 DO jj = mj0(101), mj1(101) 439 DO ji = mi0(139), mi1(139) 440 zsumt3 = zsumt3 + tn(ji, jj,jk) * e2u(ji, jj) * fse3u(ji, jj,jk) * zu1_ms(jk) 441 zsums3 = zsums3 + sn(ji, jj,jk) * e2u(ji, jj) * fse3u(ji, jj,jk) * zu1_ms(jk) 442 END DO 443 END DO 444 END DO 445 446 ! Total transport = 0. 447 zsumt4 = zsumt2 + zsumt1 - zsumt - zsumt3 448 zsums4 = zsums2 + zsums1 - zsums - zsums3 449 450 ! Temperature and Salinity at West gibraltar , Level 22 451 DO jj = mj0(102), mj1(102) 452 DO ji = mi0(140), mi1(140) 453 zt = zsumt4 / ( zu2_ms(22) * e2u(ji-1, jj) * fse3u(ji-1, jj, 22) ) 454 zs = zsums4 / ( zu2_ms(22) * e2u(ji-1, jj) * fse3u(ji-1, jj, 22) ) 455 END DO 456 END DO 457 458 ! New Temperature and Salinity trend at West Gibraltar 459 ! ---------------------------------------------------- 460 461 ! south point 462 DO jk = 1, 22 463 DO jj = mj0(101), mj1(101) 464 DO ji = mi0(139), mi1(139) 465 zvt = e1t(ji, jj) * e2t(ji, jj) * fse3t(ji, jj,jk) 466 zsu = e2u(ji, jj) * fse3u(ji, jj,jk) 467 ta(ji, jj,jk) = ta(ji, jj,jk) - ( 1. / zvt ) * zsu * zu1_ms(jk) * tn(ji, jj,jk) 468 sa(ji, jj,jk) = sa(ji, jj,jk) - ( 1. / zvt ) * zsu * zu1_ms(jk) * sn(ji, jj,jk) 469 END DO 470 END DO 471 END DO 472 473 ! north point 474 DO jk = 15, 20 475 DO jj = mj0(102), mj1(102) 476 DO ji = mi0(139), mi1(139) 477 zvt = e1t(ji, jj) * e2t(ji, jj) * fse3t(ji, jj,jk) 478 zsu = e2u(ji, jj) * fse3u(ji, jj,jk) 479 ta(ji, jj,jk) = ta(ji, jj,jk) - ( 1. / zvt ) * zsu * zu2_ms(jk) * tn(ji, jj-1,jk) 480 sa(ji, jj,jk) = sa(ji, jj,jk) - ( 1. / zvt ) * zsu * zu2_ms(jk) * sn(ji, jj-1,jk) 481 END DO 482 END DO 483 END DO 484 485 ! Gibraltar outflow, north point deeper 486 jk = 22 487 DO jj = mj0(102), mj1(102) 488 DO ji = mi0(139), mi1(139) 489 zvt = e1t(ji, jj) * e2t(ji, jj) * fse3t(ji, jj,jk) 490 zsu = e2u(ji, jj) * fse3u(ji, jj,jk) 491 ta(ji, jj,jk) = ta(ji, jj,jk) - ( 1. / zvt ) * zsu * zu2_ms(jk) * zt 492 sa(ji, jj,jk) = sa(ji, jj,jk) - ( 1. / zvt ) * zsu * zu2_ms(jk) * zs 493 END DO 494 END DO 495 496 497 ! New Temperature and Salinity at East Gibraltar 498 ! ---------------------------------------------- 499 500 ! surface 501 DO jk = 1, 14 502 DO jj = mj0(102), mj1(102) 503 DO ji = mi0(141), mi1(141) 504 zvt = e1t(ji, jj) * e2t(ji, jj) * fse3t(ji, jj,jk) 505 zsu = e2u(ji-1, jj) * fse3u(ji-1, jj,jk) 506 ta(ji, jj,jk) = ta(ji, jj,jk) + ( 1. / zvt ) * zsu * zu3_ms(jk) * tn(ji-2, jj-1,jk) 507 sa(ji, jj,jk) = sa(ji, jj,jk) + ( 1. / zvt ) * zsu * zu3_ms(jk) * sn(ji-2, jj-1,jk) 508 END DO 509 END DO 510 END DO 511 ! deeper 512 jk = 21 513 DO jj = mj0(102), mj1(102) 514 DO ji = mi0(141), mi1(141) 515 zvt = e1t(ji, jj) * e2t(ji, jj) * fse3t(ji, jj,jk) 516 zsu = e2u(ji-1, jj) * fse3u(ji-1, jj,jk) 517 ta(ji, jj,jk) = ta(ji, jj,jk) + ( 1. / zvt ) * zsu * zu3_ms(jk) * tn(ji, jj,jk) 518 sa(ji, jj,jk) = sa(ji, jj,jk) + ( 1. / zvt ) * zsu * zu3_ms(jk) * sn(ji, jj,jk) 519 END DO 520 END DO 521 522 END SUBROUTINE tra_gibraltar 523 524 525 SUBROUTINE tra_hormuz 526 !!--------------------------------------------------------------------- 527 !! *** ROUTINE tra_hormuz *** 528 !! 529 !! ** Purpose : Update the horizontal advective trend of tracers 530 !! correction in Hormuz. 531 !! 532 !! ** Method : We impose transport at Hormuz . 533 !! 534 !! ** history : 535 !! ! 02-11 (A. Bozec) Original code 536 !! 8.5 ! 02-11 (A. Bozec) F90: Free form and module 537 !!--------------------------------------------------------------------- 538 !! * Local declarations 539 INTEGER :: ji, jj, jk ! dummy loop indices 540 REAL(wp) :: zsu, zvt 541 !!--------------------------------------------------------------------- 542 543 ! New trend at Hormuz strait 544 ! -------------------------- 545 DO jk = 1, 8 546 DO jj = mj0(94), mj1(94) 652 ! 653 SELECT CASE( cd_td ) 654 ! ! ---------------- ! 655 CASE( 'ini' ) ! initialisation ! 656 ! ! ---------------- ! 657 ! !** profile of horizontal divergence due to cross-land advection 658 zio_flow = 1.e6 ! imposed in/out flow 659 ! 660 hdiv_172_94(:) = 0.e0 661 ! 662 DO jj = mj0(94), mj1(94) ! in/out flow at (i,j) = (172,94) 547 663 DO ji = mi0(172), mi1(172) 548 zvt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 549 zsu = e2u(ji-1,jj) * fse3u(ji-1,jj,jk) 550 ta(ji,jj,jk) = ta(ji,jj,jk) + ( 1. / zvt ) * zsu * zu_pg(jk) * tn(ji,jj,jk) 551 sa(ji,jj,jk) = sa(ji,jj,jk) + ( 1. / zvt ) * zsu * zu_pg(jk) * sn(ji,jj,jk) 552 END DO 553 END DO 554 END DO 555 DO jk = 16, 18 556 DO jj = mj0(94), mj1(94) 664 DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0) 665 hdiv_172_94(jk) = - ( zio_flow / 8.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 666 END DO 667 DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0) 668 hdiv_172_94(jk) = + ( zio_flow / 3.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 669 END DO 670 END DO 671 END DO 672 ! !** T & S profile in the Hormuz strait (use in deep outflow) 673 ! Temperature and Salinity 674 t_171_94_hor(:) = 0.e0 ; s_171_94_hor(:) = 0.e0 675 t_171_94_hor(16) = 18.4 ; s_171_94_hor(16) = 36.27 676 t_171_94_hor(17) = 17.8 ; s_171_94_hor(17) = 36.4 677 t_171_94_hor(18) = 16. ; s_171_94_hor(18) = 36.27 678 ! 679 ! ! ---------------- ! 680 CASE( 'div' ) ! update hdivn ! (call by divcur module) 681 ! ! ---------=====-- ! 682 ! 683 DO jj = mj0(94), mj1(94) !** 172,94 (Indian ocean side) 557 684 DO ji = mi0(172), mi1(172) 558 zvt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 559 zsu = e2u(ji-1,jj) * fse3u(ji-1,jj,jk) 560 ta(ji,jj,jk) = ta(ji,jj,jk) + ( 1. / zvt ) * zsu * zu_pg(jk) * zthor(jk) 561 sa(ji,jj,jk) = sa(ji,jj,jk) + ( 1. / zvt ) * zsu * zu_pg(jk) * zshor(jk) 562 END DO 563 END DO 564 END DO 565 566 END SUBROUTINE tra_hormuz 567 568 569 SUBROUTINE tra_cla_init 570 !!--------------------------------------------------------------------- 571 !! *** ROUTINE tra_cla_init *** 572 !! 573 !! ** Purpose : Initialization of variables 574 !! 575 !! ** history : 576 !! 9.0 ! 02-11 (A. Bozec) Original code 577 !!--------------------------------------------------------------------- 578 !! * Local declarations 579 INTEGER :: ji, jj, jk ! dummy loop indices 580 !!--------------------------------------------------------------------- 581 582 ! Control print 583 ! ------------- 584 585 IF(lwp) WRITE(numout,*) 586 IF(lwp) WRITE(numout,*) 'tra_cla_init : cross land advection on tracer ' 587 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 588 589 ! Initialization at Bab el Mandeb 590 ! ------------------------------- 591 592 ! imposed transport 593 zisw_rs = 0.4e6 ! inflow surface water 594 zurw_rs = 0.2e6 ! upper recirculation water 595 !!Alex zbrw_rs = 1.2e6 ! bottom recirculation water 596 zbrw_rs = 0.5e6 ! bottom recirculation water 597 598 ! initialization of the velocity at Bab el Mandeb 599 zu1_rs_i(:) = 0.e0 ! velocity profile at 161,88 South point 600 zu2_rs_i(:) = 0.e0 ! velocity profile at 161,87 North point 601 zu3_rs_i(:) = 0.e0 ! velocity profile at 160,88 East point 602 603 ! velocity profile at 161,88 East Bab el Mandeb North point 604 ! we imposed zisw_rs + EMP above the Red Sea 605 DO jk = 1, 8 606 DO jj = mj0(88), mj1(88) 607 DO ji = mi0(160), mi1(160) 608 zu1_rs_i(jk) = -( zisw_rs / 8. ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) 609 END DO 610 END DO 611 END DO 612 613 ! recirculation water 614 DO jj = mj0(88), mj1(88) 615 DO ji = mi0(160), mi1(160) 616 zu1_rs_i(20) = -( zurw_rs ) / ( e2u(ji,jj) * fse3u(ji,jj,20) ) 617 zu1_rs_i(21) = -( zbrw_rs - zurw_rs ) / ( e2u(ji,jj) * fse3u(ji,jj,21) ) 618 END DO 619 END DO 620 621 ! velocity profile at 161,87 East Bab el Mandeb South point 622 DO jj = mj0(87), mj1(87) 623 DO ji = mi0(160), mi1(160) 624 zu2_rs_i(21) = ( zbrw_rs + zisw_rs ) / ( e2u(ji,jj) * fse3u(ji,jj,21) ) 625 END DO 626 END DO 627 628 ! velocity profile at 161, 88 West Bab el Mandeb 629 ! we imposed zisw_rs + EMP above the Red Sea 630 DO jk = 1, 10 631 DO jj = mj0(88), mj1(88) 632 DO ji = mi0(160), mi1(160) 633 zu3_rs_i(jk) = ( zisw_rs / 10. ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) 634 END DO 635 END DO 636 END DO 637 638 ! deeper 639 DO jj = mj0(88), mj1(88) 640 DO ji = mi0(160), mi1(160) 641 zu3_rs_i(16) = - zisw_rs /( e1v(ji,jj) * fse3v(ji,jj,16) ) 642 END DO 643 END DO 644 645 646 ! Initialization at Gibraltar 647 ! --------------------------- 648 649 ! imposed transport 650 zisw_ms = 0.8e6 ! atlantic-mediterranean water 651 zmrw_ms = 0.7e6 ! middle recirculation water 652 zurw_ms = 2.5e6 ! upper recirculation water 653 zbrw_ms = 3.5e6 ! bottom recirculation water 654 655 ! initialization of the velocity 656 zu1_ms_i(:) = 0.e0 ! velocity profile at 139,101 South point 657 zu2_ms_i(:) = 0.e0 ! velocity profile at 139,102 North point 658 zu3_ms_i(:) = 0.e0 ! velocity profile at 141,102 East point 659 660 ! velocity profile at 139,101 South point 661 DO jk = 1, 14 662 DO jj = mj0(102), mj1(102) 663 DO ji = mi0(140), mi1(140) 664 zu1_ms_i(jk) = ( zisw_ms / 14. ) / ( e2u(ji-1, jj-1) * fse3u(ji-1, jj-1,jk)) 665 END DO 666 END DO 667 END DO 668 669 ! middle recirculation ( uncounting in the balance ) 670 DO jk = 15, 20 671 DO jj = mj0(102), mj1(102) 672 DO ji = mi0(140), mi1(140) 673 zu1_ms_i(jk) = ( zmrw_ms / 6. ) / ( e2u(ji-1, jj-1) * fse3u(ji-1, jj-1,jk) ) 674 END DO 675 END DO 676 END DO 677 678 DO jj = mj0(102), mj1(102) 679 DO ji = mi0(140), mi1(140) 680 zu1_ms_i(21) = ( zurw_ms ) / ( e2u(ji-1, jj-1) * fse3u(ji-1, jj-1,21) ) 681 zu1_ms_i(22) = ( zbrw_ms - zurw_ms ) / ( e2u(ji-1, jj-1) * fse3u(ji-1, jj-1,22) ) 682 END DO 683 END DO 684 685 ! velocity profile at 139,102 North point 686 ! middle recirculation ( uncounting in the balance ) 687 DO jk = 15, 20 688 DO jj = mj0(102), mj1(102) 689 DO ji = mi0(140), mi1(140) 690 zu2_ms_i(jk) = -( zmrw_ms / 6. ) / ( e2u(ji-1, jj) * fse3u(ji-1, jj,jk) ) 691 END DO 692 END DO 693 END DO 694 695 DO jj = mj0(102), mj1(102) 696 DO ji = mi0(140), mi1(140) 697 zu2_ms_i(22) = -( zisw_ms + zbrw_ms ) / ( e2u(ji-1, jj) * fse3u(ji-1, jj,22) ) 698 END DO 699 END DO 700 701 ! profile at East Gibraltar 702 ! velocity profile at 141,102 703 DO jk = 1, 14 704 DO jj = mj0(102), mj1(102) 705 DO ji = mi0(140), mi1(140) 706 zu3_ms_i(jk) = ( zisw_ms / 14. ) / ( e2u(ji, jj) * fse3u(ji, jj,jk) ) 707 END DO 708 END DO 709 END DO 710 711 ! deeper 712 DO jj = mj0(102), mj1(102) 713 DO ji = mi0(140), mi1(140) 714 zu3_ms_i(21) = -zisw_ms / ( e2u(ji, jj) * fse3u(ji, jj,21) ) 715 END DO 716 END DO 717 718 719 ! Initialization at Hormuz 720 ! ------------------------ 721 722 ! imposed transport 723 zisw_pg = 4. * 0.25e6 ! surface and bottom water 724 725 ! initialization of the velocity 726 zu_pg(:) = 0.e0 ! velocity profile at 139,101 South point 727 728 ! Velocity profile 729 DO jk = 1, 8 730 DO jj = mj0(94), mj1(94) 685 hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_172_94(:) 686 END DO 687 END DO 688 ! ! ---------------- ! 689 CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module) 690 ! ! --------=======- ! 691 ! 692 DO jj = mj0(94), mj1(94) !** 172,94 (Indian ocean side) 731 693 DO ji = mi0(172), mi1(172) 732 zu_pg(jk) = -( zisw_pg / 8. ) / ( e2u(ji-1,jj) * fse3u(ji-1,jj,jk))733 END DO734 END DO735 END DO736 DO jk = 16, 18737 DO jj = mj0(94), mj1(94)738 DO ji = mi0(172), mi1(172)739 zu_pg(jk) = ( zisw_pg / 3. ) / ( e2u(ji-1,jj) * fse3u(ji-1,jj,jk) )740 END DO 741 END DO 742 END DO743 744 ! Temperature and Salinity at Hormuz745 zthor(:) = 0.e0746 zshor(:) = 0.e0747 748 zthor(16) = 18.4749 zshor(16) = 36.27750 !751 zthor(17) = 17.8752 zshor(17) = 36.4753 !754 zthor(18) = 16.755 zshor(18) = 36.27756 757 END SUBROUTINE tra_cla_init758 694 DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0) 695 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_172_94(jk) * tn(ji,jj,jk) 696 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_172_94(jk) * sn(ji,jj,jk) 697 END DO 698 DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0) 699 ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_172_94(jk) * t_171_94_hor(jk) 700 sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_172_94(jk) * s_171_94_hor(jk) 701 END DO 702 END DO 703 END DO 704 ! ! ---------------- ! 705 CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module) 706 ! ! --------=======- ! 707 ! No barotropic flow through Hormuz strait 708 ! at this stage, (ua,va) are the after velocity, not the tendancy 709 ! compute the velocity from the divergence at T-point 710 DO jj = mj0(94), mj1(94) !** 171,94 (Indian ocean side) (171 not 172 as it is the western U-point) 711 DO ji = mi0(171), mi1(171) ! div >0 => ua >0, opposite sign 712 ua(ji,jj,:) = - hdiv_172_94(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) & 713 & * e2u(ji,jj) * fse3u(ji,jj,:) 714 END DO 715 END DO 716 ! 717 END SELECT 718 ! 719 END SUBROUTINE cla_hormuz 720 759 721 #else 760 722 !!---------------------------------------------------------------------- 761 !! Default option NO cross land advection723 !! Default key Dummy module 762 724 !!---------------------------------------------------------------------- 763 USE in_out_manager 725 USE in_out_manager ! I/O manager 764 726 CONTAINS 765 SUBROUTINE tra_cla_init 766 END SUBROUTINE tra_cla_init 767 SUBROUTINE tra_cla( kt ) 768 INTEGER, INTENT(in) :: kt ! ocean time-step indice 769 IF( kt == nit000 .AND. lwp ) THEN 770 WRITE(numout,*) 771 WRITE(numout,*) 'tra_cla : No use of cross land advection' 772 WRITE(numout,*) '~~~~~~~' 773 ENDIF 774 END SUBROUTINE tra_cla 727 SUBROUTINE cla_init 728 CALL ctl_stop( 'cla_init: Cross Land Advection hard coded for ORCA_R2 with 31 levels' ) 729 END SUBROUTINE cla_init 730 SUBROUTINE cla_div( kt ) 731 WRITE(*,*) 'cla_div: You should have not see this print! error?', kt 732 END SUBROUTINE cla_div 733 SUBROUTINE cla_traadv( kt ) 734 WRITE(*,*) 'cla_traadv: You should have not see this print! error?', kt 735 END SUBROUTINE cla_traadv 736 SUBROUTINE cla_dynspg( kt ) 737 WRITE(*,*) 'dyn_spg_cla: You should have not see this print! error?', kt 738 END SUBROUTINE cla_dynspg 775 739 #endif 776 740 777 741 !!====================================================================== 778 742 END MODULE cla
Note: See TracChangeset
for help on using the changeset viewer.