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.
trdmod.F90 in trunk/NEMOGCM/NEMO/OPA_SRC/TRD – NEMO

source: trunk/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod.F90 @ 4624

Last change on this file since 4624 was 4624, checked in by acc, 10 years ago

#1305. Fix slow start-up problems on some systems by introducing and using lwm logical to restrict output of merged namelists to the first (or only) processor. lwm is true only on the first processor regardless of ln_ctl. Small changes to all flavours of nemogcm.F90 are also required to write namctl and namcfg after the call to mynode which now opens output.namelist.dyn and writes nammpp.
  • Property svn:keywords set to Id
File size: 17.2 KB
Line 
1MODULE trdmod
2   !!======================================================================
3   !!                       ***  MODULE  trdmod  ***
4   !! Ocean diagnostics:  ocean tracers and dynamic trends
5   !!=====================================================================
6   !! History :  1.0  !  2004-08  (C. Talandier) Original code
7   !!             -   !  2005-04  (C. Deltel)    Add Asselin trend in the ML budget
8   !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
9   !!----------------------------------------------------------------------
10#if  defined key_trdtra || defined key_trddyn || defined key_trdmld || defined key_trdvor || defined key_esopa
11   !!----------------------------------------------------------------------
12   !!   trd_mod          : Call the trend to be computed
13   !!   trd_mod_init     : Initialization step
14   !!----------------------------------------------------------------------
15   USE oce                     ! ocean dynamics and tracers variables
16   USE dom_oce                 ! ocean space and time domain variables
17   USE zdf_oce                 ! ocean vertical physics variables
18   USE trdmod_oce              ! ocean variables trends
19   USE ldftra_oce              ! ocean active tracers lateral physics
20   USE sbc_oce                 ! surface boundary condition: ocean
21   USE phycst                  ! physical constants
22   USE trdvor                  ! ocean vorticity trends
23   USE trdicp                  ! ocean bassin integral constraints properties
24   USE trdmld                  ! ocean active mixed layer tracers trends
25   USE in_out_manager          ! I/O manager
26   USE lib_mpp                 ! MPP library
27   USE wrk_nemo                ! Memory allocation
28
29
30   IMPLICIT NONE
31   PRIVATE
32
33   REAL(wp) ::   r2dt          ! time-step, = 2 rdttra except at nit000 (=rdttra) if neuler=0
34
35   PUBLIC trd_mod              ! called by all dynXX or traXX modules
36   PUBLIC trd_mod_init         ! called by opa.F90 module
37
38   !! * Substitutions
39#  include "domzgr_substitute.h90"
40#  include "vectopt_loop_substitute.h90"
41   !!----------------------------------------------------------------------
42   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
43   !! $Id$
44   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
45   !!----------------------------------------------------------------------
46
47CONTAINS
48
49   SUBROUTINE trd_mod( ptrdx, ptrdy, ktrd, ctype, kt )
50      !!---------------------------------------------------------------------
51      !!                  ***  ROUTINE trd_mod  ***
52      !!
53      !! ** Purpose : Dispatch all trends computation, e.g. vorticity, mld or
54      !!              integral constraints
55      !!----------------------------------------------------------------------
56      !
57      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdx   ! Temperature or U trend
58      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdy   ! Salinity    or V trend
59      CHARACTER(len=3)          , INTENT(in   ) ::   ctype   ! momentum or tracers trends type 'DYN'/'TRA'
60      INTEGER                   , INTENT(in   ) ::   kt      ! time step
61      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
62      !!
63      INTEGER ::   ji, jj   ! dummy loop indices
64      REAL(wp), POINTER, DIMENSION(:,:)  :: ztswu, ztswv, ztbfu, ztbfv, z2dx, z2dy 
65      !!----------------------------------------------------------------------
66
67      CALL wrk_alloc( jpi, jpj, ztswu, ztswv, ztbfu, ztbfv, z2dx, z2dy )
68
69      z2dx(:,:) = 0._wp   ;   z2dy(:,:) = 0._wp                            ! initialization of workspace arrays
70
71      IF( neuler == 0 .AND. kt == nit000    ) THEN   ;   r2dt =      rdt   ! = rdtra (restart with Euler time stepping)
72      ELSEIF(               kt <= nit000 + 1) THEN   ;   r2dt = 2. * rdt   ! = 2 rdttra (leapfrog)
73      ENDIF
74
75      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
76      ! I. Integral Constraints Properties for momentum and/or tracers trends
77      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
78
79      IF( ( mod(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend) )   THEN
80         !
81         IF( lk_trdtra .AND. ctype == 'TRA' )   THEN       ! active tracer trends
82            SELECT CASE ( ktrd )
83            CASE ( jptra_trd_ldf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_ldf, ctype )   ! lateral diff
84            CASE ( jptra_trd_zdf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_zdf, ctype )   ! vertical diff (Kz)
85            CASE ( jptra_trd_bbc )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_bbc, ctype )   ! bottom boundary cond
86            CASE ( jptra_trd_bbl )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_bbl, ctype )   ! bottom boundary layer
87            CASE ( jptra_trd_npc )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_npc, ctype )   ! static instability mixing
88            CASE ( jptra_trd_dmp )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_dmp, ctype )   ! damping
89            CASE ( jptra_trd_qsr )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_qsr, ctype )   ! penetrative solar radiat.
90            CASE ( jptra_trd_nsr )   ;   z2dx(:,:) = ptrdx(:,:,1)   
91                                         z2dy(:,:) = ptrdy(:,:,1)
92                                         CALL trd_icp( z2dx , z2dy , jpicpt_nsr, ctype )   ! non solar radiation
93            CASE ( jptra_trd_xad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_xad, ctype )   ! x- horiz adv
94            CASE ( jptra_trd_yad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_yad, ctype )   ! y- horiz adv
95            CASE ( jptra_trd_zad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpt_zad, ctype )   ! z- vertical adv
96                                         CALL trd_icp( ptrdx, ptrdy, jpicpt_zad, ctype )   
97                                         ! compute the surface flux condition wn(:,:,1)*tsn(:,:,1,jp_tem)
98                                         z2dx(:,:) = wn(:,:,1)*tsn(:,:,1,jp_tem)/fse3t(:,:,1)
99                                         z2dy(:,:) = wn(:,:,1)*tsn(:,:,1,jp_sal)/fse3t(:,:,1)
100                                         CALL trd_icp( z2dx , z2dy , jpicpt_zl1, ctype )   ! 1st z- vertical adv
101            END SELECT
102         END IF
103
104         IF( lk_trddyn .AND. ctype == 'DYN' )   THEN       ! momentum trends
105            !
106            SELECT CASE ( ktrd )
107            CASE ( jpdyn_trd_hpg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_hpg, ctype )   ! hydrost. pressure grad
108            CASE ( jpdyn_trd_keg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_keg, ctype )   ! KE gradient
109            CASE ( jpdyn_trd_rvo )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_rvo, ctype )   ! relative vorticity
110            CASE ( jpdyn_trd_pvo )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_pvo, ctype )   ! planetary vorticity
111            CASE ( jpdyn_trd_ldf )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_ldf, ctype )   ! lateral diffusion
112            CASE ( jpdyn_trd_had )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_had, ctype )   ! horizontal advection
113            CASE ( jpdyn_trd_zad )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_zad, ctype )   ! vertical advection
114            CASE ( jpdyn_trd_spg )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_spg, ctype )   ! surface pressure grad.
115            CASE ( jpdyn_trd_dat )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_dat, ctype )   ! damping term
116            CASE ( jpdyn_trd_zdf )                                                         ! vertical diffusion
117               ! subtract surface forcing/bottom friction trends
118               ! from vertical diffusive momentum trends
119               ztswu(:,:) = 0._wp   ;   ztswv(:,:) = 0._wp
120               ztbfu(:,:) = 0._wp   ;   ztbfv(:,:) = 0._wp 
121               DO jj = 2, jpjm1   
122                  DO ji = fs_2, fs_jpim1   ! vector opt.
123                     ! save the surface forcing momentum fluxes
124                     ztswu(ji,jj) = utau(ji,jj) / ( fse3u(ji,jj,1)*rau0 )
125                     ztswv(ji,jj) = vtau(ji,jj) / ( fse3v(ji,jj,1)*rau0 )
126                     ! bottom friction contribution now handled explicitly
127                     ptrdx(ji,jj,1) = ptrdx(ji,jj,1) - ztswu(ji,jj)
128                     ptrdy(ji,jj,1) = ptrdy(ji,jj,1) - ztswv(ji,jj)
129                  END DO
130               END DO
131               !
132               CALL trd_icp( ptrdx, ptrdy, jpicpd_zdf, ctype )   
133               CALL trd_icp( ztswu, ztswv, jpicpd_swf, ctype )                               ! wind stress forcing term
134               ! bottom friction contribution now handled explicitly
135            CASE ( jpdyn_trd_bfr )   ;   CALL trd_icp( ptrdx, ptrdy, jpicpd_bfr, ctype )     ! bottom friction term
136            END SELECT
137            !
138         END IF
139         !
140      END IF
141
142      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
143      ! II. Vorticity trends
144      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
145
146      IF( lk_trdvor .AND. ctype == 'DYN' )   THEN
147         !
148         SELECT CASE ( ktrd ) 
149         CASE ( jpdyn_trd_hpg )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_prg )   ! Hydrostatique Pressure Gradient
150         CASE ( jpdyn_trd_keg )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_keg )   ! KE Gradient
151         CASE ( jpdyn_trd_rvo )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_rvo )   ! Relative Vorticity
152         CASE ( jpdyn_trd_pvo )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_pvo )   ! Planetary Vorticity Term
153         CASE ( jpdyn_trd_ldf )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_ldf )   ! Horizontal Diffusion
154         CASE ( jpdyn_trd_had )   ;   CALL ctl_warn('Vorticity for horizontal advection trend never checked')   
155         CASE ( jpdyn_trd_zad )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_zad )   ! Vertical Advection
156         CASE ( jpdyn_trd_spg )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_spg )   ! Surface Pressure Grad.
157         CASE ( jpdyn_trd_dat )   ;   CALL trd_vor_zint( ptrdx, ptrdy, jpvor_bev )   ! Beta V 
158         CASE ( jpdyn_trd_zdf )                                                      ! Vertical Diffusion
159            ! subtract surface forcing/bottom friction trends
160            ! from vertical diffusive momentum trends
161            ztswu(:,:) = 0.e0   ;   ztswv(:,:) = 0.e0
162            ztbfu(:,:) = 0.e0   ;   ztbfv(:,:) = 0.e0 
163            DO jj = 2, jpjm1   
164               DO ji = fs_2, fs_jpim1   ! vector opt.
165                  ! save the surface forcing momentum fluxes
166                  ztswu(ji,jj) = utau(ji,jj) / ( fse3u(ji,jj,1)*rau0 )
167                  ztswv(ji,jj) = vtau(ji,jj) / ( fse3v(ji,jj,1)*rau0 )
168                  !
169                  ptrdx(ji,jj,1     ) = ptrdx(ji,jj,1     ) - ztswu(ji,jj)
170                  ptrdy(ji,jj,1     ) = ptrdy(ji,jj,1     ) - ztswv(ji,jj)
171               END DO
172            END DO
173            !
174            CALL trd_vor_zint( ptrdx, ptrdy, jpvor_zdf )   
175            CALL trd_vor_zint( ztswu, ztswv, jpvor_swf )                               ! Wind stress forcing term
176         CASE ( jpdyn_trd_bfr )
177            CALL trd_vor_zint( ptrdx, ptrdy, jpvor_bfr )                               ! Bottom friction term
178         END SELECT
179         !
180      ENDIF
181
182      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
183      ! III. Mixed layer trends for active tracers
184      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
185
186      IF( lk_trdmld .AND. ctype == 'TRA' )   THEN
187         
188         !-----------------------------------------------------------------------------------------------
189         ! W.A.R.N.I.N.G :
190         ! jptra_trd_ldf : called by traldf.F90
191         !                 at this stage we store:
192         !                  - the lateral geopotential diffusion (here, lateral = horizontal)
193         !                  - and the iso-neutral diffusion if activated
194         ! jptra_trd_zdf : called by trazdf.F90
195         !                 * in case of iso-neutral diffusion we store the vertical diffusion component in the
196         !                   lateral trend including the K_z contrib, which will be removed later (see trd_mld)
197         !-----------------------------------------------------------------------------------------------
198
199         SELECT CASE ( ktrd )
200         CASE ( jptra_trd_xad )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_xad, '3D' )   ! merid. advection
201         CASE ( jptra_trd_yad )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_yad, '3D' )   ! zonal  advection
202         CASE ( jptra_trd_zad )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_zad, '3D' )   ! vertical advection
203         CASE ( jptra_trd_ldf )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_ldf, '3D' )   ! lateral diffusive
204         CASE ( jptra_trd_bbl )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_bbl, '3D' )   ! bottom boundary layer
205         CASE ( jptra_trd_zdf )
206            IF( ln_traldf_iso )   THEN
207               CALL trd_mld_zint( ptrdx, ptrdy, jpmld_ldf, '3D' )   ! vertical diffusion (K_z)
208            ELSE
209               CALL trd_mld_zint( ptrdx, ptrdy, jpmld_zdf, '3D' )   ! vertical diffusion (K_z)
210            ENDIF
211         CASE ( jptra_trd_dmp )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_dmp, '3D' )   ! internal 3D restoring (tradmp)
212         CASE ( jptra_trd_qsr )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_for, '3D' )   ! air-sea : penetrative sol radiat
213         CASE ( jptra_trd_nsr )
214            ptrdx(:,:,2:jpk) = 0.e0   ;   ptrdy(:,:,2:jpk) = 0.e0
215            CALL trd_mld_zint( ptrdx, ptrdy, jpmld_for, '2D' )                             ! air-sea : non penetr sol radiat
216         CASE ( jptra_trd_bbc )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_bbc, '3D' )   ! bottom bound cond (geoth flux)
217         CASE ( jptra_trd_atf )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_atf, '3D' )   ! asselin numerical
218         CASE ( jptra_trd_npc )   ;   CALL trd_mld_zint( ptrdx, ptrdy, jpmld_npc, '3D' )   ! non penetr convect adjustment
219         END SELECT
220
221      ENDIF
222      !
223      CALL wrk_dealloc( jpi, jpj, ztswu, ztswv, ztbfu, ztbfv, z2dx, z2dy )
224      !
225   END SUBROUTINE trd_mod
226
227#else
228   !!----------------------------------------------------------------------
229   !!   Default case :                                         Empty module
230   !!----------------------------------------------------------------------
231   USE trdmod_oce      ! ocean variables trends
232   USE trdvor          ! ocean vorticity trends
233   USE trdicp          ! ocean bassin integral constraints properties
234   USE trdmld          ! ocean active mixed layer tracers trends
235   !!----------------------------------------------------------------------
236CONTAINS
237   SUBROUTINE trd_mod(ptrd3dx, ptrd3dy, ktrd , ctype, kt )   ! Empty routine
238      REAL(wp) ::   ptrd3dx(:,:,:), ptrd3dy(:,:,:)
239      INTEGER  ::   ktrd, kt                           
240      CHARACTER(len=3) ::  ctype                 
241      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', ptrd3dx(1,1,1), ptrd3dy(1,1,1)
242      WRITE(*,*) ' "   ": You should not have seen this print! error ?', ktrd, ctype, kt
243   END SUBROUTINE trd_mod
244#endif
245
246   SUBROUTINE trd_mod_init
247      !!----------------------------------------------------------------------
248      !!                  ***  ROUTINE trd_mod_init  ***
249      !!
250      !! ** Purpose :   Initialization of activated trends
251      !!----------------------------------------------------------------------
252      USE in_out_manager          ! I/O manager
253      USE lib_mpp                 ! MPP library
254      !!   
255      NAMELIST/namtrd/ nn_trd, nn_ctls, cn_trdrst_in, cn_trdrst_out, ln_trdmld_restart, rn_ucf, ln_trdmld_instant
256      INTEGER  ::   ios           ! Local integer output status for namelist read
257      !!----------------------------------------------------------------------
258
259      IF( l_trdtra .OR. l_trddyn )   THEN
260 
261         REWIND( numnam_ref )              ! Namelist namtrd in reference namelist : Diagnostics: trends
262         READ  ( numnam_ref, namtrd, IOSTAT = ios, ERR = 901)
263901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrd in reference namelist', lwp )
264
265         REWIND( numnam_cfg )              ! Namelist namtrd in configuration namelist : Diagnostics: trends
266         READ  ( numnam_cfg, namtrd, IOSTAT = ios, ERR = 902 )
267902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrd in configuration namelist', lwp )
268         IF(lwm) WRITE ( numond, namtrd )
269
270         IF(lwp) THEN
271            WRITE(numout,*)
272            WRITE(numout,*) ' trd_mod_init : Momentum/Tracers trends'
273            WRITE(numout,*) ' ~~~~~~~~~~~~~'
274            WRITE(numout,*) '   Namelist namtrd : set trends parameters'
275            WRITE(numout,*) '      frequency of trends diagnostics   nn_trd             = ', nn_trd
276            WRITE(numout,*) '      control surface type              nn_ctls            = ', nn_ctls
277            WRITE(numout,*) '      restart for ML diagnostics        ln_trdmld_restart  = ', ln_trdmld_restart
278            WRITE(numout,*) '      instantaneous or mean ML T/S      ln_trdmld_instant  = ', ln_trdmld_instant
279            WRITE(numout,*) '      unit conversion factor            rn_ucf             = ', rn_ucf
280        ENDIF
281      ENDIF
282      !
283      IF( lk_trddyn .OR. lk_trdtra )    CALL trd_icp_init       ! integral constraints trends
284      IF( lk_trdmld                )    CALL trd_mld_init       ! mixed-layer trends (active  tracers) 
285      IF( lk_trdvor                )    CALL trd_vor_init       ! vorticity trends       
286      !
287   END SUBROUTINE trd_mod_init
288
289   !!======================================================================
290END MODULE trdmod
Note: See TracBrowser for help on using the repository browser.