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 NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/DIA – NEMO

source: NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/DIA/dia25h.F90 @ 11671

Last change on this file since 11671 was 11671, checked in by acc, 5 years ago

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