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.
zdfphy.F90 in branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

source: branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfphy.F90 @ 7953

Last change on this file since 7953 was 7953, checked in by gm, 7 years ago

#1880 (HPC-09): add zdfphy (the ZDF manager) + remove all key_...

  • Property svn:keywords set to Id
File size: 14.6 KB
Line 
1MODULE zdfphy
2   !!======================================================================
3   !!                      ***  MODULE  zdfphy  ***
4   !! Ocean physics :   manager of vertical mixing parametrizations
5   !!======================================================================
6   !! History :  4.0  !  2017-04  (G. Madec)  original code
7   !!----------------------------------------------------------------------
8
9   !!----------------------------------------------------------------------
10   !!   zdf_phy_init  : initialization of all vertical physics pakages
11   !!   zdf_phy       : upadate at each time-step the vertical mixing coeff.
12   !!----------------------------------------------------------------------
13   USE par_oce        ! mesh and scale factors
14   USE zdf_oce        ! TKE vertical mixing         
15   USE sbc_oce        ! surface module (only for nn_isf in the option compatibility test)
16   USE zdfbfr         ! bottom friction
17   USE zdftke         ! TKE vertical mixing
18   USE zdfgls         ! GLS vertical mixing
19   USE zdfric         ! Richardson vertical mixing   
20   USE zdfddm         ! double diffusion mixing     
21   USE zdfevd         ! enhanced vertical diffusion
22   USE zdftmx         ! internal tide-induced mixing
23   USE zdfqiao          !Qiao module wave induced mixing   (zdf_qiao routine)
24   USE zdfmxl           ! Mixed-layer depth                (zdf_mxl routine)
25   USE tranpc         ! convection: non penetrative adjustment
26   USE sbcrnf         ! surface boundary condition: runoff variables
27   !
28   USE in_out_manager ! I/O manager
29   USE iom            ! IOM library
30   USE lib_mpp        ! distribued memory computing
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC   zdf_phy_init   ! routine called by nemogcm.F90
36   PUBLIC   zdf_phy        ! routine called by step.F90
37   
38   
39   !!----------------------------------------------------------------------
40   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
41   !! $Id$
42   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46   SUBROUTINE zdf_phy_init
47      !!----------------------------------------------------------------------
48      !!                  ***  ROUTINE zdf_phy_init  ***
49      !!
50      !! ** Purpose :   initializations of the vertical ocean physics
51      !!
52      !! ** Method  :   Read namelist namzdf, control logicals
53      !!----------------------------------------------------------------------
54      INTEGER ::   ioptio, ios   ! local integers
55      !!
56      NAMELIST/namzdf/ ln_zdfcst, ln_zdfric, ln_zdftke, ln_zdfgls,   &     ! type of closure scheme
57         &             ln_zdfevd, nn_evdm, rn_evd ,                  &     ! convection : evd
58         &             ln_zdfnpc, nn_npc , nn_npcp,                  &     ! convection : npc
59         &             ln_zdfddm, rn_avts, rn_hsbfr,                 &     ! double diffusion
60         &             ln_zdftmx,                                    &     ! tidal mixing
61         &             ln_zdfqiao,                                   &     ! surface wave-induced mixing
62         &             ln_zdfexp, nn_zdfexp,                         &     ! time-stepping
63         &             rn_avm0, rn_avt0, nn_avb, nn_havtb                  ! coefficients
64
65
66!!org      NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp,  &
67!!org         &        ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp,       &
68!!org         &        ln_zdfqiao
69      !!----------------------------------------------------------------------
70
71      REWIND( numnam_ref )              ! Namelist namzdf in reference namelist : Vertical mixing parameters
72      READ  ( numnam_ref, namzdf, IOSTAT = ios, ERR = 901)
73901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf in reference namelist', lwp )
74
75      REWIND( numnam_cfg )              ! Namelist namzdf in reference namelist : Vertical mixing parameters
76      READ  ( numnam_cfg, namzdf, IOSTAT = ios, ERR = 902 )
77902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf in configuration namelist', lwp )
78      IF(lwm) WRITE ( numond, namzdf )
79
80      IF(lwp) THEN               !* Parameter print
81         WRITE(numout,*)
82         WRITE(numout,*) 'zdf_phy_init : vertical physics'
83         WRITE(numout,*) '~~~~~~~~'
84         WRITE(numout,*) '   Namelist namzdf : set vertical mixing mixing parameters'
85         WRITE(numout,*) '      vertical closure scheme'
86         WRITE(numout,*) '         constant vertical mixing coefficient    ln_zdfcst = ', ln_zdfcst
87         WRITE(numout,*) '         constant vertical mixing coefficient    ln_zdfric = ', ln_zdfric
88         WRITE(numout,*) '         constant vertical mixing coefficient    ln_zdftke = ', ln_zdftke
89         WRITE(numout,*) '         constant vertical mixing coefficient    ln_zdfgls = ', ln_zdfgls
90         WRITE(numout,*) '      convection: '
91         WRITE(numout,*) '         enhanced vertical diffusion             ln_zdfevd = ', ln_zdfevd
92         WRITE(numout,*) '            applied on momentum (=1/0)             nn_evdm = ', nn_evdm
93         WRITE(numout,*) '            vertical coefficient for evd           rn_evd  = ', rn_evd
94         WRITE(numout,*) '         non-penetrative convection (npc)        ln_zdfnpc = ', ln_zdfnpc
95         WRITE(numout,*) '            npc call  frequency                    nn_npc  = ', nn_npc
96         WRITE(numout,*) '            npc print frequency                    nn_npcp = ', nn_npcp
97         WRITE(numout,*) '      double diffusive mixing                    ln_zdfddm = ', ln_zdfddm
98         WRITE(numout,*) '         maximum avs for dd mixing                 rn_avts = ', rn_avts
99         WRITE(numout,*) '         heat/salt buoyancy flux ratio             rn_hsbfr= ', rn_hsbfr
100         WRITE(numout,*) '      surface wave-induced mixing                ln_zdfqiao= ', ln_zdfqiao                                          ! surface wave induced mixing
101         WRITE(numout,*) '      tidal mixing                               ln_zdftmx = ', ln_zdftmx
102         WRITE(numout,*) '      time splitting / backward scheme           ln_zdfexp = ', ln_zdfexp
103         WRITE(numout,*) '         number of sub-time step (ln_zdfexp=T)   nn_zdfexp = ', nn_zdfexp
104         WRITE(numout,*) '      coefficients : '
105         WRITE(numout,*) '      vertical eddy viscosity             rn_avm0   = ', rn_avm0
106         WRITE(numout,*) '      vertical eddy diffusivity           rn_avt0   = ', rn_avt0
107         WRITE(numout,*) '      constant background or profile      nn_avb    = ', nn_avb
108         WRITE(numout,*) '      horizontal variation for avtb       nn_havtb  = ', nn_havtb
109      ENDIF
110
111!!gm      IF(ln_zdfddm) THEN                    ! double diffusive mixing'
112!         avs(:,:,:) = rn_avt0 * wmask(:,:,:)
113!!gm      ENDIF
114
115      !                          !* Parameter & logical controls
116      !                          !  ----------------------------
117      !
118      !                               ! ... check of vertical mixing scheme on tracers
119      !                                              ==> will be done in trazdf module
120      !
121      !                               ! ... check of mixing coefficient
122      IF(lwp) WRITE(numout,*)
123      IF(lwp) WRITE(numout,*) '   vertical mixing option :'
124      ioptio = 0
125      IF( ln_zdfcst ) THEN
126         IF(lwp) WRITE(numout,*) '      constant eddy diffusion coefficients'
127         ioptio = ioptio+1
128      ENDIF
129      IF( ln_zdfric ) THEN
130         IF(lwp) WRITE(numout,*) '      Richardson dependent eddy coefficients'
131         ioptio = ioptio+1
132      ENDIF
133      IF( ln_zdftke ) THEN
134         IF(lwp) WRITE(numout,*) '      TKE dependent eddy coefficients'
135         ioptio = ioptio+1
136      ENDIF
137      IF( ln_zdfgls ) THEN
138         IF(lwp) WRITE(numout,*) '      GLS dependent eddy coefficients'
139         ioptio = ioptio+1
140      ENDIF
141      IF( ioptio == 0 .OR. ioptio > 1 )   &
142         &   CALL ctl_stop( ' one and only one vertical diffusion option has to be defined ' )
143      IF( ( ln_zdfric .OR. ln_zdfgls ) .AND. ln_isfcav )   &
144         &   CALL ctl_stop( ' only zdfcst and zdftke were tested with ice shelves cavities ' )
145      !
146      !                               ! ... Convection
147      IF(lwp) WRITE(numout,*)
148      IF(lwp) WRITE(numout,*) '   convection :'
149      !
150#if defined key_top
151      IF( ln_zdfnpc )   CALL ctl_stop( ' zdf_phy_init: npc scheme is not working with key_top' )
152#endif
153      !
154      ioptio = 0
155      IF( ln_zdfnpc ) THEN
156         IF(lwp) WRITE(numout,*) '      use non penetrative convective scheme'
157         ioptio = ioptio+1
158      ENDIF
159      IF( ln_zdfevd ) THEN
160         IF(lwp) WRITE(numout,*) '      use enhanced vertical dif. scheme'
161         ioptio = ioptio+1
162      ENDIF
163      IF( ln_zdftke ) THEN
164         IF(lwp) WRITE(numout,*) '      use the 1.5 turbulent closure'
165      ENDIF
166      IF( ln_zdfgls ) THEN
167         IF(lwp) WRITE(numout,*) '      use the GLS closure scheme'
168      ENDIF
169      IF ( ioptio > 1 )   CALL ctl_stop( ' chose between ln_zdfnpc and ln_zdfevd' )
170      IF( ioptio == 0 .AND. .NOT.( ln_zdftke .OR. ln_zdfgls ) )           &
171         CALL ctl_stop( ' except for TKE or GLS physics, a convection scheme is',   &
172         &              ' required: ln_zdfevd or ln_zdfnpc logicals' )
173
174      !                               !* Background eddy viscosity and diffusivity profil
175      IF( nn_avb == 0 ) THEN                ! Define avmb, avtb from namelist parameter
176         avmb(:) = rn_avm0
177         avtb(:) = rn_avt0                     
178      ELSE                                  ! Background profile of avt (fit a theoretical/observational profile (Krauss 1990)
179         avmb(:) = rn_avm0
180         avtb(:) = rn_avt0 + ( 3.e-4_wp - 2._wp * rn_avt0 ) * 1.e-4_wp * gdepw_1d(:)   ! m2/s
181         IF(ln_sco .AND. lwp)   CALL ctl_warn( 'avtb profile not valid in sco' )
182      ENDIF
183      !
184      IF( ln_rstart ) THEN                  !  Read avmb, avtb in restart (if exist)
185         ! if ln_traadv_cen, avmb, avtb have been modified in traadv_cen2 module.
186         ! To ensure the restartability, avmb & avtb are written in the restart
187         ! file in traadv_cen2 end read here.
188         IF( iom_varid( numror, 'avmb', ldstop = .FALSE. ) > 0 ) THEN
189            CALL iom_get( numror, jpdom_unknown, 'avmb', avmb )
190            CALL iom_get( numror, jpdom_unknown, 'avtb', avtb )
191         ENDIF
192      ENDIF
193      !                                     ! 2D shape of the avtb
194      avtb_2d(:,:) = 1.e0                        ! uniform
195      !
196      IF( nn_havtb == 1 ) THEN                   ! decrease avtb in the equatorial band
197           !  -15S -5S : linear decrease from avt0 to avt0/10.
198           !  -5S  +5N : cst value avt0/10.
199           !   5N  15N : linear increase from avt0/10, to avt0
200           WHERE(-15. <= gphit .AND. gphit < -5 )   avtb_2d = (1.  - 0.09 * (gphit + 15.))
201           WHERE( -5. <= gphit .AND. gphit <  5 )   avtb_2d =  0.1
202           WHERE(  5. <= gphit .AND. gphit < 15 )   avtb_2d = (0.1 + 0.09 * (gphit -  5.))
203      ENDIF
204      !
205
206!!gm moved into zdf_phy_init
207!
208                            CALL zdf_bfr_init      ! bottom friction
209
210      ioptio = 0                 !==  type of vertical turbulent closure  ==!    (set nzdfphy)
211      !
212!      IF( ln_zdfcst ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_CST   ;   ENDIF
213!      IF( ln_zdfric ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_RIC   ;   CALL zdf_ric_init   ;   ENDIF
214!      IF( ln_zdftke ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_TKE   ;   CALL zdf_tke_init   ;   ENDIF
215!      IF( ln_zdfgls ) THEN   ;   ioptio = ioptio + 1   ;    nzdf_phy = np_GLS   ;   CALL zdf_gls_init   ;   ENDIF
216
217
218      !
219      IF( ln_zdfric     )   CALL zdf_ric_init      ! Richardson number dependent Kz
220      IF( ln_zdftke     )   CALL zdf_tke_init      ! TKE closure scheme
221      IF( ln_zdfgls     )   CALL zdf_gls_init      ! GLS closure scheme
222      IF( ln_zdftmx     )   CALL zdf_tmx_init      ! tidal vertical mixing
223!!gm
224      !
225   END SUBROUTINE zdf_phy_init
226
227
228   SUBROUTINE zdf_phy( kstp )
229      !!----------------------------------------------------------------------
230      !!                     ***  ROUTINE zdf_phy  ***
231      !!
232      !! ** Purpose :  Update ocean physics at each time-step
233      !!
234      !! ** Method  :
235      !!
236      !! ** Action  :   avm, avt vertical eddy viscosity and diffusivity at w-points
237      !!                nmld ??? mixed layer depth in level and meters   <<<<====verifier !
238      !!                bottom stress.....                               <<<<====verifier !
239      !!----------------------------------------------------------------------
240      INTEGER, INTENT(in) ::   kstp   ! ocean time-step index
241      !
242      INTEGER ::   ji, jj, jk    ! dummy loop indice
243      !!----------------------------------------------------------------------
244      !
245                         CALL zdf_bfr( kstp )         ! bottom friction (if quadratic)
246      !                                               ! Vertical eddy viscosity and diffusivity coefficients
247      IF( ln_zdfric  )   CALL zdf_ric ( kstp )             ! Richardson number dependent Kz
248      IF( ln_zdftke  )   CALL zdf_tke ( kstp )             ! TKE closure scheme for Kz
249      IF( ln_zdfgls  )   CALL zdf_gls ( kstp )             ! GLS closure scheme for Kz
250      IF( ln_zdfqiao )   CALL zdf_qiao( kstp )             ! Qiao vertical mixing
251      !
252      IF( ln_zdfcst  ) THEN                                ! Constant Kz (reset avt, avm[uv] to the background value)
253         avt (:,:,:) = rn_avt0 * wmask (:,:,:)
254         avm (:,:,:) = rn_avm0 * wmask (:,:,:)
255         avmu(:,:,:) = rn_avm0 * wumask(:,:,:)
256         avmv(:,:,:) = rn_avm0 * wvmask(:,:,:)
257      ENDIF
258      !
259      IF( ln_rnf_mouth ) THEN                         ! increase diffusivity at rivers mouths
260         DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2._wp * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO
261      ENDIF
262      !
263      IF( ln_zdfevd  )   CALL zdf_evd( kstp )         ! enhanced vertical eddy diffusivity
264      !
265      IF( ln_zdfddm  ) THEN                           ! double diffusive mixing
266                         CALL zdf_ddm( kstp )
267      ELSE                                            ! avs=avt
268         DO jk = 2, jpkm1   ;   avs(:,:,jk) = avt(:,:,jk)   ;   END DO
269      ENDIF
270      !
271      IF( ln_zdftmx  )   CALL zdf_tmx( kstp )         ! tidal vertical mixing
272
273                         CALL zdf_mxl( kstp )         ! mixed layer depth
274
275                                                      ! write TKE or GLS information in the restart file
276      IF( lrst_oce .AND. ln_zdftke )   CALL tke_rst( kstp, 'WRITE' )
277      IF( lrst_oce .AND. ln_zdfgls )   CALL gls_rst( kstp, 'WRITE' )
278      !
279   END SUBROUTINE zdf_phy
280
281   !!======================================================================
282END MODULE zdfphy
Note: See TracBrowser for help on using the repository browser.