- Timestamp:
- 2020-05-14T21:46:00+02:00 (4 years ago)
- Location:
- NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser
- Property svn:externals
-
old new 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 9 # SETTE 10 ^/utils/CI/sette@HEAD sette
-
- Property svn:externals
-
NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser/src/OCE/DYN/dynzdf.F90
r12178 r12928 37 37 38 38 !! * Substitutions 39 # include " vectopt_loop_substitute.h90"39 # include "do_loop_substitute.h90" 40 40 !!---------------------------------------------------------------------- 41 41 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 45 45 CONTAINS 46 46 47 SUBROUTINE dyn_zdf( kt )47 SUBROUTINE dyn_zdf( kt, Kbb, Kmm, Krhs, puu, pvv, Kaa ) 48 48 !!---------------------------------------------------------------------- 49 49 !! *** ROUTINE dyn_zdf *** … … 54 54 !! 55 55 !! ** Method : - Leap-Frog time stepping on all trends but the vertical mixing 56 !! u a = ub + 2*dt * uavector form or linear free surf.57 !! u a = ( e3u_b*ub + 2*dt * e3u_n*ua ) / e3u_aotherwise56 !! 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) otherwise 58 58 !! - update the after velocity with the implicit vertical mixing. 59 59 !! This requires to solver the following system: 60 !! u a = ua + 1/e3u_a dk+1[ mi(avm) / e3uw_adk[ua] ]60 !! u(after) = u(after) + 1/e3u(after) dk+1[ mi(avm) / e3uw(after) dk[ua] ] 61 61 !! with the following surface/top/bottom boundary condition: 62 62 !! surface: wind stress input (averaged over kt-1/2 & kt+1/2) 63 63 !! top & bottom : top stress (iceshelf-ocean) & bottom stress (cf zdfdrg.F90) 64 64 !! 65 !! ** Action : ( ua,va) after velocity65 !! ** Action : (puu(:,:,:,Kaa),pvv(:,:,:,Kaa)) after velocity 66 66 !!--------------------------------------------------------------------- 67 INTEGER, INTENT(in) :: kt ! ocean time-step index 67 INTEGER , INTENT( in ) :: kt ! ocean time-step index 68 INTEGER , INTENT( in ) :: Kbb, Kmm, Krhs, Kaa ! ocean time level indices 69 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 68 70 ! 69 71 INTEGER :: ji, jj, jk ! dummy loop indices … … 90 92 ENDIF 91 93 ENDIF 92 ! !* set time step93 IF( neuler == 0 .AND. kt == nit000 ) THEN ; r2dt = rdt ! = rdt (restart with Euler time stepping)94 ELSEIF( kt <= nit000 + 1 ) THEN ; r2dt = 2. * rdt ! = 2 rdt (leapfrog)95 ENDIF96 !97 94 ! !* explicit top/bottom drag case 98 IF( .NOT.ln_drgimp ) CALL zdf_drg_exp( kt, ub, vb, ua, va ) ! add top/bottom friction trend to (ua,va)95 IF( .NOT.ln_drgimp ) CALL zdf_drg_exp( kt, Kmm, puu(:,:,:,Kbb), pvv(:,:,:,Kbb), puu(:,:,:,Krhs), pvv(:,:,:,Krhs) ) ! add top/bottom friction trend to (puu(Kaa),pvv(Kaa)) 99 96 ! 100 97 ! 101 98 IF( l_trddyn ) THEN !* temporary save of ta and sa trends 102 99 ALLOCATE( ztrdu(jpi,jpj,jpk), ztrdv(jpi,jpj,jpk) ) 103 ztrdu(:,:,:) = ua(:,:,:)104 ztrdv(:,:,:) = va(:,:,:)105 ENDIF 106 ! 107 ! !== RHS: Leap-Frog time stepping on all trends but the vertical mixing ==! (put in ua,va)100 ztrdu(:,:,:) = puu(:,:,:,Krhs) 101 ztrdv(:,:,:) = pvv(:,:,:,Krhs) 102 ENDIF 103 ! 104 ! !== RHS: Leap-Frog time stepping on all trends but the vertical mixing ==! (put in puu(:,:,:,Kaa),pvv(:,:,:,Kaa)) 108 105 ! 109 106 ! ! time stepping except vertical diffusion 110 107 IF( ln_dynadv_vec .OR. ln_linssh ) THEN ! applied on velocity 111 108 DO jk = 1, jpkm1 112 ua(:,:,jk) = ( ub(:,:,jk) + r2dt * ua(:,:,jk) ) * umask(:,:,jk)113 va(:,:,jk) = ( vb(:,:,jk) + r2dt * va(:,:,jk) ) * vmask(:,:,jk)109 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) 114 111 END DO 115 112 ELSE ! applied on thickness weighted velocity 116 113 DO jk = 1, jpkm1 117 ua(:,:,jk) = ( e3u_b(:,:,jk) * ub(:,:,jk) &118 & + r 2dt * e3u_n(:,:,jk) * ua(:,:,jk) ) / e3u_a(:,:,jk) * umask(:,:,jk)119 va(:,:,jk) = ( e3v_b(:,:,jk) * vb(:,:,jk) &120 & + r 2dt * e3v_n(:,:,jk) * va(:,:,jk) ) / e3v_a(:,:,jk) * vmask(:,:,jk)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) 121 118 END DO 122 119 ENDIF … … 125 122 ! J. Chanut: The bottom stress is computed considering after barotropic velocities, which does 126 123 ! not lead to the effective stress seen over the whole barotropic loop. 127 ! G. Madec : in linear free surface, e3u _a = e3u_n = e3u_0, so systematic use of e3u_a124 ! G. Madec : in linear free surface, e3u(:,:,:,Kaa) = e3u(:,:,:,Kmm) = e3u_0, so systematic use of e3u(:,:,:,Kaa) 128 125 IF( ln_drgimp .AND. ln_dynspg_ts ) THEN 129 126 DO jk = 1, jpkm1 ! remove barotropic velocities 130 ua(:,:,jk) = ( ua(:,:,jk) - ua_b(:,:) ) * umask(:,:,jk)131 va(:,:,jk) = ( va(:,:,jk) - va_b(:,:) ) * vmask(:,:,jk)127 puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kaa) - uu_b(:,:,Kaa) ) * umask(:,:,jk) 128 pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kaa) - vv_b(:,:,Kaa) ) * vmask(:,:,jk) 132 129 END DO 133 DO jj = 2, jpjm1 ! Add bottom/top stress due to barotropic component only 134 DO ji = fs_2, fs_jpim1 ! vector opt. 135 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 136 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 137 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) 138 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) 139 ua(ji,jj,iku) = ua(ji,jj,iku) + r2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * ua_b(ji,jj) / ze3ua 140 va(ji,jj,ikv) = va(ji,jj,ikv) + r2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * va_b(ji,jj) / ze3va 141 END DO 142 END DO 130 DO_2D_00_00 131 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 132 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) 135 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 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 END_2D 143 138 IF( ln_isfcav ) THEN ! Ocean cavities (ISF) 144 DO jj = 2, jpjm1 145 DO ji = fs_2, fs_jpim1 ! vector opt. 146 iku = miku(ji,jj) ! top ocean level at u- and v-points 147 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 148 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) 149 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) 150 ua(ji,jj,iku) = ua(ji,jj,iku) + r2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * ua_b(ji,jj) / ze3ua 151 va(ji,jj,ikv) = va(ji,jj,ikv) + r2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * va_b(ji,jj) / ze3va 152 END DO 153 END DO 139 DO_2D_00_00 140 iku = miku(ji,jj) ! top ocean level at u- and v-points 141 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) 144 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 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 146 END_2D 154 147 END IF 155 148 ENDIF … … 158 151 ! 159 152 ! !* Matrix construction 160 zdt = r 2dt * 0.5153 zdt = rDt * 0.5 161 154 IF( ln_zad_Aimp ) THEN !! 162 155 SELECT CASE( nldf_dyn ) 163 156 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzu) 164 DO jk = 1, jpkm1 165 DO jj = 2, jpjm1 166 DO ji = fs_2, fs_jpim1 ! vector opt. 167 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk) ! after scale factor at U-point 168 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) + akzu(ji,jj,jk ) ) & 169 & / ( ze3ua * e3uw_n(ji,jj,jk ) ) * wumask(ji,jj,jk ) 170 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) ) & 171 & / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 172 zWui = ( wi(ji,jj,jk ) + wi(ji+1,jj,jk ) ) / ze3ua 173 zWus = ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) / ze3ua 174 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWui, 0._wp ) 175 zws(ji,jj,jk) = zzws - zdt * MAX( zWus, 0._wp ) 176 zwd(ji,jj,jk) = 1._wp - zzwi - zzws + zdt * ( MAX( zWui, 0._wp ) - MIN( zWus, 0._wp ) ) 177 END DO 178 END DO 179 END DO 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 159 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) + akzu(ji,jj,jk ) ) & 160 & / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) 161 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) ) & 162 & / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 163 zWui = ( wi(ji,jj,jk ) + wi(ji+1,jj,jk ) ) / ze3ua 164 zWus = ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) / ze3ua 165 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWui, 0._wp ) 166 zws(ji,jj,jk) = zzws - zdt * MAX( zWus, 0._wp ) 167 zwd(ji,jj,jk) = 1._wp - zzwi - zzws + zdt * ( MAX( zWui, 0._wp ) - MIN( zWus, 0._wp ) ) 168 END_3D 180 169 CASE DEFAULT ! iso-level lateral mixing 181 DO jk = 1, jpkm1 182 DO jj = 2, jpjm1 183 DO ji = fs_2, fs_jpim1 ! vector opt. 184 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk) ! after scale factor at U-point 185 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) / ( ze3ua * e3uw_n(ji,jj,jk ) ) * wumask(ji,jj,jk ) 186 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 187 zWui = ( wi(ji,jj,jk ) + wi(ji+1,jj,jk ) ) / ze3ua 188 zWus = ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) / ze3ua 189 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWui, 0._wp ) 190 zws(ji,jj,jk) = zzws - zdt * MAX( zWus, 0._wp ) 191 zwd(ji,jj,jk) = 1._wp - zzwi - zzws + zdt * ( MAX( zWui, 0._wp ) - MIN( zWus, 0._wp ) ) 192 END DO 193 END DO 194 END DO 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) 174 zWui = ( wi(ji,jj,jk ) + wi(ji+1,jj,jk ) ) / ze3ua 175 zWus = ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) / ze3ua 176 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWui, 0._wp ) 177 zws(ji,jj,jk) = zzws - zdt * MAX( zWus, 0._wp ) 178 zwd(ji,jj,jk) = 1._wp - zzwi - zzws + zdt * ( MAX( zWui, 0._wp ) - MIN( zWus, 0._wp ) ) 179 END_3D 195 180 END SELECT 196 DO jj = 2, jpjm1 !* Surface boundary conditions 197 DO ji = fs_2, fs_jpim1 ! vector opt. 198 zwi(ji,jj,1) = 0._wp 199 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,1) + r_vvl * e3u_a(ji,jj,1) 200 zzws = - zdt * ( avm(ji+1,jj,2) + avm(ji ,jj,2) ) / ( ze3ua * e3uw_n(ji,jj,2) ) * wumask(ji,jj,2) 201 zWus = ( wi(ji ,jj,2) + wi(ji+1,jj,2) ) / ze3ua 202 zws(ji,jj,1 ) = zzws - zdt * MAX( zWus, 0._wp ) 203 zwd(ji,jj,1 ) = 1._wp - zzws - zdt * ( MIN( zWus, 0._wp ) ) 204 END DO 205 END DO 181 DO_2D_00_00 182 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) 185 zWus = ( wi(ji ,jj,2) + wi(ji+1,jj,2) ) / ze3ua 186 zws(ji,jj,1 ) = zzws - zdt * MAX( zWus, 0._wp ) 187 zwd(ji,jj,1 ) = 1._wp - zzws - zdt * ( MIN( zWus, 0._wp ) ) 188 END_2D 206 189 ELSE 207 190 SELECT CASE( nldf_dyn ) 208 191 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzu) 209 DO jk = 1, jpkm1 210 DO jj = 2, jpjm1 211 DO ji = fs_2, fs_jpim1 ! vector opt. 212 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk) ! after scale factor at U-point 213 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) + akzu(ji,jj,jk ) ) & 214 & / ( ze3ua * e3uw_n(ji,jj,jk ) ) * wumask(ji,jj,jk ) 215 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) ) & 216 & / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 217 zwi(ji,jj,jk) = zzwi 218 zws(ji,jj,jk) = zzws 219 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 220 END DO 221 END DO 222 END DO 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 194 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) + akzu(ji,jj,jk ) ) & 195 & / ( ze3ua * e3uw(ji,jj,jk ,Kmm) ) * wumask(ji,jj,jk ) 196 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) ) & 197 & / ( ze3ua * e3uw(ji,jj,jk+1,Kmm) ) * wumask(ji,jj,jk+1) 198 zwi(ji,jj,jk) = zzwi 199 zws(ji,jj,jk) = zzws 200 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 201 END_3D 223 202 CASE DEFAULT ! iso-level lateral mixing 224 DO jk = 1, jpkm1 225 DO jj = 2, jpjm1 226 DO ji = fs_2, fs_jpim1 ! vector opt. 227 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk) ! after scale factor at U-point 228 zzwi = - zdt * ( avm(ji+1,jj,jk ) + avm(ji,jj,jk ) ) / ( ze3ua * e3uw_n(ji,jj,jk ) ) * wumask(ji,jj,jk ) 229 zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 230 zwi(ji,jj,jk) = zzwi 231 zws(ji,jj,jk) = zzws 232 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 233 END DO 234 END DO 235 END DO 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) 207 zwi(ji,jj,jk) = zzwi 208 zws(ji,jj,jk) = zzws 209 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 210 END_3D 236 211 END SELECT 237 DO jj = 2, jpjm1 !* Surface boundary conditions 238 DO ji = fs_2, fs_jpim1 ! vector opt. 239 zwi(ji,jj,1) = 0._wp 240 zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 241 END DO 242 END DO 212 DO_2D_00_00 213 zwi(ji,jj,1) = 0._wp 214 zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 215 END_2D 243 216 ENDIF 244 217 ! … … 251 224 ! 252 225 IF ( ln_drgimp ) THEN ! implicit bottom friction 253 DO jj = 2, jpjm1 254 DO ji = 2, jpim1 255 iku = mbku(ji,jj) ! ocean bottom level at u- and v-points 256 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) ! after scale factor at T-point 257 zwd(ji,jj,iku) = zwd(ji,jj,iku) - r2dt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 258 END DO 259 END DO 226 DO_2D_00_00 227 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 229 zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) / ze3ua 230 END_2D 260 231 IF ( ln_isfcav ) THEN ! top friction (always implicit) 261 DO jj = 2, jpjm1 262 DO ji = 2, jpim1 263 !!gm top Cd is masked (=0 outside cavities) no need of test on mik>=2 ==>> it has been suppressed 264 iku = miku(ji,jj) ! ocean top level at u- and v-points 265 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,iku) + r_vvl * e3u_a(ji,jj,iku) ! after scale factor at T-point 266 zwd(ji,jj,iku) = zwd(ji,jj,iku) - r2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 267 END DO 268 END DO 232 DO_2D_00_00 233 !!gm top Cd is masked (=0 outside cavities) no need of test on mik>=2 ==>> it has been suppressed 234 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 236 zwd(ji,jj,iku) = zwd(ji,jj,iku) - rDt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3ua 237 END_2D 269 238 END IF 270 239 ENDIF … … 282 251 ! m is decomposed in the product of an upper and a lower triangular matrix 283 252 ! The 3 diagonal terms are in 2d arrays: zwd, zws, zwi 284 ! The solution (the after velocity) is in ua253 ! The solution (the after velocity) is in puu(:,:,:,Kaa) 285 254 !----------------------------------------------------------------------- 286 255 ! 287 DO jk = 2, jpkm1 !== First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 (increasing k) == 288 DO jj = 2, jpjm1 289 DO ji = fs_2, fs_jpim1 ! vector opt. 290 zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1) 291 END DO 292 END DO 293 END DO 294 ! 295 DO jj = 2, jpjm1 !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! 296 DO ji = fs_2, fs_jpim1 ! vector opt. 297 ze3ua = ( 1._wp - r_vvl ) * e3u_n(ji,jj,1) + r_vvl * e3u_a(ji,jj,1) 298 ua(ji,jj,1) = ua(ji,jj,1) + r2dt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) ) & 299 & / ( ze3ua * rau0 ) * umask(ji,jj,1) 300 END DO 301 END DO 302 DO jk = 2, jpkm1 303 DO jj = 2, jpjm1 304 DO ji = fs_2, fs_jpim1 305 ua(ji,jj,jk) = ua(ji,jj,jk) - zwi(ji,jj,jk) / zwd(ji,jj,jk-1) * ua(ji,jj,jk-1) 306 END DO 307 END DO 308 END DO 309 ! 310 DO jj = 2, jpjm1 !== thrid recurrence : SOLk = ( Lk - Uk * Ek+1 ) / Dk ==! 311 DO ji = fs_2, fs_jpim1 ! vector opt. 312 ua(ji,jj,jpkm1) = ua(ji,jj,jpkm1) / zwd(ji,jj,jpkm1) 313 END DO 314 END DO 315 DO jk = jpk-2, 1, -1 316 DO jj = 2, jpjm1 317 DO ji = fs_2, fs_jpim1 318 ua(ji,jj,jk) = ( ua(ji,jj,jk) - zws(ji,jj,jk) * ua(ji,jj,jk+1) ) / zwd(ji,jj,jk) 319 END DO 320 END DO 321 END DO 256 DO_3D_00_00( 2, jpkm1 ) 257 zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1) 258 END_3D 259 ! 260 DO_2D_00_00 261 ze3ua = ( 1._wp - r_vvl ) * e3u(ji,jj,1,Kmm) + r_vvl * e3u(ji,jj,1,Kaa) 262 puu(ji,jj,1,Kaa) = puu(ji,jj,1,Kaa) + rDt * 0.5_wp * ( utau_b(ji,jj) + utau(ji,jj) ) & 263 & / ( ze3ua * rho0 ) * umask(ji,jj,1) 264 END_2D 265 DO_3D_00_00( 2, jpkm1 ) 266 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 END_3D 268 ! 269 DO_2D_00_00 270 puu(ji,jj,jpkm1,Kaa) = puu(ji,jj,jpkm1,Kaa) / zwd(ji,jj,jpkm1) 271 END_2D 272 DO_3DS_00_00( jpk-2, 1, -1 ) 273 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 END_3D 322 275 ! 323 276 ! !== Vertical diffusion on v ==! 324 277 ! 325 278 ! !* Matrix construction 326 zdt = r 2dt * 0.5279 zdt = rDt * 0.5 327 280 IF( ln_zad_Aimp ) THEN !! 328 281 SELECT CASE( nldf_dyn ) 329 282 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzv) 330 DO jk = 1, jpkm1 331 DO jj = 2, jpjm1 332 DO ji = fs_2, fs_jpim1 ! vector opt. 333 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk) ! after scale factor at V-point 334 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) + akzv(ji,jj,jk ) ) & 335 & / ( ze3va * e3vw_n(ji,jj,jk ) ) * wvmask(ji,jj,jk ) 336 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) + akzv(ji,jj,jk+1) ) & 337 & / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 338 zWvi = ( wi(ji,jj,jk ) + wi(ji,jj+1,jk ) ) / ze3va 339 zWvs = ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) / ze3va 340 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWvi, 0._wp ) 341 zws(ji,jj,jk) = zzws - zdt * MAX( zWvs, 0._wp ) 342 zwd(ji,jj,jk) = 1._wp - zzwi - zzws - zdt * ( - MAX( zWvi, 0._wp ) + MIN( zWvs, 0._wp ) ) 343 END DO 344 END DO 345 END DO 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 285 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) + akzv(ji,jj,jk ) ) & 286 & / ( ze3va * e3vw(ji,jj,jk ,Kmm) ) * wvmask(ji,jj,jk ) 287 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) + akzv(ji,jj,jk+1) ) & 288 & / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 289 zWvi = ( wi(ji,jj,jk ) + wi(ji,jj+1,jk ) ) / ze3va 290 zWvs = ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) / ze3va 291 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWvi, 0._wp ) 292 zws(ji,jj,jk) = zzws - zdt * MAX( zWvs, 0._wp ) 293 zwd(ji,jj,jk) = 1._wp - zzwi - zzws - zdt * ( - MAX( zWvi, 0._wp ) + MIN( zWvs, 0._wp ) ) 294 END_3D 346 295 CASE DEFAULT ! iso-level lateral mixing 347 DO jk = 1, jpkm1 348 DO jj = 2, jpjm1 349 DO ji = fs_2, fs_jpim1 ! vector opt. 350 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk) ! after scale factor at V-point 351 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) / ( ze3va * e3vw_n(ji,jj,jk ) ) * wvmask(ji,jj,jk ) 352 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 353 zWvi = ( wi(ji,jj,jk ) + wi(ji,jj+1,jk ) ) / ze3va 354 zWvs = ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) / ze3va 355 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWvi, 0._wp ) 356 zws(ji,jj,jk) = zzws - zdt * MAX( zWvs, 0._wp ) 357 zwd(ji,jj,jk) = 1._wp - zzwi - zzws - zdt * ( - MAX( zWvi, 0._wp ) + MIN( zWvs, 0._wp ) ) 358 END DO 359 END DO 360 END DO 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) 300 zWvi = ( wi(ji,jj,jk ) + wi(ji,jj+1,jk ) ) / ze3va 301 zWvs = ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) / ze3va 302 zwi(ji,jj,jk) = zzwi + zdt * MIN( zWvi, 0._wp ) 303 zws(ji,jj,jk) = zzws - zdt * MAX( zWvs, 0._wp ) 304 zwd(ji,jj,jk) = 1._wp - zzwi - zzws - zdt * ( - MAX( zWvi, 0._wp ) + MIN( zWvs, 0._wp ) ) 305 END_3D 361 306 END SELECT 362 DO jj = 2, jpjm1 !* Surface boundary conditions 363 DO ji = fs_2, fs_jpim1 ! vector opt. 364 zwi(ji,jj,1) = 0._wp 365 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,1) + r_vvl * e3v_a(ji,jj,1) 366 zzws = - zdt * ( avm(ji,jj+1,2) + avm(ji,jj,2) ) / ( ze3va * e3vw_n(ji,jj,2) ) * wvmask(ji,jj,2) 367 zWvs = ( wi(ji,jj ,2) + wi(ji,jj+1,2) ) / ze3va 368 zws(ji,jj,1 ) = zzws - zdt * MAX( zWvs, 0._wp ) 369 zwd(ji,jj,1 ) = 1._wp - zzws - zdt * ( MIN( zWvs, 0._wp ) ) 370 END DO 371 END DO 307 DO_2D_00_00 308 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) 311 zWvs = ( wi(ji,jj ,2) + wi(ji,jj+1,2) ) / ze3va 312 zws(ji,jj,1 ) = zzws - zdt * MAX( zWvs, 0._wp ) 313 zwd(ji,jj,1 ) = 1._wp - zzws - zdt * ( MIN( zWvs, 0._wp ) ) 314 END_2D 372 315 ELSE 373 316 SELECT CASE( nldf_dyn ) 374 317 CASE( np_lap_i ) ! rotated lateral mixing: add its vertical mixing (akzu) 375 DO jk = 1, jpkm1 376 DO jj = 2, jpjm1 377 DO ji = fs_2, fs_jpim1 ! vector opt. 378 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk) ! after scale factor at V-point 379 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) + akzv(ji,jj,jk ) ) & 380 & / ( ze3va * e3vw_n(ji,jj,jk ) ) * wvmask(ji,jj,jk ) 381 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) + akzv(ji,jj,jk+1) ) & 382 & / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 383 zwi(ji,jj,jk) = zzwi 384 zws(ji,jj,jk) = zzws 385 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 386 END DO 387 END DO 388 END DO 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 320 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) + akzv(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) + akzv(ji,jj,jk+1) ) & 323 & / ( ze3va * e3vw(ji,jj,jk+1,Kmm) ) * wvmask(ji,jj,jk+1) 324 zwi(ji,jj,jk) = zzwi 325 zws(ji,jj,jk) = zzws 326 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 327 END_3D 389 328 CASE DEFAULT ! iso-level lateral mixing 390 DO jk = 1, jpkm1 391 DO jj = 2, jpjm1 392 DO ji = fs_2, fs_jpim1 ! vector opt. 393 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk) ! after scale factor at V-point 394 zzwi = - zdt * ( avm(ji,jj+1,jk ) + avm(ji,jj,jk ) ) / ( ze3va * e3vw_n(ji,jj,jk ) ) * wvmask(ji,jj,jk ) 395 zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 396 zwi(ji,jj,jk) = zzwi 397 zws(ji,jj,jk) = zzws 398 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 399 END DO 400 END DO 401 END DO 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) 333 zwi(ji,jj,jk) = zzwi 334 zws(ji,jj,jk) = zzws 335 zwd(ji,jj,jk) = 1._wp - zzwi - zzws 336 END_3D 402 337 END SELECT 403 DO jj = 2, jpjm1 !* Surface boundary conditions 404 DO ji = fs_2, fs_jpim1 ! vector opt. 405 zwi(ji,jj,1) = 0._wp 406 zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 407 END DO 408 END DO 338 DO_2D_00_00 339 zwi(ji,jj,1) = 0._wp 340 zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 341 END_2D 409 342 ENDIF 410 343 ! … … 416 349 ! 417 350 IF( ln_drgimp ) THEN 418 DO jj = 2, jpjm1 419 DO ji = 2, jpim1 420 ikv = mbkv(ji,jj) ! (deepest ocean u- and v-points) 421 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) ! after scale factor at T-point 422 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - r2dt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va 423 END DO 424 END DO 351 DO_2D_00_00 352 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 354 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) / ze3va 355 END_2D 425 356 IF ( ln_isfcav ) THEN 426 DO jj = 2, jpjm1 427 DO ji = 2, jpim1 428 ikv = mikv(ji,jj) ! (first wet ocean u- and v-points) 429 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,ikv) + r_vvl * e3v_a(ji,jj,ikv) ! after scale factor at T-point 430 zwd(ji,jj,iku) = zwd(ji,jj,iku) - r2dt * 0.5*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) / ze3va 431 END DO 432 END DO 357 DO_2D_00_00 358 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 360 zwd(ji,jj,ikv) = zwd(ji,jj,ikv) - rDt * 0.5*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) / ze3va 361 END_2D 433 362 ENDIF 434 363 ENDIF … … 449 378 !----------------------------------------------------------------------- 450 379 ! 451 DO jk = 2, jpkm1 !== First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 (increasing k) == 452 DO jj = 2, jpjm1 453 DO ji = fs_2, fs_jpim1 ! vector opt. 454 zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1) 455 END DO 456 END DO 457 END DO 458 ! 459 DO jj = 2, jpjm1 !== second recurrence: SOLk = RHSk - Lk / Dk-1 Lk-1 ==! 460 DO ji = fs_2, fs_jpim1 ! vector opt. 461 ze3va = ( 1._wp - r_vvl ) * e3v_n(ji,jj,1) + r_vvl * e3v_a(ji,jj,1) 462 va(ji,jj,1) = va(ji,jj,1) + r2dt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) ) & 463 & / ( ze3va * rau0 ) * vmask(ji,jj,1) 464 END DO 465 END DO 466 DO jk = 2, jpkm1 467 DO jj = 2, jpjm1 468 DO ji = fs_2, fs_jpim1 ! vector opt. 469 va(ji,jj,jk) = va(ji,jj,jk) - zwi(ji,jj,jk) / zwd(ji,jj,jk-1) * va(ji,jj,jk-1) 470 END DO 471 END DO 472 END DO 473 ! 474 DO jj = 2, jpjm1 !== third recurrence : SOLk = ( Lk - Uk * SOLk+1 ) / Dk ==! 475 DO ji = fs_2, fs_jpim1 ! vector opt. 476 va(ji,jj,jpkm1) = va(ji,jj,jpkm1) / zwd(ji,jj,jpkm1) 477 END DO 478 END DO 479 DO jk = jpk-2, 1, -1 480 DO jj = 2, jpjm1 481 DO ji = fs_2, fs_jpim1 482 va(ji,jj,jk) = ( va(ji,jj,jk) - zws(ji,jj,jk) * va(ji,jj,jk+1) ) / zwd(ji,jj,jk) 483 END DO 484 END DO 485 END DO 380 DO_3D_00_00( 2, jpkm1 ) 381 zwd(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwd(ji,jj,jk-1) 382 END_3D 383 ! 384 DO_2D_00_00 385 ze3va = ( 1._wp - r_vvl ) * e3v(ji,jj,1,Kmm) + r_vvl * e3v(ji,jj,1,Kaa) 386 pvv(ji,jj,1,Kaa) = pvv(ji,jj,1,Kaa) + rDt * 0.5_wp * ( vtau_b(ji,jj) + vtau(ji,jj) ) & 387 & / ( ze3va * rho0 ) * vmask(ji,jj,1) 388 END_2D 389 DO_3D_00_00( 2, jpkm1 ) 390 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 END_3D 392 ! 393 DO_2D_00_00 394 pvv(ji,jj,jpkm1,Kaa) = pvv(ji,jj,jpkm1,Kaa) / zwd(ji,jj,jpkm1) 395 END_2D 396 DO_3DS_00_00( jpk-2, 1, -1 ) 397 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 END_3D 486 399 ! 487 400 IF( l_trddyn ) THEN ! save the vertical diffusive trends for further diagnostics 488 ztrdu(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) / r2dt - ztrdu(:,:,:)489 ztrdv(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) / r2dt - ztrdv(:,:,:)490 CALL trd_dyn( ztrdu, ztrdv, jpdyn_zdf, kt )401 ztrdu(:,:,:) = ( puu(:,:,:,Kaa) - puu(:,:,:,Kbb) ) / rDt - ztrdu(:,:,:) 402 ztrdv(:,:,:) = ( pvv(:,:,:,Kaa) - pvv(:,:,:,Kbb) ) / rDt - ztrdv(:,:,:) 403 CALL trd_dyn( ztrdu, ztrdv, jpdyn_zdf, kt, Kmm ) 491 404 DEALLOCATE( ztrdu, ztrdv ) 492 405 ENDIF 493 406 ! ! print mean trends (used for debugging) 494 IF( ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' zdf - Ua: ', mask1=umask, &495 & tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )407 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Kaa), clinfo1=' zdf - Ua: ', mask1=umask, & 408 & tab3d_2=pvv(:,:,:,Kaa), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 496 409 ! 497 410 IF( ln_timing ) CALL timing_stop('dyn_zdf')
Note: See TracChangeset
for help on using the changeset viewer.