source: NEMO/branches/UKMO/dev_10448_WAD_SBC_BUGFIX/src/OCE/DIA/dia25h.F90 @ 10456

Last change on this file since 10456 was 10456, checked in by deazer, 22 months ago

Added option to taper sbc fluxes near very shallow water when using WAD
Corrected some IO bugs in dia25h, diatmb for WAD case.
User has control of the tapering. At what depth to start it, and at what fraction to start
the tanh tapering. At the WAD limit SBC is turned off completely.
Dry cells do not have any communication with the atmosphere
To DO: Documentation update.
Although not all sette tests are passed (AGRIF etc.)
it does no worse than the trunk at the revision the branch is made

  • Property svn:keywords set to Id
File size: 12.0 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      REWIND ( numnam_ref )              ! Read Namelist nam_dia25h in reference namelist : 25hour mean diagnostics
56      READ   ( numnam_ref, nam_dia25h, IOSTAT=ios, ERR= 901 )
57901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nam_dia25h in reference namelist', lwp )
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_si3
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,NINT(rdt) ) == 0 ) THEN
143         i_steps = 3600/NINT(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         IF( ll_wd ) THEN
215            CALL iom_put( "ssh25h", zw2d+ssh_ref )   ! sea surface
216         ELSE
217            CALL iom_put( "ssh25h", zw2d )   ! sea surface
218         ENDIF
219         ! Write velocities (instantaneous)
220         zw3d(:,:,:) = un_25h(:,:,:)*umask(:,:,:) + zmdi*(1.0-umask(:,:,:))
221         CALL iom_put("vozocrtx25h", zw3d)    ! i-current
222         zw3d(:,:,:) = vn_25h(:,:,:)*vmask(:,:,:) + zmdi*(1.0-vmask(:,:,:))
223         CALL iom_put("vomecrty25h", zw3d  )   ! j-current
224         zw3d(:,:,:) = wn_25h(:,:,:)*wmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
225         CALL iom_put("vomecrtz25h", zw3d )   ! k-current
226         ! Write vertical physics
227         zw3d(:,:,:) = avt_25h(:,:,:)*wmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
228         CALL iom_put("avt25h", zw3d )   ! diffusivity
229         zw3d(:,:,:) = avm_25h(:,:,:)*wmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
230         CALL iom_put("avm25h", zw3d)   ! viscosity
231         IF( ln_zdftke ) THEN
232            zw3d(:,:,:) = en_25h(:,:,:)*wmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
233            CALL iom_put("tke25h", zw3d)   ! tke
234         ENDIF
235         IF( ln_zdfgls ) THEN
236            zw3d(:,:,:) = en_25h(:,:,:)*wmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
237            CALL iom_put("tke25h", zw3d)   ! tke
238            zw3d(:,:,:) = rmxln_25h(:,:,:)*wmask(:,:,:) + zmdi*(1.0-tmask(:,:,:))
239            CALL iom_put( "mxln25h",zw3d)
240         ENDIF
241         !
242         ! After the write reset the values to cnt=1 and sum values equal current value
243         tn_25h  (:,:,:) = tsn (:,:,:,jp_tem)
244         sn_25h  (:,:,:) = tsn (:,:,:,jp_sal)
245         sshn_25h(:,:)   = sshn(:,:)
246         un_25h  (:,:,:) = un  (:,:,:)
247         vn_25h  (:,:,:) = vn  (:,:,:)
248         wn_25h  (:,:,:) = wn  (:,:,:)
249         avt_25h (:,:,:) = avt (:,:,:)
250         avm_25h (:,:,:) = avm (:,:,:)
251         IF( ln_zdftke ) THEN
252            en_25h(:,:,:) = en(:,:,:)
253         ENDIF
254         IF( ln_zdfgls ) THEN
255            en_25h   (:,:,:) = en    (:,:,:)
256            rmxln_25h(:,:,:) = hmxl_n(:,:,:)
257         ENDIF
258         cnt_25h = 1
259         IF(lwp)  WRITE(numout,*) 'dia_wri_tide :   &
260            &    After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average', cnt_25h
261      ENDIF !  cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps * 24) == 0 .AND. kt .NE. nn_it000
262      !
263   END SUBROUTINE dia_25h 
264
265   !!======================================================================
266END MODULE dia25h
Note: See TracBrowser for help on using the repository browser.