- Timestamp:
- 2020-09-29T12:41:06+02:00 (4 years ago)
- Location:
- NEMO/branches/2020/r12377_ticket2386
- Files:
-
- 2 edited
-
. (modified) (1 prop)
-
src/OCE/DYN/dynzdf.F90 (modified) (20 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/r12377_ticket2386
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev @HEADext/AGRIF5 ^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 8 9 9 # SETTE 10 ^/utils/CI/sette@ HEADsette10 ^/utils/CI/sette@13507 sette
-
- Property svn:externals
-
NEMO/branches/2020/r12377_ticket2386/src/OCE/DYN/dynzdf.F90
r12511 r13540 38 38 !! * Substitutions 39 39 # include "do_loop_substitute.h90" 40 # include "domzgr_substitute.h90" 40 41 !!---------------------------------------------------------------------- 41 42 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 55 56 !! ** Method : - Leap-Frog time stepping on all trends but the vertical mixing 56 57 !! u(after) = u(before) + 2*dt * u(rhs) vector form or linear free surf. 57 !! u(after) = ( e3u_b*u(before) + 2*dt * e3u_n*u(rhs) ) / e3u (after)otherwise58 !! u(after) = ( e3u_b*u(before) + 2*dt * e3u_n*u(rhs) ) / e3u_after otherwise 58 59 !! - update the after velocity with the implicit vertical mixing. 59 60 !! This requires to solver the following system: 60 !! u(after) = u(after) + 1/e3u (after) dk+1[ mi(avm) / e3uw(after)dk[ua] ]61 !! u(after) = u(after) + 1/e3u_after dk+1[ mi(avm) / e3uw_after dk[ua] ] 61 62 !! with the following surface/top/bottom boundary condition: 62 63 !! surface: wind stress input (averaged over kt-1/2 & kt+1/2) … … 106 107 ! ! time stepping except vertical diffusion 107 108 IF( ln_dynadv_vec .OR. ln_linssh ) THEN ! applied on velocity 108 DO jk = 1, jpkm1109 puu( :,:,jk,Kaa) = ( puu(:,:,jk,Kbb) + rDt * puu(:,:,jk,Krhs) ) * umask(:,:,jk)110 pvv( :,:,jk,Kaa) = ( pvv(:,:,jk,Kbb) + rDt * pvv(:,:,jk,Krhs) ) * vmask(:,:,jk)111 END DO109 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 110 puu(ji,jj,jk,Kaa) = ( puu(ji,jj,jk,Kbb) + rDt * puu(ji,jj,jk,Krhs) ) * umask(ji,jj,jk) 111 pvv(ji,jj,jk,Kaa) = ( pvv(ji,jj,jk,Kbb) + rDt * pvv(ji,jj,jk,Krhs) ) * vmask(ji,jj,jk) 112 END_3D 112 113 ELSE ! applied on thickness weighted velocity 113 DO jk = 1, jpkm1 114 puu(:,:,jk,Kaa) = ( e3u(:,:,jk,Kbb) * puu(:,:,jk,Kbb) & 115 & + rDt * e3u(:,:,jk,Kmm) * puu(:,:,jk,Krhs) ) / e3u(:,:,jk,Kaa) * umask(:,:,jk) 116 pvv(:,:,jk,Kaa) = ( e3v(:,:,jk,Kbb) * pvv(:,:,jk,Kbb) & 117 & + rDt * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Krhs) ) / e3v(:,:,jk,Kaa) * vmask(:,:,jk) 118 END DO 114 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 115 puu(ji,jj,jk,Kaa) = ( e3u(ji,jj,jk,Kbb) * puu(ji,jj,jk,Kbb ) & 116 & + rDt * e3u(ji,jj,jk,Kmm) * puu(ji,jj,jk,Krhs) ) & 117 & / e3u(ji,jj,jk,Kaa) * umask(ji,jj,jk) 118 pvv(ji,jj,jk,Kaa) = ( e3v(ji,jj,jk,Kbb) * pvv(ji,jj,jk,Kbb ) & 119 & + rDt * e3v(ji,jj,jk,Kmm) * pvv(ji,jj,jk,Krhs) ) & 120 & / e3v(ji,jj,jk,Kaa) * vmask(ji,jj,jk) 121 END_3D 119 122 ENDIF 120 123 ! ! add top/bottom friction … … 124 127 ! G. Madec : in linear free surface, e3u(:,:,:,Kaa) = e3u(:,:,:,Kmm) = e3u_0, so systematic use of e3u(:,:,:,Kaa) 125 128 IF( ln_drgimp .AND. ln_dynspg_ts ) THEN 126 DO jk = 1, jpkm1! remove barotropic velocities127 puu( :,:,jk,Kaa) = ( puu(:,:,jk,Kaa) - uu_b(:,:,Kaa) ) * umask(:,:,jk)128 pvv( :,:,jk,Kaa) = ( pvv(:,:,jk,Kaa) - vv_b(:,:,Kaa) ) * vmask(:,:,jk)129 END DO130 DO_2D _00_00129 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! remove barotropic velocities 130 puu(ji,jj,jk,Kaa) = ( puu(ji,jj,jk,Kaa) - uu_b(ji,jj,Kaa) ) * umask(ji,jj,jk) 131 pvv(ji,jj,jk,Kaa) = ( pvv(ji,jj,jk,Kaa) - vv_b(ji,jj,Kaa) ) * vmask(ji,jj,jk) 132 END_3D 133 DO_2D( 0, 0, 0, 0 ) ! Add bottom/top stress due to barotropic component only 131 134 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 132 135 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 133 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 134 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 136 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) & 137 & + r_vvl * e3u(ji,jj,iku,Kaa) 138 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) & 139 & + r_vvl * e3v(ji,jj,ikv,Kaa) 135 140 puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 136 141 pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va 137 142 END_2D 138 IF( ln_isfcav ) THEN ! Ocean cavities (ISF)139 DO_2D _00_00143 IF( ln_isfcav.OR.ln_drgice_imp ) THEN ! Ocean cavities (ISF) 144 DO_2D( 0, 0, 0, 0 ) 140 145 iku = miku(ji,jj) ! top ocean level at u- and v-points 141 146 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 142 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) 143 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) 147 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) & 148 & + r_vvl * e3u(ji,jj,iku,Kaa) 149 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) & 150 & + r_vvl * e3v(ji,jj,ikv,Kaa) 144 151 puu(ji,jj,iku,Kaa) = puu(ji,jj,iku,Kaa) + rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * uu_b(ji,jj,Kaa) / ze3ua 145 152 pvv(ji,jj,ikv,Kaa) = pvv(ji,jj,ikv,Kaa) + rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * vv_b(ji,jj,Kaa) / ze3va … … 155 162 SELECT CASE( nldf_dyn ) 156 163 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzu) 157 DO_3D_00_00( 1, jpkm1 ) 158 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 164 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 165 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) & 166 & + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 159 167 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) + akzu(ji,jj,jk ) ) & 160 168 & / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) … … 168 176 END_3D 169 177 CASE DEFAULT ! iso-level lateral mixing 170 DO_3D_00_00( 1, jpkm1 ) 171 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 172 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) 173 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 178 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 179 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) & ! after scale factor at U-point 180 & + r_vvl * e3u(ji,jj,jk,Kaa) 181 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) & 182 & / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) 183 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) & 184 & / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 174 185 zWui = ( wi(ji,jj,jk ) + wi(ji+1,jj,jk ) ) / ze3ua 175 186 zWus = ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) / ze3ua … … 179 190 END_3D 180 191 END SELECT 181 DO_2D _00_00192 DO_2D( 0, 0, 0, 0 ) !* Surface boundary conditions 182 193 zwi(ji,jj,1) = 0._wp 183 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa) 184 zzws = - zdt * ( avm(ji+1,jj,2) + avm(ji ,jj,2) ) / ( ze3ua * e3uw(ji,jj,2,Kmm) ) * wumask(ji,jj,2) 194 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) & 195 & + r_vvl * e3u(ji,jj,1,Kaa) 196 zzws = - zdt * ( avm(ji+1,jj,2) + avm(ji ,jj,2) ) & 197 & / ( ze3ua * e3uw(ji,jj,2,Kmm) ) * wumask(ji,jj,2) 185 198 zWus = ( wi(ji ,jj,2) + wi(ji+1,jj,2) ) / ze3ua 186 199 zws(ji,jj,1 ) = zzws - zdt * MAX( zWus, 0._wp ) … … 190 203 SELECT CASE( nldf_dyn ) 191 204 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzu) 192 DO_3D_00_00( 1, jpkm1 ) 193 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 205 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 206 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) & 207 & + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 194 208 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) + akzu(ji,jj,jk ) ) & 195 209 & / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) … … 201 215 END_3D 202 216 CASE DEFAULT ! iso-level lateral mixing 203 DO_3D_00_00( 1, jpkm1 ) 204 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 205 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) 206 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 217 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 218 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,jk,Kmm) & 219 & + r_vvl * e3u(ji,jj,jk,Kaa) ! after scale factor at U-point 220 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) & 221 & / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) 222 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) & 223 & / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 207 224 zwi(ji,jj,jk) = zzwi 208 225 zws(ji,jj,jk) = zzws … … 210 227 END_3D 211 228 END SELECT 212 DO_2D _00_00229 DO_2D( 0, 0, 0, 0 ) !* Surface boundary conditions 213 230 zwi(ji,jj,1) = 0._wp 214 231 zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) … … 224 241 ! 225 242 IF ( ln_drgimp ) THEN ! implicit bottom friction 226 DO_2D _00_00243 DO_2D( 0, 0, 0, 0 ) 227 244 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 228 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) ! after scale factor at T-point 245 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) & 246 & + r_vvl * e3u(ji,jj,iku,Kaa) ! after scale factor at T-point 229 247 zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 230 248 END_2D 231 IF ( ln_isfcav ) THEN ! top friction (always implicit)232 DO_2D _00_00249 IF ( ln_isfcav.OR.ln_drgice_imp ) THEN ! top friction (always implicit) 250 DO_2D( 0, 0, 0, 0 ) 233 251 !!gm top Cd is masked (=0 outside cavities) no need of test on mik>=2 ==>> it has been suppressed 234 252 iku = miku(ji,jj) ! ocean top level at u- and v-points 235 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) + r_vvl * e3u(ji,jj,iku,Kaa) ! after scale factor at T-point 253 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,iku,Kmm) & 254 & + r_vvl * e3u(ji,jj,iku,Kaa) ! after scale factor at T-point 236 255 zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 237 256 END_2D … … 254 273 !----------------------------------------------------------------------- 255 274 ! 256 DO_3D _00_00( 2, jpkm1 )275 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) !== First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 (increasing k) == 257 276 zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1) 258 277 END_3D 259 278 ! 260 DO_2D_00_00 261 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa) 279 DO_2D( 0, 0, 0, 0 ) !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! 280 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) & 281 & + r_vvl * e3u(ji,jj,1,Kaa) 262 282 puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) ) & 263 283 & / ( ze3ua * rho0 ) * umask(ji,jj,1) 264 284 END_2D 265 DO_3D _00_00(2, jpkm1 )285 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 266 286 puu(ji,jj,jk,Kaa) = puu(ji,jj,jk,Kaa) - zwi(ji,jj,jk) / zwd(ji,jj,jk-1) * puu(ji,jj,jk-1,Kaa) 267 287 END_3D 268 288 ! 269 DO_2D _00_00289 DO_2D( 0, 0, 0, 0 ) !== thrid recurrence : SOLk = ( Lk - Uk * Ek+1 ) / Dk ==! 270 290 puu(ji,jj,jpkm1,Kaa) = puu(ji,jj,jpkm1,Kaa) / zwd(ji,jj,jpkm1) 271 291 END_2D 272 DO_3DS _00_00(jpk-2, 1, -1 )292 DO_3DS( 0, 0, 0, 0, jpk-2, 1, -1 ) 273 293 puu(ji,jj,jk,Kaa) = ( puu(ji,jj,jk,Kaa) - zws(ji,jj,jk) * puu(ji,jj,jk+1,Kaa) ) / zwd(ji,jj,jk) 274 294 END_3D … … 281 301 SELECT CASE( nldf_dyn ) 282 302 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzv) 283 DO_3D_00_00( 1, jpkm1 ) 284 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 303 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 304 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) & 305 & + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 285 306 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) + akzv(ji,jj,jk ) ) & 286 307 & / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) … … 294 315 END_3D 295 316 CASE DEFAULT ! iso-level lateral mixing 296 DO_3D_00_00( 1, jpkm1 ) 297 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 298 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) 299 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 317 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 318 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) & 319 & + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 320 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) & 321 & / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) 322 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) & 323 & / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 300 324 zWvi = ( wi(ji,jj,jk ) + wi(ji,jj+1,jk ) ) / ze3va 301 325 zWvs = ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) / ze3va … … 305 329 END_3D 306 330 END SELECT 307 DO_2D _00_00331 DO_2D( 0, 0, 0, 0 ) !* Surface boundary conditions 308 332 zwi(ji,jj,1) = 0._wp 309 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa) 310 zzws = - zdt * ( avm(ji,jj+1,2) + avm(ji,jj,2) ) / ( ze3va * e3vw(ji,jj,2,Kmm) ) * wvmask(ji,jj,2) 333 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) & 334 & + r_vvl * e3v(ji,jj,1,Kaa) 335 zzws = - zdt * ( avm(ji,jj+1,2) + avm(ji,jj,2) ) & 336 & / ( ze3va * e3vw(ji,jj,2,Kmm) ) * wvmask(ji,jj,2) 311 337 zWvs = ( wi(ji,jj ,2) + wi(ji,jj+1,2) ) / ze3va 312 338 zws(ji,jj,1 ) = zzws - zdt * MAX( zWvs, 0._wp ) … … 316 342 SELECT CASE( nldf_dyn ) 317 343 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzu) 318 DO_3D_00_00( 1, jpkm1 ) 319 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 344 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 345 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) & 346 & + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 320 347 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) + akzv(ji,jj,jk ) ) & 321 348 & / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) … … 327 354 END_3D 328 355 CASE DEFAULT ! iso-level lateral mixing 329 DO_3D_00_00( 1, jpkm1 ) 330 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 331 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) 332 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 356 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 357 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,jk,Kmm) & 358 & + r_vvl * e3v(ji,jj,jk,Kaa) ! after scale factor at V-point 359 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) & 360 & / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) 361 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) & 362 & / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 333 363 zwi(ji,jj,jk) = zzwi 334 364 zws(ji,jj,jk) = zzws … … 336 366 END_3D 337 367 END SELECT 338 DO_2D _00_00368 DO_2D( 0, 0, 0, 0 ) !* Surface boundary conditions 339 369 zwi(ji,jj,1) = 0._wp 340 370 zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) … … 349 379 ! 350 380 IF( ln_drgimp ) THEN 351 DO_2D _00_00381 DO_2D( 0, 0, 0, 0 ) 352 382 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 353 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) ! after scale factor at T-point 383 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) & 384 & + r_vvl * e3v(ji,jj,ikv,Kaa) ! after scale factor at T-point 354 385 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va 355 386 END_2D 356 IF ( ln_isfcav ) THEN357 DO_2D _00_00387 IF ( ln_isfcav.OR.ln_drgice_imp ) THEN 388 DO_2D( 0, 0, 0, 0 ) 358 389 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 359 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) + r_vvl * e3v(ji,jj,ikv,Kaa) ! after scale factor at T-point 390 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,ikv,Kmm) & 391 & + r_vvl * e3v(ji,jj,ikv,Kaa) ! after scale factor at T-point 360 392 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va 361 393 END_2D … … 378 410 !----------------------------------------------------------------------- 379 411 ! 380 DO_3D _00_00( 2, jpkm1 )412 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) !== First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 (increasing k) == 381 413 zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1) 382 414 END_3D 383 415 ! 384 DO_2D_00_00 385 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa) 416 DO_2D( 0, 0, 0, 0 ) !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! 417 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) & 418 & + r_vvl * e3v(ji,jj,1,Kaa) 386 419 pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) ) & 387 420 & / ( ze3va * rho0 ) * vmask(ji,jj,1) 388 421 END_2D 389 DO_3D _00_00(2, jpkm1 )422 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 390 423 pvv(ji,jj,jk,Kaa) = pvv(ji,jj,jk,Kaa) - zwi(ji,jj,jk) / zwd(ji,jj,jk-1) * pvv(ji,jj,jk-1,Kaa) 391 424 END_3D 392 425 ! 393 DO_2D _00_00426 DO_2D( 0, 0, 0, 0 ) !== third recurrence : SOLk = ( Lk - Uk * SOLk+1 ) / Dk ==! 394 427 pvv(ji,jj,jpkm1,Kaa) = pvv(ji,jj,jpkm1,Kaa) / zwd(ji,jj,jpkm1) 395 428 END_2D 396 DO_3DS _00_00(jpk-2, 1, -1 )429 DO_3DS( 0, 0, 0, 0, jpk-2, 1, -1 ) 397 430 pvv(ji,jj,jk,Kaa) = ( pvv(ji,jj,jk,Kaa) - zws(ji,jj,jk) * pvv(ji,jj,jk+1,Kaa) ) / zwd(ji,jj,jk) 398 431 END_3D
Note: See TracChangeset
for help on using the changeset viewer.
