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.
ldfeiv.F90 in branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/LDF – NEMO

source: branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/LDF/ldfeiv.F90 @ 2287

Last change on this file since 2287 was 2287, checked in by smasson, 14 years ago

update licence of all NEMO files...
  • Property svn:keywords set to Id
File size: 9.4 KB
Line 
1MODULE ldfeiv
2   !!======================================================================
3   !!                     ***  MODULE  ldfeiv  ***
4   !! Ocean physics:  variable eddy induced velocity coefficients
5   !!======================================================================
6#if   defined key_traldf_eiv   &&   defined key_traldf_c2d
7   !!----------------------------------------------------------------------
8   !!   'key_traldf_eiv'      and                     eddy induced velocity
9   !!   'key_traldf_c2d'                    2D tracer lateral  mixing coef.
10   !!----------------------------------------------------------------------
11   !!   ldf_eiv      : compute the eddy induced velocity coefficients
12   !!----------------------------------------------------------------------
13   !! * Modules used
14   USE oce             ! ocean dynamics and tracers
15   USE dom_oce         ! ocean space and time domain
16   USE sbc_oce         ! surface boundary condition: ocean
17   USE sbcrnf          ! river runoffs
18   USE ldftra_oce      ! ocean tracer   lateral physics
19   USE phycst          ! physical constants
20   USE ldfslp          ! iso-neutral slopes
21   USE in_out_manager  ! I/O manager
22   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
23   USE prtctl          ! Print control
24   USE iom
25
26   IMPLICIT NONE
27   PRIVATE
28   
29   !! * Routine accessibility
30   PUBLIC ldf_eiv               ! routine called by step.F90
31   !!----------------------------------------------------------------------
32   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
33   !! $Id$
34   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
35   !!----------------------------------------------------------------------
36   !! * Substitutions
37#  include "domzgr_substitute.h90"
38#  include "vectopt_loop_substitute.h90"
39   !!----------------------------------------------------------------------
40
41CONTAINS
42
43   SUBROUTINE ldf_eiv( kt )
44      !!----------------------------------------------------------------------
45      !!                  ***  ROUTINE ldf_eiv  ***
46      !!
47      !! ** Purpose :   Compute the eddy induced velocity coefficient from the
48      !!      growth rate of baroclinic instability.
49      !!
50      !! ** Method :
51      !!
52      !! ** Action : - uslp(),  : i- and j-slopes of neutral surfaces
53      !!             - vslp()      at u- and v-points, resp.
54      !!             - wslpi(),  : i- and j-slopes of neutral surfaces
55      !!             - wslpj()     at w-points.
56      !!
57      !! History :
58      !!   8.1  !  99-03  (G. Madec, A. Jouzeau)  Original code
59      !!   8.5  !  02-06  (G. Madec)  Free form, F90
60      !!----------------------------------------------------------------------
61      !! * Arguments
62      INTEGER, INTENT( in ) ::   kt     ! ocean time-step inedx
63     
64      !! * Local declarations
65      INTEGER ::   ji, jj, jk           ! dummy loop indices
66      REAL(wp) ::   &
67         zfw, ze3w, zn2, zf20,       &  ! temporary scalars
68         zaht, zaht_min
69      REAL(wp), DIMENSION(jpi,jpj) ::   &
70         zn, zah, zhw, zross            ! workspace
71      !!----------------------------------------------------------------------
72     
73      IF( kt == nit000 ) THEN
74         IF(lwp) WRITE(numout,*)
75         IF(lwp) WRITE(numout,*) 'ldf_eiv : eddy induced velocity coefficients'
76         IF(lwp) WRITE(numout,*) '~~~~~~~'
77      ENDIF
78     
79      ! 0. Local initialization
80      ! -----------------------
81      zn   (:,:) = 0.e0
82      zhw  (:,:) = 5.e0
83      zah  (:,:) = 0.e0
84      zross(:,:) = 0.e0
85
86
87      ! 1. Compute lateral diffusive coefficient
88      ! ----------------------------------------
89
90      DO jk = 1, jpk
91#  if defined key_vectopt_loop 
92!CDIR NOVERRCHK
93         DO ji = 1, jpij   ! vector opt.
94            ! Take the max of N^2 and zero then take the vertical sum
95            ! of the square root of the resulting N^2 ( required to compute
96            ! internal Rossby radius Ro = .5 * sum_jpk(N) / f
97            zn2 = MAX( rn2b(ji,1,jk), 0.e0 )
98            zn(ji,1) = zn(ji,1) + SQRT( zn2 ) * fse3w(ji,1,jk)
99            ! Compute elements required for the inverse time scale of baroclinic
100            ! eddies using the isopycnal slopes calculated in ldfslp.F :
101            ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w))
102            ze3w = fse3w(ji,1,jk) * tmask(ji,1,jk)
103               zah(ji,1) = zah(ji,1) + zn2   &
104                              * ( wslpi(ji,1,jk) * wslpi(ji,1,jk)    &
105                                + wslpj(ji,1,jk) * wslpj(ji,1,jk) )   &
106                              * ze3w
107            zhw(ji,1) = zhw(ji,1) + ze3w
108         END DO
109#  else
110         DO jj = 2, jpjm1
111!CDIR NOVERRCHK
112            DO ji = 2, jpim1
113               ! Take the max of N^2 and zero then take the vertical sum
114               ! of the square root of the resulting N^2 ( required to compute
115               ! internal Rossby radius Ro = .5 * sum_jpk(N) / f
116               zn2 = MAX( rn2b(ji,jj,jk), 0.e0 )
117               zn(ji,jj) = zn(ji,jj) + SQRT( zn2 ) * fse3w(ji,jj,jk)
118               ! Compute elements required for the inverse time scale of baroclinic
119               ! eddies using the isopycnal slopes calculated in ldfslp.F :
120               ! T^-1 = sqrt(m_jpk(N^2*(r1^2+r2^2)*e3w))
121               ze3w = fse3w(ji,jj,jk) * tmask(ji,jj,jk)
122               zah(ji,jj) = zah(ji,jj) + zn2   &
123                              * ( wslpi(ji,jj,jk) * wslpi(ji,jj,jk)    &
124                                + wslpj(ji,jj,jk) * wslpj(ji,jj,jk) )  &
125                              * ze3w
126               zhw(ji,jj) = zhw(ji,jj) + ze3w
127            END DO
128         END DO 
129#  endif
130      END DO
131
132      DO jj = 2, jpjm1
133!CDIR NOVERRCHK
134         DO ji = fs_2, fs_jpim1   ! vector opt.
135            zfw = MAX( ABS( 2. * omega * SIN( rad * gphit(ji,jj) ) ) , 1.e-10 )
136            ! Rossby radius at w-point taken < 40km and  > 2km
137            zross(ji,jj) = MAX( MIN( .4 * zn(ji,jj) / zfw, 40.e3 ), 2.e3 )
138            ! Compute aeiw by multiplying Ro^2 and T^-1
139            aeiw(ji,jj) = zross(ji,jj) * zross(ji,jj) * SQRT( zah(ji,jj) / zhw(ji,jj) ) * tmask(ji,jj,1)
140         END DO
141      END DO
142
143      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA R02
144         DO jj = 2, jpjm1
145!CDIR NOVERRCHK
146            DO ji = fs_2, fs_jpim1   ! vector opt.
147               ! Take the minimum between aeiw and aeiv0 for depth levels
148               ! lower than 20 (21 in w- point)
149               IF( mbathy(ji,jj) <= 21. ) aeiw(ji,jj) = MIN( aeiw(ji,jj), 1000. )
150            END DO
151         END DO
152      ENDIF
153
154      ! Decrease the coefficient in the tropics (20N-20S)
155      zf20 = 2. * omega * sin( rad * 20. )
156      DO jj = 2, jpjm1
157         DO ji = fs_2, fs_jpim1   ! vector opt.
158            aeiw(ji,jj) = MIN( 1., ABS( ff(ji,jj) / zf20 ) ) * aeiw(ji,jj)
159         END DO
160      END DO
161
162      ! ORCA R05: Take the minimum between aeiw  and aeiv0
163      IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN
164         DO jj = 2, jpjm1
165            DO ji = fs_2, fs_jpim1   ! vector opt.
166               aeiw(ji,jj) = MIN( aeiw(ji,jj), aeiv0 )
167            END DO
168         END DO
169      ENDIF
170
171      ! lateral boundary condition on aeiw
172      CALL lbc_lnk( aeiw, 'W', 1. )
173
174      ! Average the diffusive coefficient at u- v- points
175      DO jj = 2, jpjm1
176         DO ji = fs_2, fs_jpim1   ! vector opt.
177            aeiu(ji,jj) = .5 * ( aeiw(ji,jj) + aeiw(ji+1,jj  ) )
178            aeiv(ji,jj) = .5 * ( aeiw(ji,jj) + aeiw(ji  ,jj+1) )
179         END DO
180      END DO 
181
182      ! lateral boundary condition on aeiu, aeiv
183      CALL lbc_lnk( aeiu, 'U', 1. )
184      CALL lbc_lnk( aeiv, 'V', 1. )
185
186      IF(ln_ctl)   THEN
187         CALL prt_ctl(tab2d_1=aeiu, clinfo1=' eiv  - u: ', ovlap=1)
188         CALL prt_ctl(tab2d_1=aeiv, clinfo1=' eiv  - v: ', ovlap=1)
189      ENDIF
190
191      ! ORCA R05: add a space variation on aht (=aeiv except at the equator and river mouth)
192      IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN
193         zf20     = 2. * omega * SIN( rad * 20. )
194         zaht_min = 100.                              ! minimum value for aht
195         DO jj = 1, jpj
196            DO ji = 1, jpi
197               zaht      = ( 1. -  MIN( 1., ABS( ff(ji,jj) / zf20 ) ) ) * ( aht0 - zaht_min )  &
198                  &      + aht0 * rnfmsk(ji,jj)                          ! enhanced near river mouths
199               ahtu(ji,jj) = MAX( MAX( zaht_min, aeiu(ji,jj) ) + zaht, aht0 )
200               ahtv(ji,jj) = MAX( MAX( zaht_min, aeiv(ji,jj) ) + zaht, aht0 )
201               ahtw(ji,jj) = MAX( MAX( zaht_min, aeiw(ji,jj) ) + zaht, aht0 )
202            END DO
203         END DO
204         IF(ln_ctl) THEN
205            CALL prt_ctl(tab2d_1=ahtu, clinfo1=' aht  - u: ', ovlap=1)
206            CALL prt_ctl(tab2d_1=ahtv, clinfo1=' aht  - v: ', ovlap=1)
207            CALL prt_ctl(tab2d_1=ahtw, clinfo1=' aht  - w: ', ovlap=1)
208         ENDIF
209      ENDIF
210
211      IF( aeiv0 == 0.e0 ) THEN
212         aeiu(:,:) = 0.e0
213         aeiv(:,:) = 0.e0
214         aeiw(:,:) = 0.e0
215      ENDIF
216
217      CALL iom_put( "aht2d"    , ahtw )   ! lateral eddy diffusivity
218      CALL iom_put( "aht2d_eiv", aeiw )   ! EIV lateral eddy diffusivity
219
220   END SUBROUTINE ldf_eiv
221
222#else
223   !!----------------------------------------------------------------------
224   !!   Default option                                         Dummy module
225   !!----------------------------------------------------------------------
226CONTAINS
227   SUBROUTINE ldf_eiv( kt )       ! Empty routine
228      WRITE(*,*) 'ldf_eiv: You should not have seen this print! error?', kt
229   END SUBROUTINE ldf_eiv
230#endif
231
232   !!======================================================================
233END MODULE ldfeiv
Note: See TracBrowser for help on using the repository browser.