New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
traldf_bilap.F90 in branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/TRA – NEMO

source: branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/TRA/traldf_bilap.F90 @ 2082

Last change on this file since 2082 was 2082, checked in by cetlod, 14 years ago

Improve the merge of TRA-TRC, see ticket #717
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 9.1 KB
Line 
1MODULE traldf_bilap
2   !!==============================================================================
3   !!                   ***  MODULE  traldf_bilap  ***
4   !! Ocean  tracers:  horizontal component of the lateral tracer mixing trend
5   !!==============================================================================
6   !! History :     !  91-11  (G. Madec)  Original code
7   !!               !  93-03  (M. Guyon)  symetrical conditions
8   !!               !  95-11  (G. Madec)  suppress volumetric scale factors
9   !!               !  96-01  (G. Madec)  statement function for e3
10   !!               !  96-01  (M. Imbard)  mpp exchange
11   !!               !  97-07  (G. Madec)  optimization, and ahtt
12   !!          8.5  !  02-08  (G. Madec)  F90: Free form and module
13   !!          9.0  !  04-08  (C. Talandier) New trends organization
14   !!               !  05-11  (G. Madec)  zps or sco as default option
15   !!          3.3  !  10-05  (C. Ethe, G. Madec)  merge TRC-TRA
16   !!==============================================================================
17
18   !!----------------------------------------------------------------------
19   !!   tra_ldf_bilap : update the tracer trend with the horizontal diffusion
20   !!                   using a iso-level biharmonic operator
21   !!----------------------------------------------------------------------
22   !! * Modules used
23   USE oce             ! ocean dynamics and active tracers
24   USE dom_oce         ! ocean space and time domain
25   USE ldftra_oce      ! ocean tracer   lateral physics
26   USE in_out_manager  ! I/O manager
27   USE ldfslp          ! iso-neutral slopes
28   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
29   USE diaptr          ! poleward transport diagnostics
30   USE trc_oce         ! share passive tracers/Ocean variables
31
32   IMPLICIT NONE
33   PRIVATE
34
35   !! * Routine accessibility
36   PUBLIC tra_ldf_bilap   ! routine called by step.F90
37
38   !! * Substitutions
39#  include "domzgr_substitute.h90"
40#  include "ldftra_substitute.h90"
41#  include "ldfeiv_substitute.h90"
42#  include "vectopt_loop_substitute.h90"
43   !!----------------------------------------------------------------------
44   !! NEMO/OPA 3.3 , LOCEAN-IPSL (2010)
45   !! $Id$
46   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
47   !!----------------------------------------------------------------------
48
49CONTAINS
50 
51   SUBROUTINE tra_ldf_bilap( kt, cdtype, pgu, pgv,  &
52      &                                  ptb, pta, kjpt ) 
53      !!----------------------------------------------------------------------
54      !!                  ***  ROUTINE tra_ldf_bilap  ***
55      !!
56      !! ** Purpose :   Compute the before horizontal tracer diffusive
57      !!      trend and add it to the general trend of tracer equation.
58      !!
59      !! ** Method  :   4th order diffusive operator along model level surfaces
60      !!      evaluated using before fields (forward time scheme). The hor.
61      !!      diffusive trends  is given by:
62      !!      Laplacian of tb:
63      !!         zlt   = 1/(e1t*e2t*e3t) {  di-1[ e2u*e3u/e1u di(tb) ]
64      !!                                  + dj-1[ e1v*e3v/e2v dj(tb) ]  }
65      !!      Multiply by the eddy diffusivity coef. and insure lateral bc:
66      !!        zlt   = ahtt * zlt
67      !!        call to lbc_lnk
68      !!      Bilaplacian (laplacian of zlt):
69      !!         difft = 1/(e1t*e2t*e3t) {  di-1[ e2u*e3u/e1u di(zlt) ]
70      !!                                  + dj-1[ e1v*e3v/e2v dj(zlt) ]  }
71      !!
72      !!      Add this trend to the general trend
73      !!         (pta) = (pta) + ( difft )
74      !!
75      !! ** Action : - Update pta arrays with the before iso-level
76      !!               biharmonic mixing trend.
77      !!----------------------------------------------------------------------
78      !!
79      USE oce         , ztu => ua   ! use ua as workspace
80      USE oce         , ztv => va   ! use va as workspace
81      !!
82      INTEGER         , INTENT(in   )                                ::   kt             ! ocean time-step index
83      CHARACTER(len=3), INTENT(in   )                                ::   cdtype         ! =TRA or TRC (tracer indicator)
84      INTEGER         , INTENT(in   )                                ::   kjpt            ! number of tracers
85      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,kjpt  )     ::   pgu, pgv     ! tracer gradient at pstep levels
86      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk,kjpt)   ::   ptb          ! before and now tracer fields
87      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)   ::   pta          ! tracer trend
88      !!
89      INTEGER  ::  ji, jj, jk, jn         ! dummy loop indices
90      INTEGER  ::  iku, ikv               ! temporary integers
91      REAL(wp) ::  zbtr, ztra             ! temporary scalars
92      REAL(wp), DIMENSION(jpi,jpj) ::   & 
93         zeeu, zeev, zlt               ! 2D workspace
94      !!----------------------------------------------------------------------
95
96      IF( ( cdtype == 'TRA' .AND. kt == nit000 ) .OR. ( cdtype == 'TRC' .AND. kt == nittrc000 ) )  THEN
97         IF(lwp) WRITE(numout,*)
98         IF(lwp) WRITE(numout,*) 'tra_ldf_bilap : iso-level biharmonic operator on ', cdtype
99         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~'
100      ENDIF
101      !                                                          ! ===========
102      DO jn = 1, kjpt                                            ! tracer loop
103         !                                                       ! ===========
104         !                                               
105         DO jk = 1, jpkm1                                 
106            !                                             
107           
108            ! 0. Initialization of metric arrays (for z- or s-coordinates)
109            ! ----------------------------------
110            DO jj = 1, jpjm1
111               DO ji = 1, fs_jpim1   ! vector opt.
112                  zeeu(ji,jj) = e2u(ji,jj) * fse3u(ji,jj,jk) / e1u(ji,jj) * umask(ji,jj,jk)
113                  zeev(ji,jj) = e1v(ji,jj) * fse3v(ji,jj,jk) / e2v(ji,jj) * vmask(ji,jj,jk)
114               END DO
115            END DO
116
117
118            ! 1. Laplacian
119            ! ------------
120           
121            ! First derivative (gradient)
122            DO jj = 1, jpjm1
123               DO ji = 1, fs_jpim1   ! vector opt.
124                  ztu(ji,jj,jk) = zeeu(ji,jj) * ( ptb(ji+1,jj  ,jk,jn) - ptb(ji,jj,jk,jn) )
125                  ztv(ji,jj,jk) = zeev(ji,jj) * ( ptb(ji  ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) )
126               END DO
127            END DO
128            IF( ln_zps ) THEN      ! set gradient at partial step level
129               DO jj = 1, jpjm1
130                  DO ji = 1, jpim1
131                     ! last level
132                     iku = MIN ( mbathy(ji,jj), mbathy(ji+1,jj  ) ) - 1
133                     ikv = MIN ( mbathy(ji,jj), mbathy(ji  ,jj+1) ) - 1
134                     IF( iku == jk )  ztu(ji,jj,jk) = zeeu(ji,jj) * pgu(ji,jj,jn)
135                     IF( ikv == jk )  ztv(ji,jj,jk) = zeev(ji,jj) * pgv(ji,jj,jn)
136                  END DO
137               END DO
138            ENDIF
139
140            ! Second derivative (divergence) multiply by the eddy diffusivity coefficient
141            DO jj = 2, jpjm1
142               DO ji = fs_2, fs_jpim1   ! vector opt.
143                  zbtr = 1.0 / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
144                  zlt(ji,jj) = fsahtt(ji,jj,jk) &
145                     &                * zbtr * (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  )
146               END DO
147            END DO
148
149            ! Lateral boundary conditions on the laplacian (zlt)   (unchanged sgn)
150            CALL lbc_lnk( zlt, 'T', 1. )
151
152            ! 2. Bilaplacian
153            ! --------------
154           
155            ! third derivative (gradient)
156            DO jj = 1, jpjm1
157               DO ji = 1, fs_jpim1   ! vector opt.
158                  ztu(ji,jj,jk) = zeeu(ji,jj) * ( zlt(ji+1,jj  ) - zlt(ji,jj) )
159                  ztv(ji,jj,jk) = zeev(ji,jj) * ( zlt(ji  ,jj+1) - zlt(ji,jj) )
160               END DO
161            END DO
162
163            ! fourth derivative (divergence) and add to the general tracer trend
164            DO jj = 2, jpjm1
165               DO ji = fs_2, fs_jpim1   ! vector opt.
166                  ! horizontal diffusive trends
167                  zbtr = 1.0 / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
168                  ztra = zbtr * (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  )
169                  ! add it to the general tracer trends
170                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra
171               END DO
172            END DO
173            !                                             ! ===============
174         END DO                                           ! Horizontal slab
175         !                                                ! ===============
176         ! "zonal" mean lateral diffusive heat and salt transport
177         IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 
178           IF( jn == jp_tem )  pht_ldf(:) = ptr_vj( ztv(:,:,:) )
179           IF( jn == jp_sal )  pst_ldf(:) = ptr_vj( ztv(:,:,:) )
180         ENDIF
181         !
182      END DO
183
184   END SUBROUTINE tra_ldf_bilap
185
186   !!==============================================================================
187END MODULE traldf_bilap
Note: See TracBrowser for help on using the repository browser.