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 branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/TRD – NEMO

source: branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/TRD/trdmod.F90 @ 2281

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

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