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.
trdvor.F90 in NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/TRD – NEMO

source: NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/TRD/trdvor.F90 @ 12340

Last change on this file since 12340 was 12340, checked in by acc, 4 years ago

Branch 2019/dev_r11943_MERGE_2019. This commit introduces basic do loop macro
substitution to the 2019 option 1, merge branch. These changes have been SETTE
tested. The only addition is the do_loop_substitute.h90 file in the OCE directory but
the macros defined therein are used throughout the code to replace identifiable, 2D-
and 3D- nested loop opening and closing statements with single-line alternatives. Code
indents are also adjusted accordingly.

The following explanation is taken from comments in the new header file:

This header file contains preprocessor definitions and macros used in the do-loop
substitutions introduced between version 4.0 and 4.2. The primary aim of these macros
is to assist in future applications of tiling to improve performance. This is expected
to be achieved by alternative versions of these macros in selected locations. The
initial introduction of these macros simply replaces all identifiable nested 2D- and
3D-loops with single line statements (and adjusts indenting accordingly). Do loops
are identifiable if they comform to either:

DO jk = ....

DO jj = .... DO jj = ...

DO ji = .... DO ji = ...
. OR .
. .

END DO END DO

END DO END DO

END DO

and white-space variants thereof.

Additionally, only loops with recognised jj and ji loops limits are treated; these are:
Lower limits of 1, 2 or fs_2
Upper limits of jpi, jpim1 or fs_jpim1 (for ji) or jpj, jpjm1 or fs_jpjm1 (for jj)

The macro naming convention takes the form: DO_2D_BT_LR where:

B is the Bottom offset from the PE's inner domain;
T is the Top offset from the PE's inner domain;
L is the Left offset from the PE's inner domain;
R is the Right offset from the PE's inner domain

So, given an inner domain of 2,jpim1 and 2,jpjm1, a typical example would replace:

DO jj = 2, jpj

DO ji = 1, jpim1
.
.

END DO

END DO

with:

DO_2D_01_10
.
.
END_2D

similar conventions apply to the 3D loops macros. jk loop limits are retained
through macro arguments and are not restricted. This includes the possibility of
strides for which an extra set of DO_3DS macros are defined.

In the example definition below the inner PE domain is defined by start indices of
(kIs, kJs) and end indices of (kIe, KJe)

#define DO_2D_00_00 DO jj = kJs, kJe ; DO ji = kIs, kIe
#define END_2D END DO ; END DO

TO DO:

Only conventional nested loops have been identified and replaced by this step. There are constructs such as:

DO jk = 2, jpkm1

z2d(:,:) = z2d(:,:) + e3w(:,:,jk,Kmm) * z3d(:,:,jk) * wmask(:,:,jk)

END DO

which may need to be considered.
  • Property svn:keywords set to Id
