MODULE trazdf_iso !!============================================================================== !! *** MODULE trazdf_iso *** !! Ocean active tracers: vertical component of the tracer mixing trend !!============================================================================== #if defined key_ldfslp || defined key_esopa !!---------------------------------------------------------------------- !! 'key_ldfslp' rotation of the lateral mixing tensor !!---------------------------------------------------------------------- !! tra_zdf_iso : update the tracer trend with the vertical part of !! the isopycnal or geopotential s-coord. operator and !! the vertical diffusion !!---------------------------------------------------------------------- !! * 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) IMPLICIT NONE PRIVATE !! * Accessibility PUBLIC tra_zdf_iso ! routine called by step.F90 !! * Substitutions # include "domzgr_substitute.h90" # include "ldftra_substitute.h90" # include "ldfeiv_substitute.h90" # include "zdfddm_substitute.h90" !!---------------------------------------------------------------------- !!---------------------------------------------------------------------- !! OPA 9.0 , LOCEAN-IPSL (2005) !! $Header$ !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt !!---------------------------------------------------------------------- CONTAINS SUBROUTINE tra_zdf_iso( kt ) !!---------------------------------------------------------------------- !! *** ROUTINE tra_zdf_iso *** !! !! ** 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 : !! The vertical component of the lateral diffusive trends is !! provided by a 2nd order operator rotated along neural or geopo- !! tential surfaces to which an eddy induced advection can be added !! It is computed using before fields (forward in time) and isopyc- !! nal or geopotential slopes computed in routine ldfslp. !! !! First part: vertical trends associated with the lateral mixing !! ========== (excluding the vertical flux proportional to dk[t] ) !! vertical fluxes associated with the rotated lateral mixing: !! zftw =-aht { e2t*wslpi di[ mi(mk(tb)) ] !! + e1t*wslpj dj[ mj(mk(tb)) ] } !! 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 !! add vertical Eddy Induced advective fluxes ('lk_traldf_eiv=T): !! zftw = zftw + { di[aht e2u mi(wslpi)] !! +dj[aht e1v mj(wslpj)] } mk(tb) !! take the horizontal divergence of the fluxes: !! difft = 1/(e1t*e2t*e3t) dk[ zftw ] !! Add this trend to the general trend (ta,sa): !! ta = ta + difft !! !! 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 ) !! !! 'key_trdtra' defined: trend saved for futher diagnostics. !! !! macro-tasked on vertical slab (jj-loop) !! !! ** Action : !! Update (ta,sa) arrays with the before vertical diffusion trend !! Save in (ztdta,ztdsa) arrays the trends if 'key_trdtra' defined !! !! 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 !!--------------------------------------------------------------------- !! * 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 !! * Local save REAL(wp), DIMENSION(jpk), SAVE :: & z2dt !! * Local declarations INTEGER :: ji, jj, jk ! dummy loop indices INTEGER :: ikst, ikenm2, ikstp1 ! temporary integers #if defined key_partial_steps INTEGER :: iku, ikv, ikv1 ! temporary integers #endif REAL(wp) :: zta, zsa REAL(wp) :: & zcoef0, zcoef3, & ! ??? zcoef4, zavi, & ! ??? zbtr, zmku, zmkv, & ! ztav, zsav REAL(wp), DIMENSION(jpi,jpk) :: & zwd, zws, zwi, & ! ??? zwx, zwy, zwz, zwt ! ??? REAL(wp), DIMENSION(jpi,jpk) :: & ztfw, zdit, zdjt, zdj1t, & zsfw, zdis, zdjs, zdj1s REAL(wp), DIMENSION(jpi,jpj,jpk) :: & ztavg, zsavg, & ! workspace arrays ztdta, ztdsa ! workspace arrays #if defined key_traldf_eiv || defined key_esopa REAL(wp), DIMENSION(jpi,jpj,jpk) :: & ztfwg, zsfwg REAL(wp) :: & zcoeg3, & zuwk, zvwk, & zuwki, zvwki #endif !!--------------------------------------------------------------------- !! OPA 8.5, LODYC-IPSL (2002) !!--------------------------------------------------------------------- IF( kt == nit000 ) THEN IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'tra_zdf_iso : vertical mixing (including isopycnal component)' IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' #if defined key_diaeiv w_eiv(:,:,:) = 0.e0 #endif ENDIF ! 0. Local constant initialization ! -------------------------------- ztavg(:,:,:) = 0.e0 zsavg(:,:,:) = 0.e0 ! time step = 2 rdttra ex IF( neuler == 0 .AND. kt == nit000 ) THEN z2dt(:) = rdttra(:) ! restarting with Euler time stepping ELSEIF( kt <= nit000 + 1) THEN z2dt(:) = 2. * rdttra(:) ! leapfrog ENDIF ! Save ta and sa trends IF( l_trdtra ) THEN ztdta(:,:,:) = ta(:,:,:) ztdsa(:,:,:) = sa(:,:,:) ENDIF ! ! =============== DO jj = 2, jpjm1 ! Vertical slab ! ! =============== ! I. vertical trends associated with the lateral mixing ! ===================================================== ! (excluding the vertical flux proportional to dk[t] ! I.1 horizontal tracer gradient ! ------------------------------ DO jk = 1, jpkm1 DO ji = 1, jpim1 ! i-gradient of T and S at jj zdit (ji,jk) = ( tb(ji+1,jj,jk)-tb(ji,jj,jk) ) * umask(ji,jj,jk) zdis (ji,jk) = ( sb(ji+1,jj,jk)-sb(ji,jj,jk) ) * umask(ji,jj,jk) ! j-gradient of T and S at jj zdjt (ji,jk) = ( tb(ji,jj+1,jk)-tb(ji,jj,jk) ) * vmask(ji,jj,jk) zdjs (ji,jk) = ( sb(ji,jj+1,jk)-sb(ji,jj,jk) ) * vmask(ji,jj,jk) ! j-gradient of T and S at jj+1 zdj1t(ji,jk) = ( tb(ji,jj,jk)-tb(ji,jj-1,jk) ) * vmask(ji,jj-1,jk) zdj1s(ji,jk) = ( sb(ji,jj,jk)-sb(ji,jj-1,jk) ) * vmask(ji,jj-1,jk) END DO END DO # if defined key_partial_steps ! partial steps correction at the bottom ocean level DO ji = 1, jpim1 ! last ocean level iku = MIN( mbathy(ji,jj), mbathy(ji+1,jj ) ) - 1 ikv = MIN( mbathy(ji,jj), mbathy(ji ,jj+1) ) - 1 ikv1 = MIN( mbathy(ji,jj), mbathy(ji ,jj-1) ) - 1 ! i-gradient of T and S at jj zdit (ji,iku) = gtu(ji,jj) zdis (ji,iku) = gsu(ji,jj) ! j-gradient of T and S at jj zdjt (ji,ikv) = gtv(ji,jj) zdjs (ji,ikv) = gsv(ji,jj) ! j-gradient of T and S at jj+1 zdj1t(ji,ikv1)= gtv(ji,jj-1) zdj1s(ji,ikv1)= gsv(ji,jj-1) END DO #endif ! I.2 Vertical fluxes ! ------------------- ! Surface and bottom vertical fluxes set to zero ztfw(:, 1 ) = 0.e0 zsfw(:, 1 ) = 0.e0 ztfw(:,jpk) = 0.e0 zsfw(:,jpk) = 0.e0 #if defined key_traldf_eiv ztfwg(:,:, 1 ) = 0.e0 zsfwg(:,:, 1 ) = 0.e0 ztfwg(:,:,jpk) = 0.e0 zsfwg(:,:,jpk) = 0.e0 #endif ! interior (2=