- Timestamp:
- 2016-11-28T17:04:10+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/TRA/trazdf_exp.F90
r3294 r7351 20 20 21 21 !!---------------------------------------------------------------------- 22 !! tra_zdf_exp : compute the tracer the vertical diffusion trend using a23 !! split-explicit time stepping and provide the after tracer22 !! tra_zdf_exp : compute the tracer the vertical diffusion trend using a 23 !! split-explicit time stepping and provide the after tracer 24 24 !!---------------------------------------------------------------------- 25 USE oce ! ocean dynamics and active tracers 26 USE dom_oce ! ocean space and time domain 27 USE domvvl ! variable volume levels 28 USE zdf_oce ! ocean vertical physics 29 USE zdfddm ! ocean vertical physics: double diffusion 30 USE trc_oce ! share passive tracers/Ocean variables 31 USE in_out_manager ! I/O manager 32 USE lib_mpp ! MPP library 33 USE wrk_nemo ! Memory Allocation 34 USE timing ! Timing 25 USE oce ! ocean dynamics and active tracers 26 USE dom_oce ! ocean space and time domain 27 USE domvvl ! variable volume levels 28 USE zdf_oce ! ocean vertical physics 29 USE zdfddm ! ocean vertical physics: double diffusion 30 USE trc_oce ! share passive tracers/Ocean variables 31 ! 32 USE in_out_manager ! I/O manager 33 USE lib_mpp ! MPP library 34 USE wrk_nemo ! Memory Allocation 35 USE timing ! Timing 35 36 36 37 IMPLICIT NONE … … 40 41 41 42 !! * Substitutions 42 # include "domzgr_substitute.h90"43 43 # include "zdfddm_substitute.h90" 44 44 # include "vectopt_loop_substitute.h90" … … 50 50 CONTAINS 51 51 52 SUBROUTINE tra_zdf_exp( kt, kit000, cdtype, p2dt, k n_zdfexp, &53 & ptb , pta, kjpt )52 SUBROUTINE tra_zdf_exp( kt, kit000, cdtype, p2dt, ksts, & 53 & ptb , pta , kjpt ) 54 54 !!---------------------------------------------------------------------- 55 55 !! *** ROUTINE tra_zdf_exp *** … … 60 60 !! ** Method : - The after tracer fields due to the vertical diffusion 61 61 !! of tracers alone is given by: 62 !! z wx= ptb + p2dt difft62 !! ztb = ptb + p2dt difft 63 63 !! where difft = dz( avt dz(ptb) ) = 1/e3t dk+1( avt/e3w dk(ptb) ) 64 64 !! (if lk_zdfddm=T use avs on salinity and passive tracers instead of avt) … … 67 67 !! (N.B. bottom condition is applied through the masked field avt). 68 68 !! - the after tracer fields due to the whole trend is 69 !! obtained in leap-frog environment by : 70 !! pta = zwx + p2dt pta 71 !! - in case of variable level thickness (lk_vvl=T) the 72 !! the leap-frog is applied on thickness weighted tracer. That is: 73 !! pta = [ ptb*e3tb + e3tn*( zwx - ptb + p2dt pta ) ] / e3tn 69 !! obtained in leap-frog environment applied on thickness weighted tracer by : 70 !! pta = [ ptb*e3tb + e3tn*( ztb - ptb + p2dt pta ) ] / e3tn 74 71 !! 75 72 !! ** Action : - after tracer fields pta 76 73 !!--------------------------------------------------------------------- 74 INTEGER , INTENT(in ) :: kt ! ocean time-step index 75 INTEGER , INTENT(in ) :: kit000 ! first time step index 76 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 77 INTEGER , INTENT(in ) :: kjpt ! number of tracers 78 INTEGER , INTENT(in ) :: ksts ! number of sub-time step 79 REAL(wp) , INTENT(in ) :: p2dt ! vertical profile of tracer time-step 80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! in: tracer trend ; out: after tracer field 77 82 ! 78 INTEGER , INTENT(in ) :: kt ! ocean time-step index 79 INTEGER , INTENT(in ) :: kit000 ! first time step index 80 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 81 INTEGER , INTENT(in ) :: kjpt ! number of tracers 82 INTEGER , INTENT(in ) :: kn_zdfexp ! number of sub-time step 83 REAL(wp), DIMENSION( jpk ), INTENT(in ) :: p2dt ! vertical profile of tracer time-step 84 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 85 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 86 ! 87 INTEGER :: ji, jj, jk, jn, jl ! dummy loop indices 88 REAL(wp) :: zlavmr, zave3r, ze3tr ! local scalars 89 REAL(wp) :: ztra, ze3tb ! - - 90 REAL(wp), POINTER, DIMENSION(:,:,:) :: zwx, zwy 83 INTEGER :: ji, jj, jk, jn, jl ! dummy loop indices 84 REAL(wp) :: z1_ksts, ze3tr ! local scalars 85 REAL(wp) :: ztra, ze3tb ! - - 86 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztb, zwf 91 87 !!--------------------------------------------------------------------- 92 88 ! 93 89 IF( nn_timing == 1 ) CALL timing_start('tra_zdf_exp') 94 90 ! 95 CALL wrk_alloc( jpi, jpj, jpk, zwx, zwy)91 CALL wrk_alloc( jpi,jpj,jpk, ztb, zwf ) 96 92 ! 97 98 93 IF( kt == kit000 ) THEN 99 94 IF(lwp) WRITE(numout,*) … … 104 99 ! Initializations 105 100 ! --------------- 106 zlavmr = 1. / float( kn_zdfexp ) ! Local constant 101 z1_ksts = 1._wp / REAL( ksts, wp ) 102 zwf(:,:, 1 ) = 0._wp ! no flux at the surface and at bottom level 103 zwf(:,:,jpk) = 0._wp 107 104 ! 108 105 ! 109 DO jn = 1, kjpt ! loop over tracers106 DO jn = 1, kjpt !== loop over tracers ==! 110 107 ! 111 zwy(:,:, 1 ) = 0.e0 ! surface boundary conditions: no flux 112 zwy(:,:,jpk) = 0.e0 ! bottom boundary conditions: no flux 113 ! 114 zwx(:,:,:) = ptb(:,:,:,jn) ! zwx array set to before tracer values 115 116 ! Split-explicit loop (after tracer due to the vertical diffusion alone) 117 ! ------------------- 118 ! 119 DO jl = 1, kn_zdfexp 120 ! ! first vertical derivative 121 DO jk = 2, jpk 108 ztb(:,:,:) = ptb(:,:,:,jn) ! initial before value for tracer 109 ! 110 DO jl = 1, ksts !== Split-explicit loop ==! 111 ! 112 DO jk = 2, jpk ! 1st vertical derivative (w-flux) 122 113 DO jj = 2, jpjm1 123 114 DO ji = fs_2, fs_jpim1 ! vector opt. 124 zave3r = 1.e0 / fse3w_n(ji,jj,jk)125 115 IF( cdtype == 'TRA' .AND. jn == jp_tem ) THEN ! temperature : use of avt 126 zw y(ji,jj,jk) = avt(ji,jj,jk) * ( zwx(ji,jj,jk-1) - zwx(ji,jj,jk) ) * zave3r116 zwf(ji,jj,jk) = avt(ji,jj,jk) * ( ztb(ji,jj,jk-1) - ztb(ji,jj,jk) ) / e3w_b(ji,jj,jk) 127 117 ELSE ! salinity or pass. tracer : use of avs 128 zw y(ji,jj,jk) = fsavs(ji,jj,jk) * ( zwx(ji,jj,jk-1) - zwx(ji,jj,jk) ) * zave3r118 zwf(ji,jj,jk) = fsavs(ji,jj,jk) * ( ztb(ji,jj,jk-1) - ztb(ji,jj,jk) ) / e3w_b(ji,jj,jk) 129 119 END IF 130 120 END DO … … 132 122 END DO 133 123 ! 134 DO jk = 1, jpkm1 ! second vertical derivative ==> tracer at kt+l*2*rdt/nn_zdfexp124 DO jk = 1, jpkm1 ! 2nd vertical derivative ==> tracer at kt+l*2*rdt/nn_zdfexp 135 125 DO jj = 2, jpjm1 136 126 DO ji = fs_2, fs_jpim1 ! vector opt. 137 ze3tr = zlavmr / fse3t_n(ji,jj,jk) 138 zwx(ji,jj,jk) = zwx(ji,jj,jk) + p2dt(jk) * ( zwy(ji,jj,jk) - zwy(ji,jj,jk+1) ) * ze3tr 127 ztb(ji,jj,jk) = ztb(ji,jj,jk) + p2dt * ( zwf(ji,jj,jk) - zwf(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 139 128 END DO 140 129 END DO 141 130 END DO 142 131 ! 143 END DO 132 END DO ! end sub-time stepping 144 133 145 ! After tracer due to all trends 146 ! ------------------------------ 147 IF( lk_vvl ) THEN ! variable level thickness : leap-frog on tracer*e3t 148 DO jk = 1, jpkm1 149 DO jj = 2, jpjm1 150 DO ji = fs_2, fs_jpim1 ! vector opt. 151 ze3tb = fse3t_b(ji,jj,jk) / fse3t(ji,jj,jk) ! before e3t 152 ztra = zwx(ji,jj,jk) - ptb(ji,jj,jk,jn) + p2dt(jk) * pta(ji,jj,jk,jn) ! total trends * 2*rdt 153 pta(ji,jj,jk,jn) = ( ze3tb * ptb(ji,jj,jk,jn) + ztra ) * tmask(ji,jj,jk) 154 END DO 134 DO jk = 1, jpkm1 !== After tracer due to all trends 135 DO jj = 2, jpjm1 136 DO ji = fs_2, fs_jpim1 ! vector opt. 137 ze3tb = e3t_b(ji,jj,jk) / e3t_n(ji,jj,jk) 138 ztra = ( ztb(ji,jj,jk) - ptb(ji,jj,jk,jn) ) + p2dt * pta(ji,jj,jk,jn) ! total trend * 2dt 139 pta(ji,jj,jk,jn) = ( ze3tb * ptb(ji,jj,jk,jn) + ztra ) * tmask(ji,jj,jk) ! after tracer 155 140 END DO 156 141 END DO 157 ELSE ! fixed level thickness : leap-frog on tracers 158 DO jk = 1, jpkm1 159 DO jj = 2, jpjm1 160 DO ji = fs_2, fs_jpim1 ! vector opt. 161 pta(ji,jj,jk,jn) = ( zwx(ji,jj,jk) + p2dt(jk) * pta(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 162 END DO 163 END DO 164 END DO 165 ENDIF 142 END DO 166 143 ! 167 END DO 144 END DO ! end of tracer loop 168 145 ! 169 CALL wrk_dealloc( jpi, jpj, jpk, zwx, zwy)146 CALL wrk_dealloc( jpi,jpj,jpk, ztb, zwf ) 170 147 ! 171 148 IF( nn_timing == 1 ) CALL timing_stop('tra_zdf_exp')
Note: See TracChangeset
for help on using the changeset viewer.