File size: 25.9 KB
Line 
1MODULE trdvor
2   !!======================================================================
3   !!                       ***  MODULE  trdvor  ***
4   !! Ocean diagnostics:  momentum trends
5   !!=====================================================================
6   !! History :  1.0  !  2006-01  (L. Brunier, A-M. Treguier) Original code
7   !!            2.0  !  2008-04  (C. Talandier) New trends organization
8   !!            3.5  !  2012-02  (G. Madec) regroup beta.V computation with pvo trend
9   !!----------------------------------------------------------------------
10
11   !!----------------------------------------------------------------------
12   !!   trd_vor      : momentum trends averaged over the depth
13   !!   trd_vor_zint : vorticity vertical integration
14   !!   trd_vor_init : initialization step
15   !!----------------------------------------------------------------------
16   USE oce             ! ocean dynamics and tracers variables
17   USE dom_oce         ! ocean space and time domain variables
18   USE trd_oce         ! trends: ocean variables
19   USE zdf_oce         ! ocean vertical physics
20   USE sbc_oce         ! surface boundary condition: ocean
21   USE phycst          ! Define parameters for the routines
22   USE ldfdyn          ! ocean active tracers: lateral physics
23   USE dianam          ! build the name of file (routine)
24   USE zdfmxl          ! mixed layer depth
25   !
26   USE in_out_manager  ! I/O manager
27   USE ioipsl          ! NetCDF library
28   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
29   USE lib_mpp         ! MPP library
30
31   IMPLICIT NONE
32   PRIVATE
33
34   INTERFACE trd_vor_zint
35      MODULE PROCEDURE trd_vor_zint_2d, trd_vor_zint_3d
36   END INTERFACE
37
38   PUBLIC   trd_vor        ! routine called by trddyn.F90
39   PUBLIC   trd_vor_init   ! routine called by opa.F90
40   PUBLIC   trd_vor_alloc  ! routine called by nemogcm.F90
41
42   INTEGER ::   nh_t, nmoydpvor, nidvor, nhoridvor, ndimvor1, icount   ! needs for IOIPSL output
43   INTEGER, SAVE, ALLOCATABLE, DIMENSION(:) ::   ndexvor1   ! needed for IOIPSL output
44   INTEGER ::   ndebug     ! (0/1) set it to 1 in case of problem to have more print
45
46   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   vor_avr      ! average
47   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   vor_avrb     ! before vorticity (kt-1)
48   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   vor_avrbb    ! vorticity at begining of the nn_write-1 timestep averaging period
49   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   vor_avrbn    ! after vorticity at time step after the
50   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   rotot        ! begining of the NN_WRITE-1 timesteps
51   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   vor_avrtot   !
52   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:)   ::   vor_avrres   !
53   REAL(wp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   vortrd       ! curl of trends
54         
55   CHARACTER(len=12) ::   cvort
56
57   !! * Substitutions
58#  include "vectopt_loop_substitute.h90"
59#  include "do_loop_substitute.h90"
60   !!----------------------------------------------------------------------
61   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
62   !! $Id$
63   !! Software governed by the CeCILL license (see ./LICENSE)
64   !!----------------------------------------------------------------------
65CONTAINS
66
67   INTEGER FUNCTION trd_vor_alloc()
68      !!----------------------------------------------------------------------------
69      !!                  ***  ROUTINE trd_vor_alloc  ***
70      !!----------------------------------------------------------------------------
71      ALLOCATE( vor_avr   (jpi,jpj) , vor_avrb(jpi,jpj) , vor_avrbb (jpi,jpj) ,   &
72         &      vor_avrbn (jpi,jpj) , rotot   (jpi,jpj) , vor_avrtot(jpi,jpj) ,   &
73         &      vor_avrres(jpi,jpj) , vortrd  (jpi,jpj,jpltot_vor) ,              &
74         &      ndexvor1  (jpi*jpj)                                ,   STAT= trd_vor_alloc )
75         !
76      CALL mpp_sum ( 'trdvor', trd_vor_alloc )
77      IF( trd_vor_alloc /= 0 )   CALL ctl_stop( 'STOP', 'trd_vor_alloc: failed to allocate arrays' )
78   END FUNCTION trd_vor_alloc
79
80
81   SUBROUTINE trd_vor( putrd, pvtrd, ktrd, kt, Kmm )
82      !!----------------------------------------------------------------------
83      !!                  ***  ROUTINE trd_vor  ***
84      !!
85      !! ** Purpose :  computation of cumulated trends over analysis period
86      !!               and make outputs (NetCDF format)
87      !!----------------------------------------------------------------------
88      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   putrd, pvtrd   ! U and V trends
89      INTEGER                   , INTENT(in   ) ::   ktrd           ! trend index
90      INTEGER                   , INTENT(in   ) ::   kt             ! time step
91      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
92      !
93      INTEGER ::   ji, jj   ! dummy loop indices
94      REAL(wp), DIMENSION(jpi,jpj) ::   ztswu, ztswv    ! 2D workspace
95      !!----------------------------------------------------------------------
96
97      SELECT CASE( ktrd ) 
98      CASE( jpdyn_hpg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_prg, Kmm )   ! Hydrostatique Pressure Gradient
99      CASE( jpdyn_keg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_keg, Kmm )   ! KE Gradient
100      CASE( jpdyn_rvo )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_rvo, Kmm )   ! Relative Vorticity
101      CASE( jpdyn_pvo )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_pvo, Kmm )   ! Planetary Vorticity Term
102      CASE( jpdyn_ldf )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_ldf, Kmm )   ! Horizontal Diffusion
103      CASE( jpdyn_zad )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_zad, Kmm )   ! Vertical Advection
104      CASE( jpdyn_spg )   ;   CALL trd_vor_zint( putrd, pvtrd, jpvor_spg, Kmm )   ! Surface Pressure Grad.
105      CASE( jpdyn_zdf )                                                      ! Vertical Diffusion
106         ztswu(:,:) = 0.e0   ;   ztswv(:,:) = 0.e0
107         DO_2D_00_00
108            ztswu(ji,jj) = 0.5 * ( utau_b(ji,jj) + utau(ji,jj) ) / ( e3u(ji,jj,1,Kmm) * rau0 )
109            ztswv(ji,jj) = 0.5 * ( vtau_b(ji,jj) + vtau(ji,jj) ) / ( e3v(ji,jj,1,Kmm) * rau0 )
110         END_2D
111         !
112         CALL trd_vor_zint( putrd, pvtrd, jpvor_zdf, Kmm )                             ! zdf trend including surf./bot. stresses
113         CALL trd_vor_zint( ztswu, ztswv, jpvor_swf, Kmm )                             ! surface wind stress
114      CASE( jpdyn_bfr )
115         CALL trd_vor_zint( putrd, pvtrd, jpvor_bfr, Kmm )                             ! Bottom stress
116         !
117      CASE( jpdyn_atf )       ! last trends: perform the output of 2D vorticity trends
118         CALL trd_vor_iom( kt, Kmm )
119      END SELECT
120      !
121   END SUBROUTINE trd_vor
122
123
124   SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd, Kmm )
125      !!----------------------------------------------------------------------------
126      !!                  ***  ROUTINE trd_vor_zint  ***
127      !!
128      !! ** Purpose :   computation of vertically integrated vorticity budgets
129      !!              from ocean surface down to control surface (NetCDF output)
130      !!
131      !! ** Method/usage :   integration done over nn_write-1 time steps
132      !!
133      !! ** Action :   trends :
134      !!                  vortrd (,, 1) = Pressure Gradient Trend
135      !!                  vortrd (,, 2) = KE Gradient Trend
136      !!                  vortrd (,, 3) = Relative Vorticity Trend
137      !!                  vortrd (,, 4) = Coriolis Term Trend
138      !!                  vortrd (,, 5) = Horizontal Diffusion Trend
139      !!                  vortrd (,, 6) = Vertical Advection Trend
140      !!                  vortrd (,, 7) = Vertical Diffusion Trend
141      !!                  vortrd (,, 8) = Surface Pressure Grad. Trend
142      !!                  vortrd (,, 9) = Beta V
143      !!                  vortrd (,,10) = forcing term
144      !!                  vortrd (,,11) = bottom friction term
145      !!                  rotot(,) : total cumulative trends over nn_write-1 time steps
146      !!                  vor_avrtot(,) : first membre of vrticity equation
147      !!                  vor_avrres(,) : residual = dh/dt entrainment
148      !!
149      !!      trends output in netCDF format using ioipsl
150      !!----------------------------------------------------------------------
151      INTEGER                     , INTENT(in   ) ::   ktrd       ! ocean trend index
152      INTEGER                     , INTENT(in   ) ::   Kmm        ! time level index
153      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   putrdvor   ! u vorticity trend
154      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   pvtrdvor   ! v vorticity trend
155      !
156      INTEGER ::   ji, jj       ! dummy loop indices
157      INTEGER ::   ikbu, ikbv   ! local integers
158      REAL(wp), DIMENSION(jpi,jpj) :: zudpvor, zvdpvor  ! total cmulative trends
159      !!----------------------------------------------------------------------
160
161      !
162
163      zudpvor(:,:) = 0._wp                 ;   zvdpvor(:,:) = 0._wp                    ! Initialisation
164      CALL lbc_lnk_multi( 'trdvor', putrdvor, 'U', -1. , pvtrdvor, 'V', -1. )      ! lateral boundary condition
165     
166
167      !  =====================================
168      !  I vertical integration of 2D trends
169      !  =====================================
170
171      SELECT CASE( ktrd ) 
172      !
173      CASE( jpvor_bfr )        ! bottom friction
174         DO_2D_00_00
175            ikbu = mbkv(ji,jj)
176            ikbv = mbkv(ji,jj)           
177            zudpvor(ji,jj) = putrdvor(ji,jj) * e3u(ji,jj,ikbu,Kmm) * e1u(ji,jj) * umask(ji,jj,ikbu)
178            zvdpvor(ji,jj) = pvtrdvor(ji,jj) * e3v(ji,jj,ikbv,Kmm) * e2v(ji,jj) * vmask(ji,jj,ikbv)
179         END_2D
180         !
181      CASE( jpvor_swf )        ! wind stress
182         zudpvor(:,:) = putrdvor(:,:) * e3u(:,:,1,Kmm) * e1u(:,:) * umask(:,:,1)
183         zvdpvor(:,:) = pvtrdvor(:,:) * e3v(:,:,1,Kmm) * e2v(:,:) * vmask(:,:,1)
184         !
185      END SELECT
186
187      ! Average except for Beta.V
188      zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm)
189      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm)
190   
191      ! Curl
192      DO ji = 1, jpim1
193         DO jj = 1, jpjm1
194            vortrd(ji,jj,ktrd) = (    zvdpvor(ji+1,jj) - zvdpvor(ji,jj)       &
195                 &                - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) )   ) / ( e1f(ji,jj) * e2f(ji,jj) )
196         END DO
197      END DO
198      vortrd(:,:,ktrd) = vortrd(:,:,ktrd) * fmask(:,:,1)      ! Surface mask
199
200      IF( ndebug /= 0 ) THEN
201         IF(lwp) WRITE(numout,*) ' debuging trd_vor_zint: I done'
202         CALL FLUSH(numout)
203      ENDIF
204      !
205   END SUBROUTINE trd_vor_zint_2d
206
207
208   SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd , Kmm )
209      !!----------------------------------------------------------------------------
210      !!                  ***  ROUTINE trd_vor_zint  ***
211      !!
212      !! ** Purpose :   computation of vertically integrated vorticity budgets
213      !!              from ocean surface down to control surface (NetCDF output)
214      !!
215      !! ** Method/usage :   integration done over nn_write-1 time steps
216      !!
217      !! ** Action :     trends :
218      !!                  vortrd (,,1) = Pressure Gradient Trend
219      !!                  vortrd (,,2) = KE Gradient Trend
220      !!                  vortrd (,,3) = Relative Vorticity Trend
221      !!                  vortrd (,,4) = Coriolis Term Trend
222      !!                  vortrd (,,5) = Horizontal Diffusion Trend
223      !!                  vortrd (,,6) = Vertical Advection Trend
224      !!                  vortrd (,,7) = Vertical Diffusion Trend
225      !!                  vortrd (,,8) = Surface Pressure Grad. Trend
226      !!                  vortrd (,,9) = Beta V
227      !!                  vortrd (,,10) = forcing term
228      !!      vortrd (,,11) = bottom friction term
229      !!                  rotot(,) : total cumulative trends over nn_write-1 time steps
230      !!                  vor_avrtot(,) : first membre of vrticity equation
231      !!                  vor_avrres(,) : residual = dh/dt entrainment
232      !!
233      !!      trends output in netCDF format using ioipsl
234      !!----------------------------------------------------------------------
235      !
236      INTEGER                         , INTENT(in   ) ::   ktrd       ! ocean trend index
237      INTEGER                         , INTENT(in   ) ::   Kmm        ! time level index
238      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   putrdvor   ! u vorticity trend
239      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pvtrdvor   ! v vorticity trend
240      !
241      INTEGER ::   ji, jj, jk   ! dummy loop indices
242      REAL(wp), DIMENSION(jpi,jpj) :: zubet  , zvbet    ! Beta.V   
243      REAL(wp), DIMENSION(jpi,jpj) :: zudpvor, zvdpvor  ! total cmulative trends
244      !!----------------------------------------------------------------------
245     
246      ! Initialization
247      zubet  (:,:) = 0._wp
248      zvbet  (:,:) = 0._wp
249      zudpvor(:,:) = 0._wp
250      zvdpvor(:,:) = 0._wp
251      !                            ! lateral boundary condition on input momentum trends
252      CALL lbc_lnk_multi( 'trdvor', putrdvor, 'U', -1. , pvtrdvor, 'V', -1. )
253
254      !  =====================================
255      !  I vertical integration of 3D trends
256      !  =====================================
257      ! putrdvor and pvtrdvor terms
258      DO jk = 1,jpk
259        zudpvor(:,:) = zudpvor(:,:) + putrdvor(:,:,jk) * e3u(:,:,jk,Kmm) * e1u(:,:) * umask(:,:,jk)
260        zvdpvor(:,:) = zvdpvor(:,:) + pvtrdvor(:,:,jk) * e3v(:,:,jk,Kmm) * e2v(:,:) * vmask(:,:,jk)
261      END DO
262
263      ! Planetary vorticity: 2nd computation (Beta.V term) store the vertical sum
264      ! as Beta.V term need intergration, not average
265      IF( ktrd == jpvor_pvo ) THEN
266         zubet(:,:) = zudpvor(:,:)
267         zvbet(:,:) = zvdpvor(:,:)
268         DO ji = 1, jpim1
269            DO jj = 1, jpjm1
270               vortrd(ji,jj,jpvor_bev) = (    zvbet(ji+1,jj) - zvbet(ji,jj)     &
271                  &                       - ( zubet(ji,jj+1) - zubet(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) )
272            END DO
273         END DO
274         ! Average of the Curl and Surface mask
275         vortrd(:,:,jpvor_bev) = vortrd(:,:,jpvor_bev) * r1_hu(:,:,Kmm) * fmask(:,:,1)
276      ENDIF
277      !
278      ! Average
279      zudpvor(:,:) = zudpvor(:,:) * r1_hu(:,:,Kmm)
280      zvdpvor(:,:) = zvdpvor(:,:) * r1_hv(:,:,Kmm)
281      !
282      ! Curl
283      DO ji=1,jpim1
284         DO jj=1,jpjm1
285            vortrd(ji,jj,ktrd) = (    zvdpvor(ji+1,jj) - zvdpvor(ji,jj)     &
286               &                  - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) )
287         END DO
288      END DO
289      ! Surface mask
290      vortrd(:,:,ktrd) = vortrd(:,:,ktrd) * fmask(:,:,1)
291   
292      IF( ndebug /= 0 ) THEN
293         IF(lwp) WRITE(numout,*) ' debuging trd_vor_zint: I done'
294         CALL FLUSH(numout)
295      ENDIF
296      !
297   END SUBROUTINE trd_vor_zint_3d
298
299
300   SUBROUTINE trd_vor_iom( kt , Kmm )
301      !!----------------------------------------------------------------------
302      !!                  ***  ROUTINE trd_vor  ***
303      !!
304      !! ** Purpose :  computation of cumulated trends over analysis period
305      !!               and make outputs (NetCDF format)
306      !!----------------------------------------------------------------------
307      INTEGER                   , INTENT(in   ) ::   kt             ! time step
308      INTEGER                   , INTENT(in   ) ::   Kmm            ! time level index
309      !
310      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
311      INTEGER  ::   it, itmod        ! local integers
312      REAL(wp) ::   zmean            ! local scalars
313      REAL(wp), DIMENSION(jpi,jpj) :: zuu, zvv
314      !!----------------------------------------------------------------------
315
316      !  =================
317      !  I. Initialization
318      !  =================
319     
320     
321      ! I.1 set before values of vertically average u and v
322      ! ---------------------------------------------------
323
324      IF( kt > nit000 )   vor_avrb(:,:) = vor_avr(:,:)
325
326      ! I.2 vertically integrated vorticity
327      !  ----------------------------------
328
329      vor_avr   (:,:) = 0._wp
330      zuu       (:,:) = 0._wp
331      zvv       (:,:) = 0._wp
332      vor_avrtot(:,:) = 0._wp
333      vor_avrres(:,:) = 0._wp
334     
335      ! Vertically averaged velocity
336      DO jk = 1, jpk - 1
337         zuu(:,:) = zuu(:,:) + e1u(:,:) * uu(:,:,jk,Kmm) * e3u(:,:,jk,Kmm)
338         zvv(:,:) = zvv(:,:) + e2v(:,:) * vv(:,:,jk,Kmm) * e3v(:,:,jk,Kmm)
339      END DO
340 
341      zuu(:,:) = zuu(:,:) * r1_hu(:,:,Kmm)
342      zvv(:,:) = zvv(:,:) * r1_hv(:,:,Kmm)
343
344      ! Curl
345      DO ji = 1, jpim1
346         DO jj = 1, jpjm1
347            vor_avr(ji,jj) = (  ( zvv(ji+1,jj) - zvv(ji,jj) )    &
348               &              - ( zuu(ji,jj+1) - zuu(ji,jj) ) ) / ( e1f(ji,jj) * e2f(ji,jj) ) * fmask(ji,jj,1)
349         END DO
350      END DO
351     
352      !  =================================
353      !   II. Cumulated trends
354      !  =================================
355
356      ! II.1 set `before' mixed layer values for kt = nit000+1
357      ! ------------------------------------------------------
358      IF( kt == nit000+1 ) THEN
359         vor_avrbb(:,:) = vor_avrb(:,:)
360         vor_avrbn(:,:) = vor_avr (:,:)
361      ENDIF
362
363      ! II.2 cumulated trends over analysis period (kt=2 to nn_write)
364      ! ----------------------
365      ! trends cumulated over nn_write-2 time steps
366
367      IF( kt >= nit000+2 ) THEN
368         nmoydpvor = nmoydpvor + 1
369         DO jl = 1, jpltot_vor
370            IF( jl /= 9 ) THEN
371               rotot(:,:) = rotot(:,:) + vortrd(:,:,jl)
372            ENDIF
373         END DO
374      ENDIF
375
376      !  =============================================
377      !   III. Output in netCDF + residual computation
378      !  =============================================
379     
380      ! define time axis
381      it    = kt
382      itmod = kt - nit000 + 1
383
384      IF( MOD( it, nn_trd ) == 0 ) THEN
385
386         ! III.1 compute total trend
387         ! ------------------------
388         zmean = 1._wp / (  REAL( nmoydpvor, wp ) * 2._wp * rdt  )
389         vor_avrtot(:,:) = (  vor_avr(:,:) - vor_avrbn(:,:) + vor_avrb(:,:) - vor_avrbb(:,:) ) * zmean
390
391
392         ! III.2 compute residual
393         ! ---------------------
394         zmean = 1._wp / REAL( nmoydpvor, wp )
395         vor_avrres(:,:) = vor_avrtot(:,:) - rotot(:,:) / zmean
396
397         ! Boundary conditions
398         CALL lbc_lnk_multi( 'trdvor', vor_avrtot, 'F', 1. , vor_avrres, 'F', 1. )
399
400
401         ! III.3 time evolution array swap
402         ! ------------------------------
403         vor_avrbb(:,:) = vor_avrb(:,:)
404         vor_avrbn(:,:) = vor_avr (:,:)
405         !
406         nmoydpvor = 0
407         !
408      ENDIF
409
410      ! III.4 write trends to output
411      ! ---------------------------
412
413      IF( kt >=  nit000+1 ) THEN
414
415         IF( lwp .AND. MOD( itmod, nn_trd ) == 0 ) THEN
416            WRITE(numout,*) ''
417            WRITE(numout,*) 'trd_vor : write trends in the NetCDF file at kt = ', kt
418            WRITE(numout,*) '~~~~~~~  '
419         ENDIF
420 
421         CALL histwrite( nidvor,"sovortPh",it,vortrd(:,:,jpvor_prg),ndimvor1,ndexvor1)  ! grad Ph
422         CALL histwrite( nidvor,"sovortEk",it,vortrd(:,:,jpvor_keg),ndimvor1,ndexvor1)  ! Energy
423         CALL histwrite( nidvor,"sovozeta",it,vortrd(:,:,jpvor_rvo),ndimvor1,ndexvor1)  ! rel vorticity
424         CALL histwrite( nidvor,"sovortif",it,vortrd(:,:,jpvor_pvo),ndimvor1,ndexvor1)  ! coriolis
425         CALL histwrite( nidvor,"sovodifl",it,vortrd(:,:,jpvor_ldf),ndimvor1,ndexvor1)  ! lat diff
426         CALL histwrite( nidvor,"sovoadvv",it,vortrd(:,:,jpvor_zad),ndimvor1,ndexvor1)  ! vert adv
427         CALL histwrite( nidvor,"sovodifv",it,vortrd(:,:,jpvor_zdf),ndimvor1,ndexvor1)  ! vert diff
428         CALL histwrite( nidvor,"sovortPs",it,vortrd(:,:,jpvor_spg),ndimvor1,ndexvor1)  ! grad Ps
429         CALL histwrite( nidvor,"sovortbv",it,vortrd(:,:,jpvor_bev),ndimvor1,ndexvor1)  ! beta.V
430         CALL histwrite( nidvor,"sovowind",it,vortrd(:,:,jpvor_swf),ndimvor1,ndexvor1) ! wind stress
431         CALL histwrite( nidvor,"sovobfri",it,vortrd(:,:,jpvor_bfr),ndimvor1,ndexvor1) ! bottom friction
432         CALL histwrite( nidvor,"1st_mbre",it,vor_avrtot    ,ndimvor1,ndexvor1) ! First membre
433         CALL histwrite( nidvor,"sovorgap",it,vor_avrres    ,ndimvor1,ndexvor1) ! gap between 1st and 2 nd mbre
434         !
435         IF( ndebug /= 0 ) THEN
436            WRITE(numout,*) ' debuging trd_vor: III.4 done'
437            CALL FLUSH(numout)
438         ENDIF
439         !
440      ENDIF
441      !
442      IF( MOD( it, nn_trd ) == 0 ) rotot(:,:)=0
443      !
444      IF( kt == nitend )   CALL histclo( nidvor )
445      !
446   END SUBROUTINE trd_vor_iom
447
448
449   SUBROUTINE trd_vor_init
450      !!----------------------------------------------------------------------
451      !!                  ***  ROUTINE trd_vor_init  ***
452      !!
453      !! ** Purpose :   computation of vertically integrated T and S budgets
454      !!      from ocean surface down to control surface (NetCDF output)
455      !!----------------------------------------------------------------------
456      REAL(wp) ::   zjulian, zsto, zout
457      CHARACTER (len=40) ::   clhstnam
458      CHARACTER (len=40) ::   clop
459      !!----------------------------------------------------------------------
460
461      !  ===================
462      !   I. initialization
463      !  ===================
464
465      cvort='averaged-vor'
466
467      ! Open specifier
468      ndebug = 0      ! set it to 1 in case of problem to have more Print
469
470      IF(lwp) THEN
471         WRITE(numout,*) ' '
472         WRITE(numout,*) ' trd_vor_init: vorticity trends'
473         WRITE(numout,*) ' ~~~~~~~~~~~~'
474         WRITE(numout,*) ' '
475         WRITE(numout,*) '               ##########################################################################'
476         WRITE(numout,*) '                CAUTION: The interpretation of the vorticity trends is'
477         WRITE(numout,*) '                not obvious, please contact Anne-Marie TREGUIER at: treguier@ifremer.fr '
478         WRITE(numout,*) '               ##########################################################################'
479         WRITE(numout,*) ' '
480      ENDIF
481
482      IF( trd_vor_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'trd_vor_init : unable to allocate trdvor arrays' )
483
484
485      ! cumulated trends array init
486      nmoydpvor = 0
487      rotot(:,:)=0
488      vor_avrtot(:,:)=0
489      vor_avrres(:,:)=0
490
491      IF( ndebug /= 0 ) THEN
492         WRITE(numout,*) ' debuging trd_vor_init: I. done'
493         CALL FLUSH(numout)
494      ENDIF
495
496      !  =================================
497      !   II. netCDF output initialization
498      !  =================================
499
500      !-----------------------------------------
501      ! II.1 Define frequency of output and means
502      ! -----------------------------------------
503      IF( ln_mskland )   THEN   ;   clop = "only(x)"   ! put 1.e+20 on land (very expensive!!)
504      ELSE                      ;   clop = "x"         ! no use of the mask value (require less cpu time)
505      ENDIF
506#if defined key_diainstant
507      zsto = nn_write*rdt
508      clop = "inst("//TRIM(clop)//")"
509#else
510      zsto = rdt
511      clop = "ave("//TRIM(clop)//")"
512#endif
513      zout = nn_trd*rdt
514
515      IF(lwp) WRITE(numout,*) '               netCDF initialization'
516
517      ! II.2 Compute julian date from starting date of the run
518      ! ------------------------
519      CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )
520      zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
521      IF(lwp) WRITE(numout,*)' ' 
522      IF(lwp) WRITE(numout,*)'               Date 0 used :',nit000,    &
523         &                   ' YEAR ', nyear,' MONTH '      , nmonth,   &
524         &                   ' DAY ' , nday, 'Julian day : ', zjulian
525
526      ! II.3 Define the T grid trend file (nidvor)
527      ! ---------------------------------
528      CALL dia_nam( clhstnam, nn_trd, 'vort' )                  ! filename
529      IF(lwp) WRITE(numout,*) ' Name of NETCDF file ', clhstnam
530      CALL histbeg( clhstnam, jpi, glamf, jpj, gphif,1, jpi,   &  ! Horizontal grid : glamt and gphit
531         &          1, jpj, nit000-1, zjulian, rdt, nh_t, nidvor, domain_id=nidom, snc4chunks=snc4set )
532      CALL wheneq( jpi*jpj, fmask, 1, 1., ndexvor1, ndimvor1 )    ! surface
533
534      ! Declare output fields as netCDF variables
535      CALL histdef( nidvor, "sovortPh", cvort//"grad Ph" , "s-2",        & ! grad Ph
536         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
537      CALL histdef( nidvor, "sovortEk", cvort//"Energy", "s-2",          & ! Energy
538         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
539      CALL histdef( nidvor, "sovozeta", cvort//"rel vorticity", "s-2",   & ! rel vorticity
540         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
541      CALL histdef( nidvor, "sovortif", cvort//"coriolis", "s-2",        & ! coriolis
542         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
543      CALL histdef( nidvor, "sovodifl", cvort//"lat diff ", "s-2",       & ! lat diff
544         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
545      CALL histdef( nidvor, "sovoadvv", cvort//"vert adv", "s-2",        & ! vert adv
546         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
547      CALL histdef( nidvor, "sovodifv", cvort//"vert diff" , "s-2",      & ! vert diff
548         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
549      CALL histdef( nidvor, "sovortPs", cvort//"grad Ps", "s-2",         & ! grad Ps
550         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
551      CALL histdef( nidvor, "sovortbv", cvort//"Beta V", "s-2",          & ! beta.V
552         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
553      CALL histdef( nidvor, "sovowind", cvort//"wind stress", "s-2",     & ! wind stress
554         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
555      CALL histdef( nidvor, "sovobfri", cvort//"bottom friction", "s-2", & ! bottom friction
556         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
557      CALL histdef( nidvor, "1st_mbre", cvort//"1st mbre", "s-2",        & ! First membre
558         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
559      CALL histdef( nidvor, "sovorgap", cvort//"gap", "s-2",             & ! gap between 1st and 2 nd mbre
560         &          jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout)
561      CALL histend( nidvor, snc4set )
562
563      IF( ndebug /= 0 ) THEN
564         WRITE(numout,*) ' debuging trd_vor_init: II. done'
565         CALL FLUSH(numout)
566      ENDIF
567      !
568   END SUBROUTINE trd_vor_init
569
570   !!======================================================================
571END MODULE trdvor
Note: See TracBrowser for help on using the repository browser.