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.
dia25h.F90 in branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/2017/dev_merge_2017/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90 @ 9019

Last change on this file since 9019 was 9019, checked in by timgraham, 6 years ago

Merge of dev_CNRS_2017 into branch

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