Changeset 2807 for branches/2011/dev_r2802_NOCS_vvlfix
- Timestamp:
- 2011-07-19T18:35:40+02:00 (13 years ago)
- Location:
- branches/2011/dev_r2802_NOCS_vvlfix/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 3 edited
-
DOM/domvvl.F90 (modified) (7 diffs)
-
DOM/istate.F90 (modified) (2 diffs)
-
DYN/dynnxt.F90 (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2802_NOCS_vvlfix/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r2779 r2807 24 24 PRIVATE 25 25 26 PUBLIC dom_vvl ! called by domain.F9027 PUBLIC dom_vvl_ alloc ! called by nemogcm.F9028 29 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ee_t, ee_u, ee_v, ee_f !: ??? 30 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: mut , muu , muv , muf !: ???26 PUBLIC dom_vvl ! called by domain.F90 27 PUBLIC dom_vvl_2 ! called by domain.F90 28 PUBLIC dom_vvl_alloc ! called by nemogcm.F90 29 30 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: mut , muu , muv , muf !: 1/H_0 at t-,u-,v-,f-points 31 31 32 32 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: r2dt ! vertical profile time-step, = 2 rdttra … … 49 49 ! 50 50 ALLOCATE( mut (jpi,jpj,jpk) , muu (jpi,jpj,jpk) , muv (jpi,jpj,jpk) , muf (jpi,jpj,jpk) , & 51 & ee_t(jpi,jpj) , ee_u(jpi,jpj) , ee_v(jpi,jpj) , ee_f(jpi,jpj) , &52 51 & r2dt (jpk) , STAT=dom_vvl_alloc ) 53 52 ! … … 62 61 !! *** ROUTINE dom_vvl *** 63 62 !! 64 !! ** Purpose : compute coefficients muX at T-U-V-F points to spread 65 !! ssh over the whole water column (scale factors) 63 !! ** Purpose : compute mu coefficients at t-, u-, v- and f-points to 64 !! spread ssh over the whole water column (scale factors) 65 !! set the before and now ssh at u- and v-points 66 !! (also f-point in now case) 66 67 !!---------------------------------------------------------------------- 67 68 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 68 USE wrk_nemo, ONLY: z s_t => wrk_2d_1 , zs_u_1 => wrk_2d_2 , zs_v_1 => wrk_2d_3! 2D workspace69 USE wrk_nemo, ONLY: zee_t => wrk_2d_1, zee_u => wrk_2d_2, zee_v => wrk_2d_3, zee_f => wrk_2d_4 ! 2D workspace 69 70 ! 70 71 INTEGER :: ji, jj, jk ! dummy loop indices 71 REAL(wp) :: zcoefu , zcoefv , zcoeff! local scalars72 REAL(wp) :: zv _t_ij, zv_t_ip1j, zv_t_ijp1, zv_t_ip1jp1 ! - -73 !!---------------------------------------------------------------------- 74 75 IF( wrk_in_use(2, 1,2,3 ) ) THEN72 REAL(wp) :: zcoefu, zcoefv , zcoeff ! local scalars 73 REAL(wp) :: zvt , zvt_ip1, zvt_jp1, zvt_ip1jp1 ! - - 74 !!---------------------------------------------------------------------- 75 76 IF( wrk_in_use(2, 1,2,3,4) ) THEN 76 77 CALL ctl_stop('dom_vvl: requested workspace arrays unavailable') ; RETURN 77 78 ENDIF … … 97 98 98 99 ! !== mu computation ==! 99 ee_t(:,:) = fse3t_0(:,:,1) ! Lower bound : thickness of the first model level100 ee_u(:,:) = fse3u_0(:,:,1)101 ee_v(:,:) = fse3v_0(:,:,1)102 ee_f(:,:) = fse3f_0(:,:,1)100 zee_t(:,:) = fse3t_0(:,:,1) ! Lower bound : thickness of the first model level 101 zee_u(:,:) = fse3u_0(:,:,1) 102 zee_v(:,:) = fse3v_0(:,:,1) 103 zee_f(:,:) = fse3f_0(:,:,1) 103 104 DO jk = 2, jpkm1 ! Sum of the masked vertical scale factors 104 ee_t(:,:) =ee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk)105 ee_u(:,:) =ee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk)106 ee_v(:,:) =ee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk)105 zee_t(:,:) = zee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk) 106 zee_u(:,:) = zee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk) 107 zee_v(:,:) = zee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk) 107 108 DO jj = 1, jpjm1 ! f-point : fmask=shlat at coasts, use the product of umask 108 ee_f(:,jj) =ee_f(:,jj) + fse3f_0(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk)109 zee_f(:,jj) = zee_f(:,jj) + fse3f_0(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk) 109 110 END DO 110 111 END DO 111 112 ! ! Compute and mask the inverse of the local depth at T, U, V and F points 112 ee_t(:,:) = 1. /ee_t(:,:) * tmask(:,:,1)113 ee_u(:,:) = 1. /ee_u(:,:) * umask(:,:,1)114 ee_v(:,:) = 1. /ee_v(:,:) * vmask(:,:,1)113 zee_t(:,:) = 1._wp / zee_t(:,:) * tmask(:,:,1) 114 zee_u(:,:) = 1._wp / zee_u(:,:) * umask(:,:,1) 115 zee_v(:,:) = 1._wp / zee_v(:,:) * vmask(:,:,1) 115 116 DO jj = 1, jpjm1 ! f-point case fmask cannot be used 116 ee_f(:,jj) = 1. /ee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1)117 END DO 118 CALL lbc_lnk( ee_f, 'F', 1. )! lateral boundary condition on ee_f117 zee_f(:,jj) = 1._wp / zee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1) 118 END DO 119 CALL lbc_lnk( zee_f, 'F', 1. ) ! lateral boundary condition on ee_f 119 120 ! 120 121 DO jk = 1, jpk ! mu coefficients 121 mut(:,:,jk) = ee_t(:,:) * tmask(:,:,jk) ! T-point at T levels122 muu(:,:,jk) = ee_u(:,:) * umask(:,:,jk) ! U-point at T levels123 muv(:,:,jk) = ee_v(:,:) * vmask(:,:,jk) ! V-point at T levels122 mut(:,:,jk) = zee_t(:,:) * tmask(:,:,jk) ! T-point at T levels 123 muu(:,:,jk) = zee_u(:,:) * umask(:,:,jk) ! U-point at T levels 124 muv(:,:,jk) = zee_v(:,:) * vmask(:,:,jk) ! V-point at T levels 124 125 END DO 125 126 DO jk = 1, jpk ! F-point : fmask=shlat at coasts, use the product of umask 126 127 DO jj = 1, jpjm1 127 muf(:,jj,jk) = ee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk) ! at T levels128 END DO 129 muf(:,jpj,jk) = 0. e0128 muf(:,jj,jk) = zee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk) ! at T levels 129 END DO 130 muf(:,jpj,jk) = 0._wp 130 131 END DO 131 132 CALL lbc_lnk( muf, 'F', 1. ) ! lateral boundary condition … … 139 140 END DO 140 141 141 ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations142 ! for ssh and scale factors143 zs_t (:,:) = e1t(:,:) * e2t(:,:)144 zs_u_1(:,:) = 0.5 / ( e1u(:,:) * e2u(:,:) )145 zs_v_1(:,:) = 0.5 / ( e1v(:,:) * e2v(:,:) )146 147 142 DO jj = 1, jpjm1 ! initialise before and now Sea Surface Height at u-, v-, f-points 148 143 DO ji = 1, jpim1 ! NO vector opt. 149 zcoefu = umask(ji,jj,1) * zs_u_1(ji,jj)150 zcoefv = vmask(ji,jj,1) * zs_v_1(ji,jj)151 zcoeff = 0. 5 * umask(ji,jj,1) * umask(ji,jj+1,1) / ( e1f(ji,jj) * e2f(ji,jj))152 ! before fields153 zv _t_ij = zs_t(ji ,jj ) * sshb(ji ,jj )154 zv _t_ip1j = zs_t(ji+1,jj ) * sshb(ji+1,jj )155 zv _t_ijp1 = zs_t(ji ,jj+1) * sshb(ji ,jj+1)156 sshu_b(ji,jj) = zcoefu * ( zv _t_ij + zv_t_ip1j)157 sshv_b(ji,jj) = zcoefv * ( zv _t_ij + zv_t_ijp1 )158 ! now fields159 zv _t_ij = zs_t(ji ,jj ) * sshn(ji ,jj )160 zv _t_ip1j = zs_t(ji+1,jj ) * sshn(ji+1,jj )161 zv _t_ijp1 = zs_t(ji ,jj+1) * sshn(ji ,jj+1)162 zv _t_ip1jp1 = zs_t(ji ,jj+1) * sshn(ji,jj+1)163 sshu_n(ji,jj) = zcoefu * ( zv _t_ij + zv_t_ip1j)164 sshv_n(ji,jj) = zcoefv * ( zv _t_ij + zv_t_ijp1 )165 sshf_n(ji,jj) = zcoeff * ( zv _t_ij + zv_t_ip1j + zv_t_ijp1 + zv_t_ip1jp1 )144 zcoefu = 0.50_wp / ( e1u(ji,jj) * e2u(ji,jj) ) * umask(ji,jj,1) 145 zcoefv = 0.50_wp / ( e1v(ji,jj) * e2v(ji,jj) ) * vmask(ji,jj,1) 146 zcoeff = 0.25_wp / ( e1f(ji,jj) * e2f(ji,jj) ) * umask(ji,jj,1) * umask(ji,jj+1,1) 147 ! 148 zvt = e1e2t(ji ,jj ) * sshb(ji ,jj ) ! before fields 149 zvt_ip1 = e1e2t(ji+1,jj ) * sshb(ji+1,jj ) 150 zvt_jp1 = e1e2t(ji ,jj+1) * sshb(ji ,jj+1) 151 sshu_b(ji,jj) = zcoefu * ( zvt + zvt_ip1 ) 152 sshv_b(ji,jj) = zcoefv * ( zvt + zvt_jp1 ) 153 ! 154 zvt = e1e2t(ji ,jj ) * sshn(ji ,jj ) ! now fields 155 zvt_ip1 = e1e2t(ji+1,jj ) * sshn(ji+1,jj ) 156 zvt_jp1 = e1e2t(ji ,jj+1) * sshn(ji ,jj+1) 157 zvt_ip1jp1 = e1e2t(ji+1,jj+1) * sshn(ji+1,jj+1) 158 sshu_n(ji,jj) = zcoefu * ( zvt + zvt_ip1 ) 159 sshv_n(ji,jj) = zcoefv * ( zvt + zvt_jp1 ) 160 sshf_n(ji,jj) = zcoeff * ( zvt + zvt_ip1 + zvt_jp1 + zvt_ip1jp1 ) 166 161 END DO 167 162 END DO … … 169 164 CALL lbc_lnk( sshv_n, 'V', 1. ) ; CALL lbc_lnk( sshv_b, 'V', 1. ) 170 165 CALL lbc_lnk( sshf_n, 'F', 1. ) 171 172 ! initialise before scale factors at (u/v)-points 173 ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 174 DO jk = 1, jpkm1 175 DO jj = 1, jpjm1 176 DO ji = 1, jpim1 177 zv_t_ij = zs_t(ji ,jj ) * fse3t_b(ji ,jj ,jk) 178 zv_t_ip1j = zs_t(ji+1,jj ) * fse3t_b(ji+1,jj ,jk) 179 zv_t_ijp1 = zs_t(ji ,jj+1) * fse3t_b(ji ,jj+1,jk) 180 fse3u_b(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 181 fse3v_b(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 166 ! 167 IF( wrk_not_released(2, 1,2,3,4) ) CALL ctl_stop('dom_vvl: failed to release workspace arrays') 168 ! 169 END SUBROUTINE dom_vvl 170 171 172 SUBROUTINE dom_vvl_2( kt, pe3u_b, pe3v_b ) 173 !!---------------------------------------------------------------------- 174 !! *** ROUTINE dom_vvl_2 *** 175 !! 176 !! ** Purpose : compute the vertical scale factors at u- and v-points 177 !! in variable volume case. 178 !! 179 !! ** Method : In variable volume case (non linear sea surface) the 180 !! the vertical scale factor at velocity points is computed 181 !! as the average of the cell surface weighted e3t. 182 !! It uses the sea surface heigth so it have to be initialized 183 !! after ssh is read/set 184 !!---------------------------------------------------------------------- 185 INTEGER , INTENT(in ) :: kt ! ocean time-step index 186 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pe3u_b, pe3v_b ! before vertical scale factor at u- & v-pts 187 ! 188 INTEGER :: ji, jj, jk ! dummy loop indices 189 INTEGER :: iku, ikv ! local integers 190 REAL(wp) :: zvt ! local scalars 191 !!---------------------------------------------------------------------- 192 193 IF( lwp .AND. kt == nit000 ) THEN 194 WRITE(numout,*) 195 WRITE(numout,*) 'dom_vvl_2 : Variable volume, fse3t_b initialization' 196 WRITE(numout,*) '~~~~~~~~~ ' 197 pe3u_b(:,:,jpk) = fse3u_0(:,:,jpk) 198 pe3v_b(:,:,jpk) = fse3u_0(:,:,jpk) 199 ENDIF 200 201 DO jk = 1, jpkm1 ! set the before scale factors at u- & v-points 202 DO jj = 2, jpjm1 203 DO ji = fs_2, fs_jpim1 204 zvt = fse3t_b(ji,jj,jk) * e1e2t(ji,jj) 205 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1e2t(ji+1,jj) ) / ( e1u(ji,jj) * e2u(ji,jj) ) 206 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e1e2t(ji,jj+1) ) / ( e1v(ji,jj) * e2v(ji,jj) ) 182 207 END DO 183 208 END DO 184 209 END DO 185 CALL lbc_lnk( fse3u_b(:,:,:), 'U', 1. ) ! lateral boundary conditions 186 CALL lbc_lnk( fse3v_b(:,:,:), 'V', 1. ) 187 ! Add initial scale factor to scale factor anomaly 188 fse3u_b(:,:,:) = fse3u_b(:,:,:) + fse3u_0(:,:,:) 189 fse3v_b(:,:,:) = fse3v_b(:,:,:) + fse3v_0(:,:,:) 190 ! 191 IF( wrk_not_released(2, 1,2,3) ) CALL ctl_stop('dom_vvl: failed to release workspace arrays') 192 ! 193 END SUBROUTINE dom_vvl 194 210 IF( ln_zps ) THEN ! minimum of the e3t at partial cell level 211 DO jj = 2, jpjm1 212 DO ji = fs_2, fs_jpim1 213 iku = mbku(ji,jj) 214 ikv = mbkv(ji,jj) 215 pe3u_b(ji,jj,iku) = MIN( fse3t_b(ji,jj,iku), fse3t_b(ji+1,jj ,iku) ) 216 pe3v_b(ji,jj,ikv) = MIN( fse3t_b(ji,jj,ikv), fse3t_b(ji ,jj+1,ikv) ) 217 END DO 218 END DO 219 ENDIF 220 pe3u_b(:,:,:) = pe3u_b(:,:,:) - fse3u_0(:,:,:) ! anomaly to avoid zero along closed boundary/extra halos 221 pe3v_b(:,:,:) = pe3v_b(:,:,:) - fse3v_0(:,:,:) 222 CALL lbc_lnk( pe3u_b(:,:,:), 'U', 1. ) ! lateral boundary conditions 223 CALL lbc_lnk( pe3v_b(:,:,:), 'V', 1. ) 224 pe3u_b(:,:,:) = pe3u_b(:,:,:) + fse3u_0(:,:,:) ! recover the full scale factor 225 pe3v_b(:,:,:) = pe3v_b(:,:,:) + fse3v_0(:,:,:) 226 ! 227 END SUBROUTINE dom_vvl_2 228 195 229 #else 196 230 !!---------------------------------------------------------------------- … … 200 234 SUBROUTINE dom_vvl 201 235 END SUBROUTINE dom_vvl 236 SUBROUTINE dom_vvl_2 237 END SUBROUTINE dom_vvl_2 202 238 #endif 203 239 -
branches/2011/dev_r2802_NOCS_vvlfix/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90
r2777 r2807 83 83 neuler = 1 ! Set time-step indicator at nit000 (leap-frog) 84 84 CALL rst_read ! Read the restart file 85 ! ! define e3u_b, e3v_b from e3t_b read in restart file 86 CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 85 87 CALL tra_swap ! swap 3D arrays (t,s) in a 4D array (ts) 86 88 CALL day_init ! model calendar (using both namelist and restart infos) … … 92 94 CALL day_init ! model calendar (using both namelist and restart infos) 93 95 ! ! Initialization of ocean to zero 94 ! before fields ! now fields 95 sshb (:,:) = 0.e0 ; sshn (:,:) = 0.e0 96 ub (:,:,:) = 0.e0 ; un (:,:,:) = 0.e0 97 vb (:,:,:) = 0.e0 ; vn (:,:,:) = 0.e0 98 rotb (:,:,:) = 0.e0 ; rotn (:,:,:) = 0.e0 99 hdivb(:,:,:) = 0.e0 ; hdivn(:,:,:) = 0.e0 96 ! before fields ! now fields 97 sshb (:,:) = 0._wp ; sshn (:,:) = 0._wp 98 ub (:,:,:) = 0._wp ; un (:,:,:) = 0._wp 99 vb (:,:,:) = 0._wp ; vn (:,:,:) = 0._wp 100 rotb (:,:,:) = 0._wp ; rotn (:,:,:) = 0._wp 101 hdivb(:,:,:) = 0._wp ; hdivn(:,:,:) = 0._wp 102 ! 103 ! ! define e3u_b, e3v_b from e3t_b initialized in domzgr 104 CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 100 105 ! 101 106 IF( cp_cfg == 'eel' ) THEN -
branches/2011/dev_r2802_NOCS_vvlfix/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r2779 r2807 92 92 !! un,vn now horizontal velocity of next time-step 93 93 !!---------------------------------------------------------------------- 94 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released95 94 USE oce , ONLY: ze3u_f => ta , ze3v_f => sa ! (ta,sa) used as 3D workspace 96 USE wrk_nemo, ONLY: zs_t => wrk_2d_1 , zs_u_1 => wrk_2d_2 , zs_v_1 => wrk_2d_397 95 ! 98 96 INTEGER, INTENT( in ) :: kt ! ocean time-step index 99 97 ! 100 98 INTEGER :: ji, jj, jk ! dummy loop indices 99 INTEGER :: iku, ikv ! local integers 101 100 #if ! defined key_dynspg_flt 102 101 REAL(wp) :: z2dt ! temporary scalar 103 102 #endif 104 REAL(wp) :: zue3a, zue3n, zue3b, zuf ! local scalars 105 REAL(wp) :: zve3a, zve3n, zve3b, zvf ! - - 106 REAL(wp) :: zec, zv_t_ij, zv_t_ip1j, zv_t_ijp1 103 REAL(wp) :: zue3a, zue3n, zue3b, zuf, zec ! local scalars 104 REAL(wp) :: zve3a, zve3n, zve3b, zvf ! - - 107 105 !!---------------------------------------------------------------------- 108 109 IF( wrk_in_use(2, 1,2,3) ) THEN110 CALL ctl_stop('dyn_nxt: requested workspace arrays unavailable') ; RETURN111 ENDIF112 106 113 107 IF( kt == nit000 ) THEN … … 238 232 ELSE ! Variable volume ! 239 233 ! ! ================! 240 ! Before scale factor at t-points 241 ! ------------------------------- 242 DO jk = 1, jpkm1 234 ! 235 DO jk = 1, jpkm1 ! Before scale factor at t-points 243 236 fse3t_b(:,:,jk) = fse3t_n(:,:,jk) & 244 237 & + atfp * ( fse3t_b(:,:,jk) + fse3t_a(:,:,jk) & 245 & - 2.e0 * fse3t_n(:,:,jk) ) 246 ENDDO 247 ! Add volume filter correction only at the first level of t-point scale factors 248 zec = atfp * rdt / rau0 238 & - 2._wp * fse3t_n(:,:,jk) ) 239 END DO 240 zec = atfp * rdt / rau0 ! Add filter correction only at the 1st level of t-point scale factors 249 241 fse3t_b(:,:,1) = fse3t_b(:,:,1) - zec * ( emp_b(:,:) - emp(:,:) ) * tmask(:,:,1) 250 ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations251 zs_t (:,:) = e1t(:,:) * e2t(:,:)252 zs_u_1(:,:) = 0.5 / ( e1u(:,:) * e2u(:,:) )253 zs_v_1(:,:) = 0.5 / ( e1v(:,:) * e2v(:,:) )254 242 ! 255 IF( ln_dynadv_vec ) THEN 256 ! Before scale factor at (u/v)-points 257 ! ----------------------------------- 258 ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 259 DO jk = 1, jpkm1 260 DO jj = 1, jpjm1 261 DO ji = 1, jpim1 262 zv_t_ij = zs_t(ji ,jj ) * fse3t_b(ji ,jj ,jk) 263 zv_t_ip1j = zs_t(ji+1,jj ) * fse3t_b(ji+1,jj ,jk) 264 zv_t_ijp1 = zs_t(ji ,jj+1) * fse3t_b(ji ,jj+1,jk) 265 fse3u_b(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 266 fse3v_b(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 267 END DO 268 END DO 269 END DO 270 ! lateral boundary conditions 271 CALL lbc_lnk( fse3u_b(:,:,:), 'U', 1. ) 272 CALL lbc_lnk( fse3v_b(:,:,:), 'V', 1. ) 273 ! Add initial scale factor to scale factor anomaly 274 fse3u_b(:,:,:) = fse3u_b(:,:,:) + fse3u_0(:,:,:) 275 fse3v_b(:,:,:) = fse3v_b(:,:,:) + fse3v_0(:,:,:) 276 ! Leap-Frog - Asselin filter and swap: applied on velocity 277 ! ----------------------------------- 278 DO jk = 1, jpkm1 279 DO jj = 1, jpj 243 IF( ln_dynadv_vec ) THEN ! vector invariant form (no thickness weighted calulation) 244 ! 245 ! ! before scale factors at u- & v-pts (computed from fse3t_b) 246 CALL dom_vvl_2( kt, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 247 ! 248 DO jk = 1, jpkm1 ! Leap-Frog - Asselin filter and swap: applied on velocity 249 DO jj = 1, jpj ! -------- 280 250 DO ji = 1, jpi 281 251 zuf = un(ji,jj,jk) + atfp * ( ub(ji,jj,jk) - 2.e0 * un(ji,jj,jk) + ua(ji,jj,jk) ) … … 290 260 END DO 291 261 ! 292 ELSE 293 ! Temporary filered scale factor at (u/v)-points (will become before scale factor) 294 !----------------------------------------------- 295 ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 296 DO jk = 1, jpkm1 297 DO jj = 1, jpjm1 298 DO ji = 1, jpim1 299 zv_t_ij = zs_t(ji ,jj ) * fse3t_b(ji ,jj ,jk) 300 zv_t_ip1j = zs_t(ji+1,jj ) * fse3t_b(ji+1,jj ,jk) 301 zv_t_ijp1 = zs_t(ji ,jj+1) * fse3t_b(ji ,jj+1,jk) 302 ze3u_f(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 303 ze3v_f(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 304 END DO 305 END DO 306 END DO 307 ! lateral boundary conditions 308 CALL lbc_lnk( ze3u_f, 'U', 1. ) 309 CALL lbc_lnk( ze3v_f, 'V', 1. ) 310 ! Add initial scale factor to scale factor anomaly 311 ze3u_f(:,:,:) = ze3u_f(:,:,:) + fse3u_0(:,:,:) 312 ze3v_f(:,:,:) = ze3v_f(:,:,:) + fse3v_0(:,:,:) 313 ! Leap-Frog - Asselin filter and swap: applied on thickness weighted velocity 314 ! ----------------------------------- =========================== 315 DO jk = 1, jpkm1 316 DO jj = 1, jpj 317 DO ji = 1, jpim1 262 ELSE ! flux form (thickness weighted calulation) 263 ! 264 CALL dom_vvl_2( kt, ze3u_f, ze3v_f ) ! before scale factors at u- & v-pts (computed from fse3t_b) 265 ! 266 DO jk = 1, jpkm1 ! Leap-Frog - Asselin filter and swap: 267 DO jj = 1, jpj ! applied on thickness weighted velocity 268 DO ji = 1, jpim1 ! --------------------------- 318 269 zue3a = ua(ji,jj,jk) * fse3u_a(ji,jj,jk) 319 270 zve3a = va(ji,jj,jk) * fse3v_a(ji,jj,jk) … … 323 274 zve3b = vb(ji,jj,jk) * fse3v_b(ji,jj,jk) 324 275 ! 325 zuf = ( zue3n + atfp * ( zue3b - 2.e0* zue3n + zue3a ) ) / ze3u_f(ji,jj,jk)326 zvf = ( zve3n + atfp * ( zve3b - 2.e0* zve3n + zve3a ) ) / ze3v_f(ji,jj,jk)276 zuf = ( zue3n + atfp * ( zue3b - 2._wp * zue3n + zue3a ) ) / ze3u_f(ji,jj,jk) 277 zvf = ( zve3n + atfp * ( zve3b - 2._wp * zve3n + zve3a ) ) / ze3v_f(ji,jj,jk) 327 278 ! 328 ub(ji,jj,jk) = zuf 279 ub(ji,jj,jk) = zuf ! ub <-- filtered velocity 329 280 vb(ji,jj,jk) = zvf 330 un(ji,jj,jk) = ua(ji,jj,jk) 281 un(ji,jj,jk) = ua(ji,jj,jk) ! un <-- ua 331 282 vn(ji,jj,jk) = va(ji,jj,jk) 332 283 END DO 333 284 END DO 334 285 END DO 335 fse3u_b(:,:, :) = ze3u_f(:,:,:)! e3u_b <-- filtered scale factor336 fse3v_b(:,:, :) = ze3v_f(:,:,:)337 CALL lbc_lnk( ub, 'U', -1. ) 286 fse3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor 287 fse3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) 288 CALL lbc_lnk( ub, 'U', -1. ) ! lateral boundary conditions 338 289 CALL lbc_lnk( vb, 'V', -1. ) 339 290 ENDIF … … 346 297 & tab3d_2=vn, clinfo2=' Vn: ' , mask2=vmask ) 347 298 ! 348 IF( wrk_not_released(2, 1,2,3) ) CALL ctl_stop('dyn_nxt: failed to release workspace arrays')349 !350 299 END SUBROUTINE dyn_nxt 351 300
Note: See TracChangeset
for help on using the changeset viewer.
