Changeset 1110 for trunk/NEMO/OPA_SRC/TRA/trazdf_exp.F90
- Timestamp:
- 2008-06-13T16:10:35+02:00 (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/TRA/trazdf_exp.F90
r789 r1110 3 3 !! *** MODULE trazdf_exp *** 4 4 !! Ocean active tracers: vertical component of the tracer mixing trend using 5 !! a n explicit time-stepping (time spllitting scheme)5 !! a split-explicit time-stepping 6 6 !!============================================================================== 7 !! History : 8 !! 6.0 ! 90-10 (B. Blanke) Original code9 !! 7.0 ! 91-11 (G. Madec)10 !! ! 92-06 (M. Imbard) correction on tracer trend loops11 !! ! 96-01 (G. Madec) statement function for e312 !! ! 97-05 (G. Madec) vertical component of isopycnal13 !! ! 97-07 (G. Madec) geopotential diffusion in s-coord14 !! ! 00-08 (G. Madec) double diffusive mixing15 !! 8.5 ! 02-08 (G. Madec) F90: Free form and module16 !! 9.0 ! 04-08 (C. Talandier) New trendsorganisation17 !! ! 05-11 (G. Madec) New organisation7 !! History : OPA ! 1990-10 (B. Blanke) Original code 8 !! 7.0 ! 1991-11 (G. Madec) 9 !! ! 1992-06 (M. Imbard) correction on tracer trend loops 10 !! ! 1996-01 (G. Madec) statement function for e3 11 !! ! 1997-05 (G. Madec) vertical component of isopycnal 12 !! ! 1997-07 (G. Madec) geopotential diffusion in s-coord 13 !! ! 2000-08 (G. Madec) double diffusive mixing 14 !! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module 15 !! - ! 2004-08 (C. Talandier) New trends organisation 16 !! - ! 2005-11 (G. Madec) New organisation 17 !! 3.0 ! 2008-04 (G. Madec) leap-frog time stepping done in trazdf 18 18 !!---------------------------------------------------------------------- 19 !! tra_zdf_exp : update the tracer trend with the vertical diffusion 20 !! using an explicit time stepping 19 21 20 !!---------------------------------------------------------------------- 22 !! * Modules used 21 !! tra_zdf_exp : compute the tracer the vertical diffusion trend using a 22 !! split-explicit time stepping and provide the after tracer 23 !!---------------------------------------------------------------------- 23 24 USE oce ! ocean dynamics and active tracers 24 25 USE dom_oce ! ocean space and time domain 26 USE domvvl ! variablevolume levels 25 27 USE trdmod ! ocean active tracers trends 26 28 USE trdmod_oce ! ocean variables trends … … 33 35 PRIVATE 34 36 35 !! * Routine accessibility 36 PUBLIC tra_zdf_exp ! routine called by step.F90 37 PUBLIC tra_zdf_exp ! routine called by step.F90 37 38 38 39 !! * Substitutions 39 40 # include "domzgr_substitute.h90" 40 41 # include "zdfddm_substitute.h90" 42 # include "vectopt_loop_substitute.h90" 41 43 !!---------------------------------------------------------------------- 42 !! OPA 9.0 , LOCEAN-IPSL (2005)43 !! $ Header$44 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt44 !! NEMO/OPA 3.0 , LOCEAN-IPSL (2008) 45 !! $Id:$ 46 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 45 47 !!---------------------------------------------------------------------- 46 48 … … 51 53 !! *** ROUTINE tra_zdf_exp *** 52 54 !! 53 !! ** Purpose : Compute the trend due to the vertical tracer mixing 54 !! using an explicit time stepping and add it to the general trend 55 !! of the tracer equations. 55 !! ** Purpose : Compute the after tracer fields due to the vertical 56 !! tracer mixing alone, and then due to the whole tracer trend. 56 57 !! 57 !! ** Method : The vertical diffusion of tracers (t & s) is given by: 58 !! difft = dz( avt dz(tb) ) = 1/e3t dk+1( avt/e3w dk(tb) ) 59 !! It is evaluated with an Euler scheme, using a time splitting 60 !! technique. 61 !! Surface and bottom boundary conditions: no diffusive flux on 62 !! both tracers (bottom, applied through the masked field avt). 63 !! Add this trend to the general trend ta,sa : 64 !! ta = ta + dz( avt dz(t) ) 65 !! (sa = sa + dz( avs dz(t) ) if lk_zdfddm= T) 58 !! ** Method : - The after tracer fields due to the vertical diffusion 59 !! of tracers alone is given by: 60 !! zwx = tb + p2dt difft 61 !! where difft = dz( avt dz(tb) ) = 1/e3t dk+1( avt/e3w dk(tb) ) 62 !! (if lk_zdfddm=T use avs on salinity instead of avt) 63 !! difft is evaluated with an Euler split-explit scheme using a 64 !! no flux boundary condition at both surface and bottomi boundaries. 65 !! (N.B. bottom condition is applied through the masked field avt). 66 !! - the after tracer fields due to the whole trend is 67 !! obtained in leap-frog environment by : 68 !! ta = zwx + p2dt ta 69 !! - in case of variable level thickness (lk_vvl=T) the 70 !! the leap-frog is applied on thickness weighted tracer. That is: 71 !! ta = [ tb*e3tb + e3tn*( zwx - tb + p2dt ta ) ] / e3tn 66 72 !! 67 !! ** Action : - Update (ta,sa) with the before vertical diffusion trend 73 !! ** Action : - after tracer fields (ta,sa) 74 !!--------------------------------------------------------------------- 75 INTEGER , INTENT(in) :: kt ! ocean time-step index 76 REAL(wp), INTENT(in), DIMENSION(jpk) :: p2dt ! vertical profile of tracer time-step 68 77 !! 69 !!--------------------------------------------------------------------- 70 !! * Arguments 71 INTEGER, INTENT( in ) :: kt ! ocean time-step index 72 REAL(wp), DIMENSION(jpk), INTENT( in ) :: & 73 p2dt ! vertical profile of tracer time-step 74 75 !! * Local declarations 76 INTEGER :: ji, jj, jk, jl ! dummy loop indices 77 REAL(wp) :: & 78 zlavmr, & ! temporary scalars 79 zave3r, ze3tr, & ! " " 80 zta, zsa ! " " 81 REAL(wp), DIMENSION(jpi,jpk) :: & 82 zwx, zwy, zwz, zww 78 INTEGER :: ji, jj, jk, jl ! dummy loop indices 79 REAL(wp) :: zlavmr, zave3r, ze3tr ! temporary scalars 80 REAL(wp) :: zta, zsa, ze3tb, zcoef ! temporary scalars 81 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwx, zwy, zwz, zww ! 3D workspace 83 82 !!--------------------------------------------------------------------- 84 83 … … 89 88 ENDIF 90 89 90 ! Initializations 91 ! --------------- 92 zlavmr = 1. / float( n_zdfexp ) ! Local constant 93 ! 94 zwy(:,:, 1 ) = 0.e0 ; zww(:,:, 1 ) = 0.e0 ! surface boundary conditions: no flux 95 zwy(:,:,jpk) = 0.e0 ; zww(:,:,jpk) = 0.e0 ! bottom boundary conditions: no flux 96 ! 97 zwx(:,:,:) = tb(:,:,:) ; zwz(:,:,:) = sb(:,:,:) ! zwx and zwz arrays set to before tracer values 91 98 92 ! 0. Local constant initialization 93 ! -------------------------------- 94 zlavmr = 1. / float( n_zdfexp ) 95 96 ! ! =============== 97 DO jj = 2, jpjm1 ! Vertical slab 98 ! ! =============== 99 ! 1. Initializations 100 ! ------------------ 101 102 ! Surface & bottom boundary conditions: no flux 103 DO ji = 2, jpim1 104 zwy(ji, 1 ) = 0.e0 105 zwy(ji,jpk) = 0.e0 106 zww(ji, 1 ) = 0.e0 107 zww(ji,jpk) = 0.e0 108 END DO 109 110 ! zwx and zwz arrays set to before tracer values 111 DO jk = 1, jpk 112 DO ji = 2, jpim1 113 zwx(ji,jk) = tb(ji,jj,jk) 114 zwz(ji,jk) = sb(ji,jj,jk) 115 END DO 116 END DO 117 118 119 ! 2. Time splitting loop 120 ! ---------------------- 121 122 DO jl = 1, n_zdfexp 123 124 ! first vertical derivative 125 IF( lk_zdfddm ) THEN ! double diffusion: avs /= avt 126 DO jk = 2, jpk 127 DO ji = 2, jpim1 128 zave3r = 1.e0 / fse3w(ji,jj,jk) 129 zwy(ji,jk) = avt(ji,jj,jk) * ( zwx(ji,jk-1) - zwx(ji,jk) ) * zave3r 130 zww(ji,jk) = fsavs(ji,jj,jk) * ( zwz(ji,jk-1) - zwz(ji,jk) ) * zave3r 131 END DO 132 END DO 133 ELSE ! default : avs = avt 134 DO jk = 2, jpk 135 DO ji = 2, jpim1 136 zave3r = avt(ji,jj,jk) / fse3w(ji,jj,jk) 137 zwy(ji,jk) = zave3r *(zwx(ji,jk-1) - zwx(ji,jk) ) 138 zww(ji,jk) = zave3r *(zwz(ji,jk-1) - zwz(ji,jk) ) 139 END DO 140 END DO 141 ENDIF 142 143 ! trend estimation at kt+l*2*rdt/n_zdfexp 144 DO jk = 1, jpkm1 145 DO ji = 2, jpim1 146 ze3tr = zlavmr / fse3t(ji,jj,jk) 147 ! 2nd vertical derivative 148 zta = ( zwy(ji,jk) - zwy(ji,jk+1) ) * ze3tr 149 zsa = ( zww(ji,jk) - zww(ji,jk+1) ) * ze3tr 150 ! update the tracer trends 151 ta(ji,jj,jk) = ta(ji,jj,jk) + zta 152 sa(ji,jj,jk) = sa(ji,jj,jk) + zsa 153 ! update tracer fields at kt+l*2*rdt/n_zdfexp 154 zwx(ji,jk) = zwx(ji,jk) + p2dt(jk) * zta * tmask(ji,jj,jk) 155 zwz(ji,jk) = zwz(ji,jk) + p2dt(jk) * zsa * tmask(ji,jj,jk) 99 ! Split-explicit loop (after tracer due to the vertical diffusion alone) 100 ! ------------------- 101 ! 102 DO jl = 1, n_zdfexp 103 ! ! first vertical derivative 104 DO jk = 2, jpk 105 DO jj = 2, jpjm1 106 DO ji = fs_2, fs_jpim1 ! vector opt. 107 zave3r = 1.e0 / fse3w(ji,jj,jk) 108 zwy(ji,jj,jk) = avt(ji,jj,jk) * ( zwx(ji,jj,jk-1) - zwx(ji,jj,jk) ) * zave3r 109 zww(ji,jj,jk) = fsavs(ji,jj,jk) * ( zwz(ji,jj,jk-1) - zwz(ji,jj,jk) ) * zave3r 156 110 END DO 157 111 END DO 158 112 END DO 159 ! ! =============== 160 END DO ! End of slab 161 ! ! =============== 113 ! 114 DO jk = 1, jpkm1 ! second vertical derivative ==> tracer at kt+l*2*rdt/n_zdfexp 115 DO jj = 2, jpjm1 116 DO ji = fs_2, fs_jpim1 ! vector opt. 117 ze3tr = zlavmr / fse3t(ji,jj,jk) 118 zwx(ji,jj,jk) = zwx(ji,jj,jk) + p2dt(jk) * ( zwy(ji,jj,jk) - zwy(ji,jj,jk+1) ) * ze3tr 119 zwz(ji,jj,jk) = zwz(ji,jj,jk) + p2dt(jk) * ( zww(ji,jj,jk) - zww(ji,jj,jk+1) ) * ze3tr 120 END DO 121 END DO 122 END DO 123 ! 124 END DO 162 125 126 ! After tracer due to all trends 127 ! ------------------------------ 128 IF( lk_vvl ) THEN ! variable level thickness : leap-frog on tracer*e3t 129 DO jk = 1, jpkm1 130 DO jj = 2, jpjm1 131 DO ji = fs_2, fs_jpim1 ! vector opt. 132 ze3tb = fsve3t(ji,jj,jk) * ( 1 + sshb(ji,jj) * mut(ji,jj,jk) ) ! before e3t 133 zta = zwx(ji,jj,jk) - tb(ji,jj,jk) + p2dt(jk) * ta(ji,jj,jk) ! total trends * 2*rdt 134 zsa = zwz(ji,jj,jk) - sb(ji,jj,jk) + p2dt(jk) * sa(ji,jj,jk) 135 zcoef = 1.e0 / fse3t(ji,jj,jk) * tmask(ji,jj,jk) 136 ta(ji,jj,jk) = ( ze3tb * tb(ji,jj,jk) + fse3t(ji,jj,jk) * zta ) * zcoef 137 sa(ji,jj,jk) = ( ze3tb * sb(ji,jj,jk) + fse3t(ji,jj,jk) * zsa ) * zcoef 138 END DO 139 END DO 140 END DO 141 ELSE ! fixed level thickness : leap-frog on tracers 142 DO jk = 1, jpkm1 143 DO jj = 2, jpjm1 144 DO ji = fs_2, fs_jpim1 ! vector opt. 145 ta(ji,jj,jk) = ( zwx(ji,jj,jk) + p2dt(jk) * ta(ji,jj,jk) ) *tmask(ji,jj,jk) 146 sa(ji,jj,jk) = ( zwz(ji,jj,jk) + p2dt(jk) * sa(ji,jj,jk) ) *tmask(ji,jj,jk) 147 END DO 148 END DO 149 END DO 150 ENDIF 151 ! 163 152 END SUBROUTINE tra_zdf_exp 164 153
Note: See TracChangeset
for help on using the changeset viewer.