Changeset 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/TRA/trazdf_imp.F90
- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/TRA/trazdf_imp.F90
r3294 r6225 16 16 !! 3.3 ! 2010-06 (C. Ethe, G. Madec) Merge TRA-TRC 17 17 !! - ! 2011-02 (A. Coward, C. Ethe, G. Madec) improvment of surface boundary condition 18 !! 3.7 ! 2015-11 (G. Madec, A. Coward) non linear free surface by default 18 19 !!---------------------------------------------------------------------- 19 20 20 21 !!---------------------------------------------------------------------- 21 !! tra_zdf_imp : Update the tracer trend with the diagonal vertical 22 !! part of the mixing tensor. 23 !!---------------------------------------------------------------------- 24 USE oce ! ocean dynamics and tracers variables 25 USE dom_oce ! ocean space and time domain variables 26 USE zdf_oce ! ocean vertical physics variables 27 USE trc_oce ! share passive tracers/ocean variables 28 USE domvvl ! variable volume 29 USE ldftra_oce ! ocean active tracers: lateral physics 30 USE ldftra ! lateral mixing type 31 USE ldfslp ! lateral physics: slope of diffusion 32 USE zdfddm ! ocean vertical physics: double diffusion 33 USE traldf_iso_grif ! active tracers: Griffies operator 34 USE in_out_manager ! I/O manager 35 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 36 USE lib_mpp ! MPP library 37 USE wrk_nemo ! Memory Allocation 38 USE timing ! Timing 22 !! tra_zdf_imp : Update the tracer trend with vertical mixing, nad compute the after tracer field 23 !!---------------------------------------------------------------------- 24 USE oce ! ocean dynamics and tracers variables 25 USE dom_oce ! ocean space and time domain variables 26 USE zdf_oce ! ocean vertical physics variables 27 USE trc_oce ! share passive tracers/ocean variables 28 USE domvvl ! variable volume 29 USE ldftra ! lateral mixing type 30 USE ldfslp ! lateral physics: slope of diffusion 31 USE zdfddm ! ocean vertical physics: double diffusion 32 USE traldf_triad ! active tracers: Method of Stabilizing Correction 33 ! 34 USE in_out_manager ! I/O manager 35 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 36 USE lib_mpp ! MPP library 37 USE wrk_nemo ! Memory Allocation 38 USE timing ! Timing 39 39 40 40 IMPLICIT NONE … … 43 43 PUBLIC tra_zdf_imp ! routine called by step.F90 44 44 45 REAL(wp) :: r_vvl ! variable volume indicator, =1 if lk_vvl=T, =0 otherwise46 47 45 !! * Substitutions 48 # include "domzgr_substitute.h90"49 # include "ldftra_substitute.h90"50 46 # include "zdfddm_substitute.h90" 51 47 # include "vectopt_loop_substitute.h90" 52 48 !!---------------------------------------------------------------------- 53 !! NEMO/OPA 3. 3 , NEMO Consortium (2010)49 !! NEMO/OPA 3.7 , NEMO Consortium (2015) 54 50 !! $Id$ 55 51 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 66 62 !! it is already computed and add to the general trend in traldf) 67 63 !! 68 !! ** Method : The vertical diffusion of the tracer t is given by: 69 !! difft = dz( avt dz(t) ) = 1/e3t dk+1( avt/e3w dk(t) ) 70 !! It is computed using a backward time scheme (t=ta). 64 !! ** Method : The vertical diffusion of a tracer ,t , is given by: 65 !! difft = dz( avt dz(t) ) = 1/e3t dk+1( avt/e3w dk(t) ) 66 !! It is computed using a backward time scheme (t=after field) 67 !! which provide directly the after tracer field. 71 68 !! If lk_zdfddm=T, use avs for salinity or for passive tracers 72 69 !! Surface and bottom boundary conditions: no diffusive flux on … … 76 73 !! ** Action : - pta becomes the after tracer 77 74 !!--------------------------------------------------------------------- 78 USE oce , ONLY: zwd => ua , zws => va ! (ua,va) used as 3D workspace79 !80 75 INTEGER , INTENT(in ) :: kt ! ocean time-step index 81 INTEGER , INTENT(in ) :: kit000 76 INTEGER , INTENT(in ) :: kit000 ! first time step index 82 77 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 83 78 INTEGER , INTENT(in ) :: kjpt ! number of tracers 84 REAL(wp) , DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile oftracer time-step79 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 85 80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 86 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! in: tracer trend ; out: after tracer field 87 82 ! 88 83 INTEGER :: ji, jj, jk, jn ! dummy loop indices 89 REAL(wp) :: zrhs , ze3tb, ze3tn, ze3ta! local scalars90 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwt 84 REAL(wp) :: zrhs ! local scalars 85 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwi, zwt, zwd, zws 91 86 !!--------------------------------------------------------------------- 92 87 ! 93 88 IF( nn_timing == 1 ) CALL timing_start('tra_zdf_imp') 94 89 ! 95 CALL wrk_alloc( jpi, jpj, jpk, zwi, zwt)90 CALL wrk_alloc( jpi,jpj,jpk, zwi, zwt, zwd, zws ) 96 91 ! 97 92 IF( kt == kit000 ) THEN … … 99 94 IF(lwp)WRITE(numout,*) 'tra_zdf_imp : implicit vertical mixing on ', cdtype 100 95 IF(lwp)WRITE(numout,*) '~~~~~~~~~~~ ' 101 !102 IF( lk_vvl ) THEN ; r_vvl = 1._wp ! Variable volume indicator103 ELSE ; r_vvl = 0._wp104 ENDIF105 96 ENDIF 106 !107 97 ! ! ============= ! 108 98 DO jn = 1, kjpt ! tracer loop ! 109 99 ! ! ============= ! 110 !111 100 ! Matrix construction 112 101 ! -------------------- … … 122 111 zwt(:,:,1) = 0._wp 123 112 ! 124 #if defined key_ldfslp 125 ! isoneutral diffusion: add the contribution126 IF( ln_traldf_grif ) THEN ! Griffies isoneutral diff127 DO jk = 2, jpkm1128 DO jj = 2, jpjm1129 DO ji = fs_2, fs_jpim1 ! vector opt.130 zwt(ji,jj,jk) = zwt(ji,jj,jk) + ah_wslp2(ji,jj,jk)113 IF( l_ldfslp ) THEN ! isoneutral diffusion: add the contribution 114 IF( ln_traldf_msc ) THEN ! MSC iso-neutral operator 115 DO jk = 2, jpkm1 116 DO jj = 2, jpjm1 117 DO ji = fs_2, fs_jpim1 ! vector opt. 118 zwt(ji,jj,jk) = zwt(ji,jj,jk) + akz(ji,jj,jk) 119 END DO 131 120 END DO 132 121 END DO 133 END DO 134 ELSE IF( l_traldf_rot ) THEN ! standard isoneutral diff 135 DO jk = 2, jpkm1 136 DO jj = 2, jpjm1 137 DO ji = fs_2, fs_jpim1 ! vector opt. 138 zwt(ji,jj,jk) = zwt(ji,jj,jk) + fsahtw(ji,jj,jk) & 139 & * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk) & 140 & + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) ) 122 ELSE ! standard or triad iso-neutral operator 123 DO jk = 2, jpkm1 124 DO jj = 2, jpjm1 125 DO ji = fs_2, fs_jpim1 ! vector opt. 126 zwt(ji,jj,jk) = zwt(ji,jj,jk) + ah_wslp2(ji,jj,jk) 127 END DO 141 128 END DO 142 129 END DO 143 END DO130 ENDIF 144 131 ENDIF 145 #endif 132 ! 146 133 ! Diagonal, lower (i), upper (s) (including the bottom boundary condition since avt is masked) 147 134 DO jk = 1, jpkm1 148 135 DO jj = 2, jpjm1 149 136 DO ji = fs_2, fs_jpim1 ! vector opt. 150 ze3ta = ( 1. - r_vvl ) + r_vvl * fse3t_a(ji,jj,jk) ! after scale factor at T-point 151 ze3tn = r_vvl + ( 1. - r_vvl ) * fse3t_n(ji,jj,jk) ! now scale factor at T-point 152 zwi(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk ) / ( ze3tn * fse3w(ji,jj,jk ) ) 153 zws(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk+1) / ( ze3tn * fse3w(ji,jj,jk+1) ) 154 zwd(ji,jj,jk) = ze3ta - zwi(ji,jj,jk) - zws(ji,jj,jk) 137 !!gm BUG I think, use e3w_a instead of e3w_n 138 zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk ) / e3w_n(ji,jj,jk ) 139 zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w_n(ji,jj,jk+1) 140 zwd(ji,jj,jk) = e3t_a(ji,jj,jk) - zwi(ji,jj,jk) - zws(ji,jj,jk) 155 141 END DO 156 142 END DO … … 176 162 ! used as a work space array: its value is modified. 177 163 ! 178 ! first recurrence: Tk = Dk - Ik Sk-1 / Tk-1 (increasing k) 179 ! done once for all passive tracers (so included in the IF instruction) 180 DO jj = 2, jpjm1 181 DO ji = fs_2, fs_jpim1 164 DO jj = 2, jpjm1 !* 1st recurrence: Tk = Dk - Ik Sk-1 / Tk-1 (increasing k) 165 DO ji = fs_2, fs_jpim1 ! done one for all passive tracers (so included in the IF instruction) 182 166 zwt(ji,jj,1) = zwd(ji,jj,1) 183 167 END DO … … 186 170 DO jj = 2, jpjm1 187 171 DO ji = fs_2, fs_jpim1 188 zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwt(ji,jj,jk-1)172 zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwt(ji,jj,jk-1) 189 173 END DO 190 174 END DO 191 175 END DO 192 176 ! 193 END 177 ENDIF 194 178 ! 195 ! second recurrence: Zk = Yk - Ik / Tk-1 Zk-1 196 DO jj = 2, jpjm1 179 DO jj = 2, jpjm1 !* 2nd recurrence: Zk = Yk - Ik / Tk-1 Zk-1 197 180 DO ji = fs_2, fs_jpim1 198 ze3tb = ( 1. - r_vvl ) + r_vvl * fse3t_b(ji,jj,1) 199 ze3tn = ( 1. - r_vvl ) + r_vvl * fse3t(ji,jj,1) 200 pta(ji,jj,1,jn) = ze3tb * ptb(ji,jj,1,jn) + p2dt(1) * ze3tn * pta(ji,jj,1,jn) 181 pta(ji,jj,1,jn) = e3t_b(ji,jj,1) * ptb(ji,jj,1,jn) + p2dt * e3t_n(ji,jj,1) * pta(ji,jj,1,jn) 201 182 END DO 202 183 END DO … … 204 185 DO jj = 2, jpjm1 205 186 DO ji = fs_2, fs_jpim1 206 ze3tb = ( 1. - r_vvl ) + r_vvl * fse3t_b(ji,jj,jk) 207 ze3tn = ( 1. - r_vvl ) + r_vvl * fse3t (ji,jj,jk) 208 zrhs = ze3tb * ptb(ji,jj,jk,jn) + p2dt(jk) * ze3tn * pta(ji,jj,jk,jn) ! zrhs=right hand side 187 zrhs = e3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * e3t_n(ji,jj,jk) * pta(ji,jj,jk,jn) ! zrhs=right hand side 209 188 pta(ji,jj,jk,jn) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * pta(ji,jj,jk-1,jn) 210 189 END DO 211 190 END DO 212 191 END DO 213 214 ! third recurrence: Xk = (Zk - Sk Xk+1 ) / Tk (result is the after tracer) 215 DO jj = 2, jpjm1 192 ! 193 DO jj = 2, jpjm1 !* 3d recurrence: Xk = (Zk - Sk Xk+1 ) / Tk (result is the after tracer) 216 194 DO ji = fs_2, fs_jpim1 217 195 pta(ji,jj,jpkm1,jn) = pta(ji,jj,jpkm1,jn) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1) … … 230 208 ! ! ================= ! 231 209 ! 232 CALL wrk_dealloc( jpi, jpj, jpk, zwi, zwt)210 CALL wrk_dealloc( jpi,jpj,jpk, zwi, zwt, zwd, zws ) 233 211 ! 234 212 IF( nn_timing == 1 ) CALL timing_stop('tra_zdf_imp')
Note: See TracChangeset
for help on using the changeset viewer.
