MODULE trazdf_imp_jki !!============================================================================== !! *** MODULE trazdf_imp_jki *** !! Ocean active tracers: vertical component of the tracer mixing trend !!============================================================================== !! History : !! 7.0 ! 91-11 (G. Madec) Original code !! ! 92-06 (M. Imbard) correction on tracer trend loops !! ! 96-01 (G. Madec) statement function for e3 !! ! 97-05 (G. Madec) vertical component of isopycnal !! ! 97-07 (G. Madec) geopotential diffusion in s-coord !! ! 00-08 (G. Madec) double diffusive mixing !! 8.5 ! 02-08 (G. Madec) F90: Free form and module !! 9.0 ! 04-08 (C. Talandier) New trends organization !! ! 05-11 (G. Madec) New step reorganisation !!---------------------------------------------------------------------- !! tra_zdf_imp_jki : Update the tracer trend with the diagonal vertical !! part of the mixing tensor. !! use j-k-i loops. !!---------------------------------------------------------------------- !! * Modules used USE oce ! ocean dynamics and tracers variables USE dom_oce ! ocean space and time domain variables USE ldfslp ! Make iso-neutral slopes available USE ldftra_oce ! ocean active tracers: lateral physics USE zdf_oce ! ocean vertical physics USE zdfddm ! ocean vertical physics: double diffusion USE trdmod ! ocean active tracers trends USE trdmod_oce ! ocean variables trends USE in_out_manager ! I/O manager USE lbclnk ! ocean lateral boundary conditions (or mpp link) USE prtctl ! Print control IMPLICIT NONE PRIVATE !! * Accessibility PUBLIC tra_zdf_imp_jki ! routine called by step.F90 !! * Substitutions # include "domzgr_substitute.h90" # include "ldftra_substitute.h90" # include "zdfddm_substitute.h90" !!---------------------------------------------------------------------- !! OPA 9.0 , LOCEAN-IPSL (2005) !!---------------------------------------------------------------------- CONTAINS SUBROUTINE tra_zdf_imp_jki( kt, p2dt ) !!---------------------------------------------------------------------- !! *** ROUTINE tra_zdf_imp_jki *** !! !! ** Purpose : !! Compute the trend due to the vertical tracer diffusion inclu- !! ding the vertical component of lateral mixing (only for second !! order operator, for fourth order it is already computed and !! add to the general trend in traldf.F) and add it to the general !! trend of the tracer equations. !! !! ** Method : !! save avt coef. resulting from vertical physics alone in zavt: !! zavt = avt !! update and save in zavt the vertical eddy viscosity coefficient: !! avt = avt + wslpi^2+wslj^2 !! !! Second part: vertical trend associated with the vertical physics !! =========== (including the vertical flux proportional to dk[t] !! associated with the lateral mixing, through the !! update of avt) !! The vertical diffusion of tracers (t & s) is given by: !! difft = dz( avt dz(t) ) = 1/e3t dk+1( avt/e3w dk(t) ) !! It is computed using a backward time scheme, t=ta. !! Surface and bottom boundary conditions: no diffusive flux on !! both tracers (bottom, applied through the masked field avt). !! Add this trend to the general trend ta,sa : !! ta = ta + dz( avt dz(t) ) !! (sa = sa + dz( avs dz(t) ) if lk_zdfddm=T ) !! !! Third part: recover avt resulting from the vertical physics !! ========== alone, for further diagnostics (for example to !! compute the turbocline depth in zdfmxl.F90). !! avt = zavt !! (avs = zavs if lk_zdfddm=T ) !! !! ** Action : !! Update (ta,sa) arrays with the before vertical diffusion trend !! !!--------------------------------------------------------------------- !! * Modules used USE oce , & # if defined key_zdfddm zavs => va, & ! use va as workspace # endif zavt => ua ! use ua as workspace !! * Arguments INTEGER, INTENT( in ) :: kt ! ocean time-step index REAL(wp), DIMENSION(jpk), INTENT( in ) :: & p2dt ! vertical profile of tracer time-step !! * Local declarations INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ikst, ikenm2, ikstp1 ! temporary integers REAL(wp) :: zavi, zavip1 ! temporary scalar REAL(wp), DIMENSION(jpi,jpk) :: & zwd, zws, zwi, & ! ??? zwx, zwy, zwz, zwt ! ??? !!--------------------------------------------------------------------- IF( kt == nit000 ) THEN IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'tra_zdf_imp_jki : implicit vertical mixing (j-k-i loops)' IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~' ENDIF !CDIR PARALLEL DO PRIVATE( zwd, zws, zwi, zwx, zwy, zwz, zwt ) !$OMP PARALLEL DO PRIVATE( zwd, zws, zwi, zwx, zwy, zwz, zwt ) ! ! =============== DO jj = 2, jpjm1 ! Vertical slab ! ! =============== ! II. Vertical trend associated with the vertical physics ! ======================================================= ! (including the vertical flux proportional to dk[t] associated ! with the lateral mixing, through the avt update) ! dk[ avt dk[ (t,s) ] ] diffusive trends ! II.0 Matrix construction ! ------------------------ ! Diagonal, inferior, superior ! (including the bottom boundary condition via avt masked) DO jk = 1, jpkm1 DO ji = 2, jpim1 #if defined key_ldfslp zavi = avt(ji,jj,jk ) + fsahtw(ji,jj,jk )*( wslpi(ji,jj,jk )*wslpi(ji,jj,jk ) & & +wslpj(ji,jj,jk )*wslpj(ji,jj,jk ) ) zavip1 = avt(ji,jj,jk+1) + fsahtw(ji,jj,jk+1)*( wslpi(ji,jj,jk+1)*wslpi(ji,jj,jk+1) & & +wslpj(ji,jj,jk+1)*wslpj(ji,jj,jk+1) ) #else zavi = avt(ji,jj,jk ) zavip1 = avt(ji,jj,jk+1) #endif zwi(ji,jk) = - p2dt(jk) * zavi / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk ) ) zws(ji,jk) = - p2dt(jk) * zavip1 / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) ) zwd(ji,jk) = 1. - zwi(ji,jk) - zws(ji,jk) END DO END DO ! Surface boudary conditions DO ji = 2, jpim1 zwi(ji,1) = 0.e0 zwd(ji,1) = 1. - zws(ji,1) END DO ! II.1. Vertical diffusion on t ! --------------------------- ! Second member construction DO jk = 1, jpkm1 DO ji = 2, jpim1 zwy(ji,jk) = tb(ji,jj,jk) + p2dt(jk) * ta(ji,jj,jk) END DO END DO ! Matrix inversion from the first level ikst = 1 # include "zdf.matrixsolver.h90" ! Save the masked temperature after in ta ! (c a u t i o n: temperature not its trend, Leap-frog scheme done ! it will not be done in tranxt) DO jk = 1, jpkm1 DO ji = 2, jpim1 ta(ji,jj,jk) = zwx(ji,jk) * tmask(ji,jj,jk) END DO END DO ! II.2 Vertical diffusion on s ! --------------------------- #if defined key_zdfddm ! Rebuild the Matrix as avt /= avs ! Diagonal, inferior, superior ! (including the bottom boundary condition via avs masked) DO jk = 1, jpkm1 DO ji = 2, jpim1 #if defined key_ldfslp zavi = fsavs(ji,jj,jk ) + fsahtw(ji,jj,jk )*( wslpi(ji,jj,jk )*wslpi(ji,jj,jk ) & & +wslpj(ji,jj,jk )*wslpj(ji,jj,jk ) ) zavip1 = fsavs(ji,jj,jk+1) + fsahtw(ji,jj,jk+1)*( wslpi(ji,jj,jk+1)*wslpi(ji,jj,jk+1) & & +wslpj(ji,jj,jk+1)*wslpj(ji,jj,jk+1) ) #else zavi = fsavs(ji,jj,jk ) zavip1 = fsavs(ji,jj,jk+1) #endif zwi(ji,jk) = - p2dt(jk) * zavi / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk ) ) zws(ji,jk) = - p2dt(jk) * zavip1 / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) ) zwd(ji,jk) = 1. - zwi(ji,jk) - zws(ji,jk) END DO END DO ! Surface boudary conditions DO ji = 2, jpim1 zwi(ji,1) = 0.e0 zwd(ji,1) = 1. - zws(ji,1) END DO #endif ! Second member construction DO jk = 1, jpkm1 DO ji = 2, jpim1 zwy(ji,jk) = sb(ji,jj,jk) + p2dt(jk) * sa(ji,jj,jk) END DO END DO ! Matrix inversion from the first level ikst = 1 # include "zdf.matrixsolver.h90" ! Save the masked salinity after in sa ! (c a u t i o n: salinity not its trend, Leap-frog scheme done ! it will not be done in tranxt) DO jk = 1, jpkm1 DO ji = 2, jpim1 sa(ji,jj,jk) = zwx(ji,jk) * tmask(ji,jj,jk) END DO END DO #if defined key_ldfslp ! III. recover the avt (avs) resulting from vertical physics only ! =============================================================== DO jk = 2, jpkm1 DO ji = 2, jpim1 avt(ji,jj,jk) = zavt(ji,jj,jk) # if defined key_zdfddm fsavs(ji,jj,jk) = zavs(ji,jj,jk) # endif END DO END DO #endif ! ! =============== END DO ! End of slab ! ! =============== END SUBROUTINE tra_zdf_imp_jki !!============================================================================== END MODULE trazdf_imp_jki