source: branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90 @ 7953

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

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

File size: 11.9 KB
Line 
1MODULE dia25h 
2   !!======================================================================
3   !!                       ***  MODULE  diaharm  ***
4   !! Harmonic analysis of tidal constituents
5   !!======================================================================
6   !! History :  3.6  !  2014  (E O'Dea)  Original code
7   !!----------------------------------------------------------------------
8   USE oce             ! ocean dynamics and tracers variables
9   USE dom_oce         ! ocean space and time domain
10   USE zdf_oce         ! ocean vertical physics   
11   USE zdfgls, ONLY: mxln
12   USE in_out_manager  ! I/O units
13   USE iom             ! I/0 library
14   USE wrk_nemo        ! work arrays
15
16   IMPLICIT NONE
17   PRIVATE
18
19   PUBLIC   dia_25h_init               ! routine called by nemogcm.F90
20   PUBLIC   dia_25h                    ! routine called by diawri.F90
21
22   LOGICAL, PUBLIC ::   ln_dia25h      !:  25h mean output
23
24   ! variables for calculating 25-hourly means
25   INTEGER , SAVE ::   cnt_25h     ! Counter for 25 hour means
26   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   tn_25h  , sn_25h
27   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:)   ::   sshn_25h 
28   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   un_25h  , vn_25h  , wn_25h
29   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   avt_25h , avm_25h
30   REAL(wp), SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   en_25h  , rmxln_25h
31
32   !!----------------------------------------------------------------------
33   !! NEMO/OPA 3.6 , NEMO Consortium (2014)
34   !! $Id:$
35   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
36   !!----------------------------------------------------------------------
37CONTAINS
38
39   SUBROUTINE dia_25h_init 
40      !!---------------------------------------------------------------------------
41      !!                  ***  ROUTINE dia_25h_init  ***
42      !!     
43      !! ** Purpose: Initialization of 25h mean namelist
44      !!       
45      !! ** Method : Read namelist
46      !!---------------------------------------------------------------------------
47      INTEGER ::   ios                 ! Local integer output status for namelist read
48      INTEGER ::   ierror              ! Local integer for memory allocation
49      !
50      NAMELIST/nam_dia25h/ ln_dia25h
51      !!----------------------------------------------------------------------
52      !
53      REWIND ( numnam_ref )              ! Read Namelist nam_dia25h in reference namelist : 25hour mean diagnostics
54      READ   ( numnam_ref, nam_dia25h, IOSTAT=ios, ERR= 901 )
55901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_dia25h in reference namelist', lwp )
56
57      REWIND( numnam_cfg )              ! Namelist nam_dia25h in configuration namelist  25hour diagnostics
58      READ  ( numnam_cfg, nam_dia25h, IOSTAT = ios, ERR = 902 )
59902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_dia25h in configuration namelist', lwp )
60      IF(lwm) WRITE ( numond, nam_dia25h )
61
62      IF(lwp) THEN                   ! Control print
63         WRITE(numout,*)
64         WRITE(numout,*) 'dia_25h_init : Output 25 hour mean diagnostics'
65         WRITE(numout,*) '~~~~~~~~~~~~'
66         WRITE(numout,*) '   Namelist nam_dia25h : set 25h outputs '
67         WRITE(numout,*) '      Switch for 25h diagnostics (T) or not (F)  ln_dia25h  = ', ln_dia25h
68      ENDIF
69      IF( .NOT. ln_dia25h )   RETURN
70      ! ------------------- !
71      ! 1 - Allocate memory !
72      ! ------------------- !
73      !                                ! ocean arrays
74      ALLOCATE( tn_25h (jpi,jpj,jpk), sn_25h (jpi,jpj,jpk), sshn_25h(jpi,jpj)  ,     &
75         &      un_25h (jpi,jpj,jpk), vn_25h (jpi,jpj,jpk), wn_25h(jpi,jpj,jpk),     &
76         &      avt_25h(jpi,jpj,jpk), avm_25h(jpi,jpj,jpk),                      STAT=ierror )
77      IF( ierror > 0 ) THEN
78         CALL ctl_stop( 'dia_25h: unable to allocate ocean arrays' )   ;   RETURN
79      ENDIF
80      IF( ln_zdftke ) THEN             ! TKE physics
81         ALLOCATE( en_25h(jpi,jpj,jpk), STAT=ierror )
82         IF( ierror > 0 ) THEN
83            CALL ctl_stop( 'dia_25h: unable to allocate en_25h' )   ;   RETURN
84         ENDIF
85      ENDIF
86      IF( ln_zdfgls ) THEN             ! GLS physics
87         ALLOCATE( en_25h(jpi,jpj,jpk), rmxln_25h(jpi,jpj,jpk), STAT=ierror )
88         IF( ierror > 0 ) THEN
89            CALL ctl_stop( 'dia_25h: unable to allocate en_25h and rmxln_25h' )   ;   RETURN
90         ENDIF
91      ENDIF
92      ! ------------------------- !
93      ! 2 - Assign Initial Values !
94      ! ------------------------- !
95      cnt_25h = 1  ! sets the first value of sum at timestep 1 (note - should strictly be at timestep zero so before values used where possible)
96      tn_25h(:,:,:) = tsb(:,:,:,jp_tem)
97      sn_25h(:,:,:) = tsb(:,:,:,jp_sal)
98      sshn_25h(:,:) = sshb(:,:)
99      un_25h(:,:,:) = ub(:,:,:)
100      vn_25h(:,:,:) = vb(:,:,:)
101      wn_25h(:,:,:) = wn(:,:,:)
102      avt_25h(:,:,:) = avt(:,:,:)
103      avm_25h(:,:,:) = avm(:,:,:)
104      IF( ln_zdftke ) THEN
105         en_25h(:,:,:) = en(:,:,:)
106      ENDIF
107      IF( ln_zdfgls ) THEN
108         en_25h(:,:,:) = en(:,:,:)
109         rmxln_25h(:,:,:) = mxln(:,:,:)
110      ENDIF
111#if defined key_lim3 || defined key_lim2
112         CALL ctl_stop('STOP', 'dia_25h not setup yet to do tidemean ice')
113#endif 
114      !
115   END SUBROUTINE dia_25h_init
116
117
118   SUBROUTINE dia_25h( kt ) 
119      !!----------------------------------------------------------------------
120      !!                 ***  ROUTINE dia_25h  ***
121      !!         
122      !! ** Purpose :   Write diagnostics with M2/S2 tide removed
123      !!
124      !! ** Method  :   25hr mean outputs for shelf seas
125      !!----------------------------------------------------------------------
126      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
127      !!
128      INTEGER ::   ji, jj, jk
129      INTEGER                          ::   iyear0, nimonth0,iday0            ! start year,imonth,day
130      LOGICAL ::   ll_print = .FALSE.    ! =T print and flush numout
131      REAL(wp)                         ::   zsto, zout, zmax, zjulian, zmdi   ! local scalars
132      INTEGER                          ::   i_steps                           ! no of timesteps per hour
133      REAL(wp), DIMENSION(jpi,jpj    ) ::   zw2d, un_dm, vn_dm                ! workspace
134      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zw3d                              ! workspace
135      REAL(wp), DIMENSION(jpi,jpj,3)   ::   zwtmb                             ! workspace
136      !!----------------------------------------------------------------------
137
138      ! 0. Initialisation
139      ! -----------------
140      ! Define frequency of summing to create 25 h mean
141      IF( MOD( 3600,INT(rdt) ) == 0 ) THEN
142         i_steps = 3600/INT(rdt)
143      ELSE
144         CALL ctl_stop('STOP', 'dia_wri_tide: timestep must give MOD(3600,rdt) = 0 otherwise no hourly values are possible')
145      ENDIF
146
147      ! local variable for debugging
148      ll_print = ll_print .AND. lwp
149
150      ! Sum of 25 hourly instantaneous values to give a 25h mean from 24hours every day
151      IF( MOD( kt, i_steps ) == 0  .AND. kt /= nn_it000 ) THEN
152
153         IF (lwp) THEN
154              WRITE(numout,*) 'dia_wri_tide : Summing instantaneous hourly diagnostics at timestep ',kt
155              WRITE(numout,*) '~~~~~~~~~~~~ '
156         ENDIF
157
158         tn_25h(:,:,:)        = tn_25h(:,:,:) + tsn(:,:,:,jp_tem)
159         sn_25h(:,:,:)        = sn_25h(:,:,:) + tsn(:,:,:,jp_sal)
160         sshn_25h(:,:)        = sshn_25h(:,:) + sshn (:,:)
161         un_25h(:,:,:)        = un_25h(:,:,:) + un(:,:,:)
162         vn_25h(:,:,:)        = vn_25h(:,:,:) + vn(:,:,:)
163         wn_25h(:,:,:)        = wn_25h(:,:,:) + wn(:,:,:)
164         avt_25h(:,:,:)       = avt_25h(:,:,:) + avt(:,:,:)
165         avm_25h(:,:,:)       = avm_25h(:,:,:) + avm(:,:,:)
166         IF( ln_zdftke ) THEN
167            en_25h(:,:,:)        = en_25h(:,:,:) + en(:,:,:)
168         ENDIF
169         IF( ln_zdfgls ) THEN
170            en_25h(:,:,:)        = en_25h(:,:,:) + en(:,:,:)
171            rmxln_25h(:,:,:)      = rmxln_25h(:,:,:) + mxln(:,:,:)
172         ENDIF
173         cnt_25h = cnt_25h + 1
174         !
175         IF (lwp) THEN
176            WRITE(numout,*) 'dia_tide : Summed the following number of hourly values so far',cnt_25h
177         ENDIF
178         !
179      ENDIF ! MOD( kt, i_steps ) == 0
180
181      ! Write data for 25 hour mean output streams
182      IF( cnt_25h == 25 .AND.  MOD( kt, i_steps*24) == 0 .AND. kt /= nn_it000 ) THEN
183         !
184         IF(lwp) THEN
185            WRITE(numout,*) 'dia_wri_tide : Writing 25 hour mean tide diagnostics at timestep', kt
186            WRITE(numout,*) '~~~~~~~~~~~~ '
187         ENDIF
188         !
189         tn_25h(:,:,:)        = tn_25h(:,:,:) / 25.0_wp
190         sn_25h(:,:,:)        = sn_25h(:,:,:) / 25.0_wp
191         sshn_25h(:,:)        = sshn_25h(:,:) / 25.0_wp
192         un_25h(:,:,:)        = un_25h(:,:,:) / 25.0_wp
193         vn_25h(:,:,:)        = vn_25h(:,:,:) / 25.0_wp
194         wn_25h(:,:,:)        = wn_25h(:,:,:) / 25.0_wp
195         avt_25h(:,:,:)       = avt_25h(:,:,:) / 25.0_wp
196         avm_25h(:,:,:)       = avm_25h(:,:,:) / 25.0_wp
197         IF( ln_zdftke ) THEN
198            en_25h(:,:,:)        = en_25h(:,:,:) / 25.0_wp
199         ENDIF
200         IF( ln_zdfgls ) THEN
201            en_25h(:,:,:)        = en_25h(:,:,:) / 25.0_wp
202            rmxln_25h(:,:,:)       = rmxln_25h(:,:,:) / 25.0_wp
203         ENDIF
204         !
205         IF(lwp)  WRITE(numout,*) 'dia_wri_tide : Mean calculated by dividing 25 hour sums and writing output'
206         zmdi=1.e+20 !missing data indicator for masking
207         ! write tracers (instantaneous)
208         zw3d(:,:,:) = tn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
209         CALL iom_put("temper25h", zw3d)   ! potential temperature
210         zw3d(:,:,:) = sn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
211         CALL iom_put( "salin25h", zw3d  )   ! salinity
212         zw2d(:,:) = sshn_25h(:,:)*tmask(:,:,1) + zmdi*(1.0-tmask(:,:,1))
213         CALL iom_put( "ssh25h", zw2d )   ! sea surface
214         ! Write velocities (instantaneous)
215         zw3d(:,:,:) = un_25h(:,:,:)*umask(:,:,:) + zmdi*(1.0-umask(:,:,:))
216         CALL iom_put("vozocrtx25h", zw3d)    ! i-current
217         zw3d(:,:,:) = vn_25h(:,:,:)*vmask(:,:,:) + zmdi*(1.0-vmask(:,:,:))
218         CALL iom_put("vomecrty25h", zw3d  )   ! j-current
219         zw3d(:,:,:) = wn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
220         CALL iom_put("vomecrtz25h", zw3d )   ! k-current
221         ! Write vertical physics
222         zw3d(:,:,:) = avt_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
223         CALL iom_put("avt25h", zw3d )   ! diffusivity
224         zw3d(:,:,:) = avm_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
225         CALL iom_put("avm25h", zw3d)   ! viscosity
226         IF( ln_zdftke ) THEN
227            zw3d(:,:,:) = en_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
228            CALL iom_put("tke25h", zw3d)   ! tke
229         ENDIF
230         IF( ln_zdfgls ) THEN
231            zw3d(:,:,:) = en_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
232            CALL iom_put("tke25h", zw3d)   ! tke
233            zw3d(:,:,:) = rmxln_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
234            CALL iom_put( "mxln25h",zw3d)
235         ENDIF
236         !
237         ! After the write reset the values to cnt=1 and sum values equal current value
238         tn_25h(:,:,:) = tsn(:,:,:,jp_tem)
239         sn_25h(:,:,:) = tsn(:,:,:,jp_sal)
240         sshn_25h(:,:) = sshn (:,:)
241         un_25h(:,:,:) = un(:,:,:)
242         vn_25h(:,:,:) = vn(:,:,:)
243         wn_25h(:,:,:) = wn(:,:,:)
244         avt_25h(:,:,:) = avt(:,:,:)
245         avm_25h(:,:,:) = avm(:,:,:)
246         IF( ln_zdftke ) THEN
247            en_25h(:,:,:) = en(:,:,:)
248         ENDIF
249         IF( ln_zdfgls ) THEN
250            en_25h(:,:,:) = en(:,:,:)
251            rmxln_25h(:,:,:) = mxln(:,:,:)
252         ENDIF
253         cnt_25h = 1
254         IF (lwp)  WRITE(numout,*) 'dia_wri_tide : After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average',cnt_25h
255         !
256      ENDIF !  cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps * 24) == 0 .AND. kt .NE. nn_it000
257      !
258   END SUBROUTINE dia_25h 
259
260   !!======================================================================
261END MODULE dia25h
Note: See TracBrowser for help on using the repository browser.