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.
dynldf.F90 in branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/DYN – NEMO

source: branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf.F90 @ 11394

Last change on this file since 11394 was 11394, checked in by mattmartin, 5 years ago

First implementation of STOPACK in the GO6 package branch.
File size: 14.2 KB
Line 
1MODULE dynldf
2   !!======================================================================
3   !!                       ***  MODULE  dynldf  ***
4   !! Ocean physics:  lateral diffusivity trends
5   !!=====================================================================
6   !! History :  9.0  !  05-11  (G. Madec)  Original code (new step architecture)
7   !!----------------------------------------------------------------------
8
9   !!----------------------------------------------------------------------
10   !!   dyn_ldf      : update the dynamics trend with the lateral diffusion
11   !!   dyn_ldf_init : initialization, namelist read, and parameters control
12   !!----------------------------------------------------------------------
13   USE oce            ! ocean dynamics and tracers
14   USE dom_oce        ! ocean space and time domain
15   USE phycst         ! physical constants
16   USE ldfdyn_oce     ! ocean dynamics lateral physics
17   USE ldftra_oce     ! ocean tracers  lateral physics
18   USE ldfslp         ! lateral mixing: slopes of mixing orientation
19   USE dynldf_bilapg  ! lateral mixing            (dyn_ldf_bilapg routine)
20   USE dynldf_bilap   ! lateral mixing            (dyn_ldf_bilap  routine)
21   USE dynldf_iso     ! lateral mixing            (dyn_ldf_iso    routine)
22   USE dynldf_lap     ! lateral mixing            (dyn_ldf_lap    routine)
23   USE trd_oce        ! trends: ocean variables
24   USE trddyn         ! trend manager: dynamics   (trd_dyn        routine)
25   !
26   USE prtctl         ! Print control
27   USE in_out_manager ! I/O manager
28   USE lib_mpp        ! distribued memory computing library
29   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
30   USE wrk_nemo        ! Memory Allocation
31   USE timing          ! Timing
32   USE stopack
33
34   IMPLICIT NONE
35   PRIVATE
36
37   PUBLIC   dyn_ldf       ! called by step module
38   PUBLIC   dyn_ldf_init  ! called by opa  module
39
40   INTEGER ::   nldf = -2   ! type of lateral diffusion used defined from ln_dynldf_... namlist logicals)
41
42#if defined key_dynldf_c3d
43   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ahm10,ahm20,ahm30,ahm40
44#endif
45#if defined key_dynldf_c2d
46   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:  ) :: ahm10,ahm20,ahm30,ahm40
47#endif
48
49   !! * Substitutions
50#  include "domzgr_substitute.h90"
51#  include "vectopt_loop_substitute.h90"
52   !!----------------------------------------------------------------------
53   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
54   !! $Id$
55   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
56   !!----------------------------------------------------------------------
57CONTAINS
58
59   SUBROUTINE dyn_ldf( kt )
60      !!----------------------------------------------------------------------
61      !!                  ***  ROUTINE dyn_ldf  ***
62      !!
63      !! ** Purpose :   compute the lateral ocean dynamics physics.
64      !!----------------------------------------------------------------------
65      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
66      !
67      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrdu, ztrdv
68      !!----------------------------------------------------------------------
69      !
70      IF( nn_timing == 1 )  CALL timing_start('dyn_ldf')
71      !
72      IF( l_trddyn )   THEN                      ! temporary save of ta and sa trends
73         CALL wrk_alloc( jpi, jpj, jpk, ztrdu, ztrdv )
74         ztrdu(:,:,:) = ua(:,:,:) 
75         ztrdv(:,:,:) = va(:,:,:) 
76      ENDIF
77
78#if defined key_dynldf_c3d
79         IF( kt .eq. nit000 .and. (nn_spp_ahm1+nn_spp_ahm2) .gt. 0 ) THEN
80             ALLOCATE ( ahm10(jpi,jpj,jpk), ahm20(jpi,jpj,jpk) )
81             ALLOCATE ( ahm30(jpi,jpj,jpk), ahm40(jpi,jpj,jpk) )
82             ahm10 = ahm1
83             ahm20 = ahm2
84             ahm30 = ahm3
85             ahm40 = ahm4
86         ENDIF
87#endif
88#if defined key_dynldf_c2d
89         IF( kt .eq. nit000 .and. (nn_spp_ahm1+nn_spp_ahm2) .gt. 0 ) THEN
90             ALLOCATE ( ahm10(jpi,jpj), ahm20(jpi,jpj) )
91             ALLOCATE ( ahm30(jpi,jpj), ahm40(jpi,jpj) )
92             ahm10 = ahm1
93             ahm20 = ahm2
94             ahm30 = ahm3
95             ahm40 = ahm4
96         ENDIF
97#endif
98
99#if defined key_traldf_c3d || defined key_traldf_c2d
100         IF( nn_spp_ahm1 .GT. 0) THEN
101           IF( ln_dynldf_lap ) THEN
102             ahm1 = ahm10
103             CALL spp_ahm(kt,ahm1,nn_spp_ahm1,rn_ahm1_sd,jk_spp_ahm1)
104           ENDIF
105           IF( ln_dynldf_bilap ) THEN
106             ahm3 = ahm30
107             CALL spp_ahm(kt,ahm3,nn_spp_ahm1,rn_ahm1_sd,jk_spp_ahm3)
108           ENDIF
109         ENDIF
110         IF( nn_spp_ahm2 .GT. 0) THEN
111           IF( ln_dynldf_lap ) THEN
112             ahm2 = ahm20
113             CALL spp_ahm(kt,ahm2,nn_spp_ahm2,rn_ahm2_sd,jk_spp_ahm2)
114           ENDIF
115           IF( ln_dynldf_bilap ) THEN
116             ahm4 = ahm40
117             CALL spp_ahm(kt,ahm4,nn_spp_ahm2,rn_ahm2_sd,jk_spp_ahm4)
118           ENDIF
119         ENDIF
120#endif
121
122      SELECT CASE ( nldf )                       ! compute lateral mixing trend and add it to the general trend
123      !
124      CASE ( 0 )    ;   CALL dyn_ldf_lap    ( kt )      ! iso-level laplacian
125      CASE ( 1 )    ;   CALL dyn_ldf_iso    ( kt )      ! rotated laplacian (except dk[ dk[.] ] part)
126      CASE ( 2 )    ;   CALL dyn_ldf_bilap  ( kt )      ! iso-level bilaplacian
127      CASE ( 3 )    ;   CALL dyn_ldf_bilapg ( kt )      ! s-coord. horizontal bilaplacian
128      CASE ( 4 )                                        ! iso-level laplacian + bilaplacian
129         CALL dyn_ldf_lap    ( kt )
130         CALL dyn_ldf_bilap  ( kt )
131      CASE ( 5 )                                        ! rotated laplacian + bilaplacian (s-coord)
132         CALL dyn_ldf_iso    ( kt )
133         CALL dyn_ldf_bilapg ( kt )
134      !
135      CASE ( -1 )                                       ! esopa: test all possibility with control print
136                        CALL dyn_ldf_lap    ( kt )
137                        CALL prt_ctl( tab3d_1=ua, clinfo1=' ldf0 - Ua: ', mask1=umask,   &
138            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
139                        CALL dyn_ldf_iso    ( kt )
140                        CALL prt_ctl( tab3d_1=ua, clinfo1=' ldf1 - Ua: ', mask1=umask,   &
141            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
142                        CALL dyn_ldf_bilap  ( kt )
143                        CALL prt_ctl( tab3d_1=ua, clinfo1=' ldf2 - Ua: ', mask1=umask,   &
144            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
145                        CALL dyn_ldf_bilapg ( kt )
146                        CALL prt_ctl( tab3d_1=ua, clinfo1=' ldf3 - Ua: ', mask1=umask,   &
147            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
148      !
149      CASE ( -2 )                                       ! neither laplacian nor bilaplacian schemes used
150         IF( kt == nit000 ) THEN
151            IF(lwp) WRITE(numout,*)
152            IF(lwp) WRITE(numout,*) 'dyn_ldf : no lateral diffusion on momentum setup'
153            IF(lwp) WRITE(numout,*) '~~~~~~~ '
154         ENDIF
155      END SELECT
156
157      IF( l_trddyn ) THEN                        ! save the horizontal diffusive trends for further diagnostics
158         ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
159         ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
160         CALL trd_dyn( ztrdu, ztrdv, jpdyn_ldf, kt )
161         CALL wrk_dealloc( jpi, jpj, jpk, ztrdu, ztrdv )
162      ENDIF
163      !                                          ! print sum trends (used for debugging)
164      IF(ln_ctl)   CALL prt_ctl( tab3d_1=ua, clinfo1=' ldf  - Ua: ', mask1=umask,   &
165         &                       tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
166      !
167      IF( nn_timing == 1 )  CALL timing_stop('dyn_ldf')
168      !
169   END SUBROUTINE dyn_ldf
170
171
172   SUBROUTINE dyn_ldf_init
173      !!----------------------------------------------------------------------
174      !!                  ***  ROUTINE dyn_ldf_init  ***
175      !!
176      !! ** Purpose :   initializations of the horizontal ocean dynamics physics
177      !!----------------------------------------------------------------------
178      INTEGER ::   ioptio, ierr         ! temporary integers
179      !!----------------------------------------------------------------------
180   
181      !                                   ! Namelist nam_dynldf: already read in ldfdyn module
182
183      IF(lwp) THEN                        ! Namelist print
184         WRITE(numout,*)
185         WRITE(numout,*) 'dyn_ldf_init : Choice of the lateral diffusive operator on dynamics'
186         WRITE(numout,*) '~~~~~~~~~~~'
187         WRITE(numout,*) '       Namelist nam_dynldf : set lateral mixing parameters (type, direction, coefficients)'
188         WRITE(numout,*) '          laplacian operator          ln_dynldf_lap   = ', ln_dynldf_lap
189         WRITE(numout,*) '          bilaplacian operator        ln_dynldf_bilap = ', ln_dynldf_bilap
190         WRITE(numout,*) '          iso-level                   ln_dynldf_level = ', ln_dynldf_level
191         WRITE(numout,*) '          horizontal (geopotential)   ln_dynldf_hor   = ', ln_dynldf_hor
192         WRITE(numout,*) '          iso-neutral                 ln_dynldf_iso   = ', ln_dynldf_iso
193      ENDIF
194
195      !                                   ! control the consistency
196      ioptio = 0
197      IF( ln_dynldf_lap   )   ioptio = ioptio + 1
198      IF( ln_dynldf_bilap )   ioptio = ioptio + 1
199      IF( ioptio <  1 ) CALL ctl_warn( '          neither laplacian nor bilaplacian operator set for dynamics' )
200      ioptio = 0
201      IF( ln_dynldf_level )   ioptio = ioptio + 1
202      IF( ln_dynldf_hor   )   ioptio = ioptio + 1
203      IF( ln_dynldf_iso   )   ioptio = ioptio + 1
204      IF( ioptio >  1 ) CALL ctl_stop( '          use only ONE direction (level/hor/iso)' )
205
206      IF( ln_dynldf_iso .AND. ln_traldf_hor ) CALL ctl_stop &
207      &   ( 'Not sensible to use geopotential diffusion for tracers with isoneutral diffusion for dynamics' )
208
209      !                                   ! Set nldf, the type of lateral diffusion, from ln_dynldf_... logicals
210      ierr = 0
211      IF ( ln_dynldf_lap ) THEN      ! laplacian operator
212         IF ( ln_zco ) THEN                ! z-coordinate
213            IF ( ln_dynldf_level )   nldf = 0      ! iso-level  (no rotation)
214            IF ( ln_dynldf_hor   )   nldf = 0      ! horizontal (no rotation)
215            IF ( ln_dynldf_iso   )   nldf = 1      ! isoneutral (   rotation)
216         ENDIF
217         IF ( ln_zps ) THEN             ! z-coordinate
218            IF ( ln_dynldf_level )   ierr = 1      ! iso-level not allowed
219            IF ( ln_dynldf_hor   )   nldf = 0      ! horizontal (no rotation)
220            IF ( ln_dynldf_iso   )   nldf = 1      ! isoneutral (   rotation)
221         ENDIF
222         IF ( ln_sco ) THEN             ! s-coordinate
223            IF ( ln_dynldf_level )   nldf = 0      ! iso-level  (no rotation)
224            IF ( ln_dynldf_hor   )   nldf = 1      ! horizontal (   rotation)
225            IF ( ln_dynldf_iso   )   nldf = 1      ! isoneutral (   rotation)
226         ENDIF
227      ENDIF
228
229      IF( ln_dynldf_bilap ) THEN      ! bilaplacian operator
230         IF ( ln_zco ) THEN                ! z-coordinate
231            IF ( ln_dynldf_level )   nldf = 2      ! iso-level  (no rotation)
232            IF ( ln_dynldf_hor   )   nldf = 2      ! horizontal (no rotation)
233            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
234         ENDIF
235         IF ( ln_zps ) THEN             ! z-coordinate
236            IF ( ln_dynldf_level )   ierr = 1      ! iso-level not allowed
237            IF ( ln_dynldf_hor   )   nldf = 2      ! horizontal (no rotation)
238            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
239         ENDIF
240         IF ( ln_sco ) THEN             ! s-coordinate
241            IF ( ln_dynldf_level )   nldf = 2      ! iso-level  (no rotation)
242            IF ( ln_dynldf_hor   )   nldf = 3      ! horizontal (   rotation)
243            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
244         ENDIF
245      ENDIF
246     
247      IF( ln_dynldf_lap .AND. ln_dynldf_bilap ) THEN  ! mixed laplacian and bilaplacian operators
248         IF ( ln_zco ) THEN                ! z-coordinate
249            IF ( ln_dynldf_level )   nldf = 4      ! iso-level  (no rotation)
250            IF ( ln_dynldf_hor   )   nldf = 4      ! horizontal (no rotation)
251            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
252         ENDIF
253         IF ( ln_zps ) THEN             ! z-coordinate
254            IF ( ln_dynldf_level )   ierr = 1      ! iso-level not allowed
255            IF ( ln_dynldf_hor   )   nldf = 4      ! horizontal (no rotation)
256            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
257         ENDIF
258         IF ( ln_sco ) THEN             ! s-coordinate
259            IF ( ln_dynldf_level )   nldf = 4      ! iso-level  (no rotation)
260            IF ( ln_dynldf_hor   )   nldf = 5      ! horizontal (   rotation)
261            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
262         ENDIF
263      ENDIF
264
265      IF( lk_esopa )                 nldf = -1     ! esopa test
266
267      IF( ierr == 1 )   CALL ctl_stop( 'iso-level in z-coordinate - partial step, not allowed' )
268      IF( ierr == 2 )   CALL ctl_stop( 'isoneutral bilaplacian operator does not exist' )
269      IF( nldf == 1 .OR. nldf == 3 ) THEN      ! rotation
270         IF( .NOT.lk_ldfslp )   CALL ctl_stop( 'the rotation of the diffusive tensor require key_ldfslp' )
271      ENDIF
272
273      IF(lwp) THEN
274         WRITE(numout,*)
275         IF( nldf == -2 )   WRITE(numout,*) '              neither laplacian nor bilaplacian schemes used'
276         IF( nldf == -1 )   WRITE(numout,*) '              ESOPA test All scheme used'
277         IF( nldf ==  0 )   WRITE(numout,*) '              laplacian operator'
278         IF( nldf ==  1 )   WRITE(numout,*) '              rotated laplacian operator'
279         IF( nldf ==  2 )   WRITE(numout,*) '              bilaplacian operator'
280         IF( nldf ==  3 )   WRITE(numout,*) '              rotated bilaplacian'
281         IF( nldf ==  4 )   WRITE(numout,*) '              laplacian and bilaplacian operators'
282         IF( nldf ==  5 )   WRITE(numout,*) '              rotated laplacian and bilaplacian operators'
283      ENDIF
284      !
285   END SUBROUTINE dyn_ldf_init
286
287   !!======================================================================
288END MODULE dynldf
Note: See TracBrowser for help on using the repository browser.