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.
ldfdyn.F90 in trunk/NEMO/OPA_SRC/LDF – NEMO

source: trunk/NEMO/OPA_SRC/LDF/ldfdyn.F90 @ 1152

Last change on this file since 1152 was 1152, checked in by rblod, 16 years ago

Convert cvs header to svn Id, step II
  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Id
File size: 12.8 KB
Line 
1MODULE ldfdyn
2   !!======================================================================
3   !!                       ***  MODULE  ldfdyn  ***
4   !! Ocean physics:  lateral viscosity coefficient
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !!   ldf_dyn_init : initialization, namelist read, and parameters control
9   !!   ldf_dyn_c3d   : 3D eddy viscosity coefficient initialization
10   !!   ldf_dyn_c2d   : 2D eddy viscosity coefficient initialization
11   !!   ldf_dyn_c1d   : 1D eddy viscosity coefficient initialization
12   !!----------------------------------------------------------------------
13   !! * Modules used
14   USE oce             ! ocean dynamics and tracers   
15   USE dom_oce         ! ocean space and time domain
16   USE ldfdyn_oce      ! ocean dynamics lateral physics
17   USE phycst          ! physical constants
18   USE ldfslp          ! ???
19   USE in_out_manager  ! I/O manager
20   USE lib_mpp         ! distribued memory computing library
21   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
22
23   IMPLICIT NONE
24   PRIVATE
25
26   !! *  Routine accessibility
27   PUBLIC ldf_dyn_init   ! called by opa.F90
28
29  INTERFACE ldf_zpf
30     MODULE PROCEDURE ldf_zpf_1d, ldf_zpf_1d_3d, ldf_zpf_3d
31  END INTERFACE
32
33   !! * Substitutions
34#  include "domzgr_substitute.h90"
35   !!----------------------------------------------------------------------
36   !!   OPA 9.0 , LOCEAN-IPSL (2005)
37   !! $Id$
38   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
39   !!----------------------------------------------------------------------
40
41CONTAINS
42
43   SUBROUTINE ldf_dyn_init
44      !!----------------------------------------------------------------------
45      !!                  ***  ROUTINE ldf_dyn_init  ***
46      !!                   
47      !! ** Purpose :   set the horizontal ocean dynamics physics
48      !!
49      !! ** Method  : 
50      !!      Eddy viscosity coefficients:
51      !!         default option   : constant coef. ahm0 (namelist)
52      !!        'key_dynldf_c1d': depth dependent coef. defined in
53      !!                        in ldf_dyn_c1d routine
54      !!        'key_dynldf_c2d': latitude and longitude dependent coef.
55      !!                        defined in ldf_dyn_c2d routine
56      !!        'key_dynldf_c3d': latitude, longitude, depth dependent coef.
57      !!                        defined in ldf_dyn_c3d routine
58      !!      N.B. User defined include files.  By default, 3d and 2d coef.
59      !!      are set to a constant value given in the namelist and the 1d
60      !!      coefficients are initialized to a hyperbolic tangent vertical
61      !!      profile.
62      !!
63      !! Reference :
64      !!      Madec, G. and M. Imbard, 1996, A global ocean mesh to overcome
65      !!      the North Pole singularity, Climate Dynamics, 12, 381-388.
66      !!
67      !! History :
68      !!        !  07-97  (G. Madec)  from inimix.F split in 2 routines
69      !!        !  08-97  (G. Madec)  multi dimensional coefficients
70      !!   8.5  !  02-09  (G. Madec)  F90: Free form and module
71      !!----------------------------------------------------------------------
72      !! * Modules used
73      USE ioipsl
74
75      !! * Local declarations
76      INTEGER ::   ioptio         ! ???
77      LOGICAL :: ll_print = .FALSE.    ! Logical flag for printing viscosity coef.
78
79       
80      NAMELIST/nam_dynldf/ ln_dynldf_lap  , ln_dynldf_bilap,                &
81         &                 ln_dynldf_level, ln_dynldf_hor, ln_dynldf_iso,   &
82         &                 ahm0, ahmb0
83      !!----------------------------------------------------------------------
84
85
86      ! Define the lateral physics parameters
87      ! ======================================
88   
89      ! Read Namelist nam_dynldf : Lateral physics
90      REWIND( numnam )
91      READ  ( numnam, nam_dynldf )
92
93      ! Parameter print
94      IF(lwp) THEN
95         WRITE(numout,*)
96         WRITE(numout,*) 'ldf_dyn : lateral momentum physics'
97         WRITE(numout,*) '~~~~~~~'
98         WRITE(numout,*) '          Namelist nam_dynldf : set lateral mixing parameters'
99         WRITE(numout,*) '             laplacian operator          ln_dynldf_lap   = ', ln_dynldf_lap
100         WRITE(numout,*) '             bilaplacian operator        ln_dynldf_bilap = ', ln_dynldf_bilap
101         WRITE(numout,*) '             iso-level                   ln_dynldf_level = ', ln_dynldf_level
102         WRITE(numout,*) '             horizontal (geopotential)   ln_dynldf_hor   = ', ln_dynldf_hor
103         WRITE(numout,*) '             iso-neutral                 ln_dynldf_iso   = ', ln_dynldf_iso
104         WRITE(numout,*) '             horizontal eddy viscosity            ahm0   = ', ahm0
105         WRITE(numout,*) '             background viscosity                 ahmb0  = ', ahmb0
106      ENDIF
107
108      ! ... check of lateral diffusive operator on tracers
109      !           ==> will be done in trazdf module
110
111      ! ... Space variation of eddy coefficients
112      ioptio = 0
113#if defined key_dynldf_c3d
114      IF(lwp) WRITE(numout,*) '          momentum mixing coef. = F( latitude, longitude, depth)'
115      ioptio = ioptio+1
116#endif
117#if defined key_dynldf_c2d
118      IF(lwp) WRITE(numout,*) '          momentum mixing coef. = F( latitude, longitude)'
119      ioptio = ioptio+1
120#endif
121#if defined key_dynldf_c1d
122      IF(lwp) WRITE(numout,*) '          momentum mixing coef. = F( depth )'
123      ioptio = ioptio+1
124      IF( ln_sco ) CALL ctl_stop( '          key_dynldf_c1d cannot be used in s-coordinate (ln_sco)' )
125#endif
126      IF( ioptio == 0 ) THEN
127          IF(lwp) WRITE(numout,*) '          momentum mixing coef. = constant  (default option)'
128        ELSEIF( ioptio > 1 ) THEN
129           CALL ctl_stop( '          use only one of the following keys:',   &
130                &         ' key_dynldf_c3d, key_dynldf_c2d, key_dynldf_c1d' )
131      ENDIF
132
133
134      IF( ln_dynldf_bilap ) THEN
135         IF(lwp) WRITE(numout,*) '          biharmonic momentum diffusion'
136         IF( ahm0 > 0 .AND. .NOT. lk_esopa )   &
137              &   CALL ctl_stop( 'The horizontal viscosity coef. ahm0 must be negative' )
138      ELSE
139         IF(lwp) WRITE(numout,*) '          harmonic momentum diff. (default)'
140         IF( ahm0 < 0 .AND. .NOT. lk_esopa ) &
141              &   CALL ctl_stop( '          The horizontal viscosity coef. ahm0 must be positive' )
142      ENDIF
143
144
145      ! Lateral eddy viscosity
146      ! ======================
147
148#if defined key_dynldf_c3d
149      CALL ldf_dyn_c3d( ll_print )   ! ahm = 3D coef. = F( longitude, latitude, depth )
150#elif defined key_dynldf_c2d
151      CALL ldf_dyn_c2d( ll_print )   ! ahm = 1D coef. = F( longitude, latitude )
152#elif defined key_dynldf_c1d
153      CALL ldf_dyn_c1d( ll_print )   ! ahm = 1D coef. = F( depth )
154#else
155                                     ! Constant coefficients
156      IF(lwp) WRITE(numout,*)
157      IF(lwp) WRITE(numout,*) 'inildf: constant eddy viscosity coef. '
158      IF(lwp) WRITE(numout,*) '~~~~~~'
159      IF(lwp) WRITE(numout,*) '        ahm1 = ahm2 = ahm0 =  ',ahm0
160#endif
161
162   END SUBROUTINE ldf_dyn_init
163
164#if defined key_dynldf_c3d
165#   include "ldfdyn_c3d.h90"
166#elif defined key_dynldf_c2d
167#   include "ldfdyn_c2d.h90"
168#elif defined key_dynldf_c1d
169#   include "ldfdyn_c1d.h90"
170#endif
171
172
173   SUBROUTINE ldf_zpf_1d( ld_print, pdam, pwam, pbot, pdep, pah )
174      !!----------------------------------------------------------------------
175      !!                  ***  ROUTINE ldf_zpf  ***
176      !!                   
177      !! ** Purpose :   vertical adimensional profile for eddy coefficient
178      !!
179      !! ** Method  :   1D eddy viscosity coefficients ( depth )
180      !!
181      !!----------------------------------------------------------------------
182      !! * Arguments
183      LOGICAL , INTENT (in   ) :: ld_print   ! If true, output arrays on numout
184      REAL(wp), INTENT (in   ) ::   &
185          pdam,     &  ! depth of the inflection point
186          pwam,     &  ! width of inflection
187          pbot         ! battom value (0<pbot<= 1)
188      REAL(wp), INTENT (in   ), DIMENSION(jpk) ::   &
189          pdep         ! depth of the gridpoint (T, U, V, F)
190      REAL(wp), INTENT (inout), DIMENSION(jpk) ::   &
191          pah          ! adimensional vertical profile
192
193      !! * Local variables
194      INTEGER  ::   jk           ! dummy loop indices
195      REAL(wp) ::   zm00, zm01, zmhb, zmhs       ! temporary scalars
196      !!----------------------------------------------------------------------
197
198      zm00 = TANH( ( pdam - gdept_0(1    ) ) / pwam )
199      zm01 = TANH( ( pdam - gdept_0(jpkm1) ) / pwam )
200      zmhs = zm00 / zm01
201      zmhb = ( 1.e0 - pbot ) / ( 1.e0 - zmhs ) / zm01
202
203      DO jk = 1, jpk
204         pah(jk) = 1.e0 + zmhb * ( zm00 - TANH( ( pdam - pdep(jk) ) / pwam )  )
205      END DO
206
207      ! Control print
208      IF(lwp .AND. ld_print ) THEN
209         WRITE(numout,*)
210         WRITE(numout,*) '         ahm profile : '
211         WRITE(numout,*)
212         WRITE(numout,'("  jk      ahm       ","  depth t-level " )')
213         DO jk = 1, jpk
214            WRITE(numout,'(i6,2f12.4,3x,2f12.4)') jk, pah(jk), pdep(jk)
215         END DO
216      ENDIF
217
218   END SUBROUTINE ldf_zpf_1d
219
220
221   SUBROUTINE ldf_zpf_1d_3d( ld_print, pdam, pwam, pbot, pdep, pah )
222      !!----------------------------------------------------------------------
223      !!                  ***  ROUTINE ldf_zpf  ***
224      !!
225      !! ** Purpose :   vertical adimensional profile for eddy coefficient
226      !!
227      !! ** Method  :   1D eddy viscosity coefficients ( depth )
228      !!
229      !!----------------------------------------------------------------------
230      !! * Arguments
231      LOGICAL , INTENT (in   ) :: ld_print   ! If true, output arrays on numout
232      REAL(wp), INTENT (in   ) ::   &
233          pdam,     &  ! depth of the inflection point
234          pwam,     &  ! width of inflection
235          pbot         ! battom value (0<pbot<= 1)
236      REAL(wp), INTENT (in   ), DIMENSION(jpk) ::   &
237          pdep         ! depth of the gridpoint (T, U, V, F)
238      REAL(wp), INTENT (inout), DIMENSION(jpi,jpj,jpk) ::   &
239          pah          ! adimensional vertical profile
240
241      !! * Local variables
242      INTEGER  ::   jk           ! dummy loop indices
243      REAL(wp) ::   zm00, zm01, zmhb, zmhs, zcf  ! temporary scalars
244      !!----------------------------------------------------------------------
245
246      zm00 = TANH( ( pdam - gdept_0(1    ) ) / pwam )
247      zm01 = TANH( ( pdam - gdept_0(jpkm1) ) / pwam )
248      zmhs = zm00 / zm01
249      zmhb = ( 1.e0 - pbot ) / ( 1.e0 - zmhs ) / zm01
250
251      DO jk = 1, jpk
252         zcf = 1.e0 + zmhb * ( zm00 - TANH( ( pdam - pdep(jk) ) / pwam )  )
253         pah(:,:,jk) = zcf
254      END DO
255
256      ! Control print
257      IF(lwp .AND. ld_print ) THEN
258         WRITE(numout,*)
259         WRITE(numout,*) '         ahm profile : '
260         WRITE(numout,*)
261         WRITE(numout,'("  jk      ahm       ","  depth t-level " )')
262         DO jk = 1, jpk
263            WRITE(numout,'(i6,2f12.4,3x,2f12.4)') jk, pah(1,1,jk), pdep(jk)
264         END DO
265      ENDIF
266
267   END SUBROUTINE ldf_zpf_1d_3d
268
269
270   SUBROUTINE ldf_zpf_3d( ld_print, pdam, pwam, pbot, pdep, pah )
271      !!----------------------------------------------------------------------
272      !!                  ***  ROUTINE ldf_zpf  ***
273      !!
274      !! ** Purpose :   vertical adimensional profile for eddy coefficient
275      !!
276      !! ** Method  :   3D for partial step or s-coordinate
277      !!
278      !!----------------------------------------------------------------------
279      !! * Arguments
280      LOGICAL , INTENT (in   ) :: ld_print   ! If true, output arrays on numout
281      REAL(wp), INTENT (in   ) ::   &
282          pdam,     &  ! depth of the inflection point
283          pwam,     &  ! width of inflection
284          pbot         ! reduction factor (surface value / bottom value)
285      REAL(wp), INTENT (in   ), DIMENSION(jpi,jpj,jpk) ::   &
286          pdep         ! dep of the gridpoint (T, U, V, F)
287      REAL(wp), INTENT (inout), DIMENSION(jpi,jpj,jpk) ::   &
288          pah          ! adimensional vertical profile
289
290      !! * Local variables
291      INTEGER  ::   jk           ! dummy loop indices
292      REAL(wp) ::   zm00, zm01, zmhb, zmhs       ! temporary scalars
293      !!----------------------------------------------------------------------
294
295      zm00 = TANH( ( pdam - gdept_0(1    ) ) / pwam )   
296      zm01 = TANH( ( pdam - gdept_0(jpkm1) ) / pwam )
297      zmhs = zm00 / zm01
298      zmhb = ( 1.e0 - pbot ) / ( 1.e0 - zmhs ) / zm01
299
300      DO jk = 1, jpk
301         pah(:,:,jk) = 1.e0 + zmhb * ( zm00 - TANH( ( pdam - pdep(:,:,jk) ) / pwam )  )
302      END DO
303
304      ! Control print
305      IF(lwp .AND. ld_print ) THEN
306         WRITE(numout,*)
307         WRITE(numout,*) '         ahm profile : '
308         WRITE(numout,*)
309         WRITE(numout,'("  jk      ahm       ","  depth t-level " )')
310         DO jk = 1, jpk
311            WRITE(numout,'(i6,2f12.4,3x,2f12.4)') jk, pah(1,1,jk), pdep(1,1,jk)
312         END DO
313      ENDIF
314
315   END SUBROUTINE ldf_zpf_3d
316   !!======================================================================
317END MODULE ldfdyn
Note: See TracBrowser for help on using the repository browser.