source: trunk/NEMOGCM/NEMO/OPA_SRC/C1D/dtauvd.F90 @ 4624

Last change on this file since 4624 was 4624, checked in by acc, 6 years ago

#1305. Fix slow start-up problems on some systems by introducing and using lwm logical to restrict output of merged namelists to the first (or only) processor. lwm is true only on the first processor regardless of ln_ctl. Small changes to all flavours of nemogcm.F90 are also required to write namctl and namcfg after the call to mynode which now opens output.namelist.dyn and writes nammpp.

File size: 13.0 KB
Line 
1MODULE dtauvd
2   !!======================================================================
3   !!                     ***  MODULE  dtauvd  ***
4   !! Ocean data  :  read ocean U & V current data from gridded data
5   !!======================================================================
6   !! History :  3.5   ! 2013-08  (D. Calvert)  Original code
7   !!----------------------------------------------------------------------
8
9   !!----------------------------------------------------------------------
10   !!   dta_uvd_init   : read namelist and allocate data structures
11   !!   dta_uvd        : read and time-interpolate ocean U & V current data
12   !!----------------------------------------------------------------------
13   USE oce             ! ocean dynamics and tracers
14   USE dom_oce         ! ocean space and time domain
15   USE fldread         ! read input fields
16   USE in_out_manager  ! I/O manager
17   USE phycst          ! physical constants
18   USE lib_mpp         ! MPP library
19   USE wrk_nemo        ! Memory allocation
20   USE timing          ! Timing
21
22   IMPLICIT NONE
23   PRIVATE
24
25   PUBLIC   dta_uvd_init   ! called by nemogcm.F90
26   PUBLIC   dta_uvd        ! called by istate.F90 and dyndmp.90
27
28   LOGICAL , PUBLIC ::   ln_uvd_init         ! Flag to initialise with U & V current data
29   LOGICAL , PUBLIC ::   ln_uvd_dyndmp       ! Flag for Newtonian damping toward U & V current data
30
31   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_uvd   ! structure for input U & V current (file information and data)
32
33   !! * Substitutions
34#  include "domzgr_substitute.h90"
35   !!----------------------------------------------------------------------
36   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
37   !! $Id: dtauvd.F90 2392 2010-11-15 21:20:05Z gm $
38   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
39   !!----------------------------------------------------------------------
40CONTAINS
41
42   SUBROUTINE dta_uvd_init( ld_dyndmp )
43      !!----------------------------------------------------------------------
44      !!                   ***  ROUTINE dta_uvd_init  ***
45      !!                   
46      !! ** Purpose :   initialization of U & V current input data
47      !!
48      !! ** Method  : - read namc1d_uvd namelist
49      !!              - allocate U & V current data structure
50      !!              - fld_fill data structure with namelist information
51      !!----------------------------------------------------------------------
52      LOGICAL, INTENT(in), OPTIONAL ::   ld_dyndmp         ! force the initialization when dyndmp is used
53      !
54      INTEGER ::   ierr0, ierr1, ierr2, ierr3              ! temporary integers
55      !
56      CHARACTER(len=100)            ::   cn_dir            ! Root directory for location of files to be used
57      TYPE(FLD_N), DIMENSION(2)     ::   suv_i             ! Combined U & V namelist information
58      TYPE(FLD_N)                   ::   sn_ucur, sn_vcur  ! U & V data namelist information
59      !!
60      NAMELIST/namc1d_uvd/ ln_uvd_init, ln_uvd_dyndmp, cn_dir, sn_ucur, sn_vcur
61      INTEGER  ::   ios
62      !!----------------------------------------------------------------------
63      !
64      IF( nn_timing == 1 )  CALL timing_start('dta_uvd_init')
65      !
66      ierr0 = 0  ;  ierr1 = 0  ;  ierr2 = 0  ;  ierr3 = 0
67
68      REWIND( numnam_ref )              ! Namelist namc1d_uvd in reference namelist :
69      READ  ( numnam_ref, namc1d_uvd, IOSTAT = ios, ERR = 901)
70901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namc1d_uvd in reference namelist', lwp )
71
72      REWIND( numnam_cfg )              ! Namelist namc1d_uvd in configuration namelist : Parameters of the run
73      READ  ( numnam_cfg, namc1d_uvd, IOSTAT = ios, ERR = 902 )
74902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namc1d_uvd in configuration namelist', lwp )
75      IF(lwm) WRITE ( numond, namc1d_uvd )
76
77      !                             ! force the initialization when dyndmp is used
78      IF( PRESENT( ld_dyndmp ) )   ln_uvd_dyndmp = .TRUE.
79     
80      IF(lwp) THEN                  ! control print
81         WRITE(numout,*)
82         WRITE(numout,*) 'dta_uvd_init : U & V current data '
83         WRITE(numout,*) '~~~~~~~~~~~~ '
84         WRITE(numout,*) '   Namelist namc1d_uvd : Set flags'
85         WRITE(numout,*) '      Initialization of ocean U & V current with input data   ln_uvd_init   = ', ln_uvd_init
86         WRITE(numout,*) '      Damping of ocean U & V current toward input data        ln_uvd_dyndmp = ', ln_uvd_dyndmp
87         WRITE(numout,*)
88         IF( .NOT. ln_uvd_init .AND. .NOT. ln_uvd_dyndmp ) THEN
89            WRITE(numout,*)
90            WRITE(numout,*) '   U & V current data not used'
91         ENDIF
92      ENDIF
93      !                             ! no initialization when restarting
94      IF( ln_rstart .AND. ln_uvd_init ) THEN
95         CALL ctl_warn( 'dta_uvd_init: ocean restart and U & V current data initialization, ',   &
96            &           'we keep the restart U & V current values and set ln_uvd_init to FALSE' )
97         ln_uvd_init = .FALSE.
98      ENDIF
99
100      !
101      IF(  ln_uvd_init .OR. ln_uvd_dyndmp  ) THEN
102         !                          !==   allocate the data arrays   ==!
103         ALLOCATE( sf_uvd(2), STAT=ierr0 )
104         IF( ierr0 > 0 ) THEN
105            CALL ctl_stop( 'dta_uvd_init: unable to allocate sf_uvd structure' )             ;   RETURN
106         ENDIF
107         !
108                                 ALLOCATE( sf_uvd(1)%fnow(jpi,jpj,jpk)   , STAT=ierr0 )
109         IF( sn_ucur%ln_tint )   ALLOCATE( sf_uvd(1)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 )
110                                 ALLOCATE( sf_uvd(2)%fnow(jpi,jpj,jpk)   , STAT=ierr2 )
111         IF( sn_vcur%ln_tint )   ALLOCATE( sf_uvd(2)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 )
112         !
113         IF( ierr0 + ierr1 + ierr2 + ierr3 > 0 ) THEN
114            CALL ctl_stop( 'dta_uvd_init : unable to allocate U & V current data arrays' )   ;   RETURN
115         ENDIF
116         !                          !==   fill sf_uvd with sn_ucur, sn_vcur and control print   ==!
117         suv_i(1) = sn_ucur   ;   suv_i(2) = sn_vcur
118         CALL fld_fill( sf_uvd, suv_i, cn_dir, 'dta_uvd', 'U & V current data', 'namc1d_uvd' )
119         !
120      ENDIF
121      !
122      IF( nn_timing == 1 )  CALL timing_stop('dta_uvd_init')
123      !
124   END SUBROUTINE dta_uvd_init
125
126
127   SUBROUTINE dta_uvd( kt, puvd )
128      !!----------------------------------------------------------------------
129      !!                   ***  ROUTINE dta_uvd  ***
130      !!                   
131      !! ** Purpose :   provides U & V current data at time step kt
132      !!
133      !! ** Method  : - call fldread routine
134      !!              - ORCA_R2: make some hand made alterations to the data (EMPTY)
135      !!              - s- or mixed s-zps coordinate: vertical interpolation onto model mesh
136      !!              - zps coordinate: vertical interpolation onto last partial level
137      !!              - ln_uvd_dyndmp=False: deallocate the U & V current data structure,
138      !!                                     as the data is no longer used
139      !!
140      !! ** Action  :   puvd,  U & V current data interpolated onto model mesh at time-step kt
141      !!----------------------------------------------------------------------
142      INTEGER                           , INTENT(in   ) ::   kt     ! ocean time-step
143      REAL(wp), DIMENSION(jpi,jpj,jpk,2), INTENT(  out) ::   puvd   ! U & V current data
144      !
145      INTEGER ::   ji, jj, jk, jl, jkk               ! dummy loop indicies
146      INTEGER ::   ik, il0, il1, ii0, ii1, ij0, ij1  ! local integers
147      REAL(wp)::   zl, zi                            ! local floats
148      REAL(wp), POINTER, DIMENSION(:) ::  zup, zvp   ! 1D workspace
149      !!----------------------------------------------------------------------
150      !
151      IF( nn_timing == 1 )  CALL timing_start('dta_uvd')
152      !
153      CALL fld_read( kt, 1, sf_uvd )      !==   read U & V current data at time step kt   ==!
154      !
155      !
156      !                                   !==   ORCA_R2 configuration and U & V current damping   ==!
157      IF( cp_cfg == "orca" .AND. jp_cfg == 2 .AND. ln_uvd_dyndmp ) THEN    ! some hand made alterations
158         !!! EMPTY- to be added for running in 3D context !!!
159      ENDIF
160      !
161      puvd(:,:,:,1) = sf_uvd(1)%fnow(:,:,:)                 ! NO mask
162      puvd(:,:,:,2) = sf_uvd(2)%fnow(:,:,:) 
163      !
164      IF( ln_sco ) THEN                   !==   s- or mixed s-zps-coordinate   ==!
165         !
166         CALL wrk_alloc( jpk, zup, zvp )
167         !
168         IF( kt == nit000 .AND. lwp )THEN
169            WRITE(numout,*)
170            WRITE(numout,*) 'dta_uvd: interpolate U & V current data onto the s- or mixed s-z-coordinate mesh'
171         ENDIF
172         !
173         DO jj = 1, jpj                   ! vertical interpolation of U & V current:
174            DO ji = 1, jpi                ! determines the interpolated U & V current profiles at each (i,j) point
175               DO jk = 1, jpk
176                  zl = fsdept(ji,jj,jk)
177                  IF    ( zl < gdept_1d(1  ) ) THEN          ! extrapolate above the first level of data
178                     zup(jk) =  puvd(ji,jj,1    ,1)
179                     zvp(jk) =  puvd(ji,jj,1    ,2)
180                  ELSEIF( zl > gdept_1d(jpk) ) THEN          ! extrapolate below the last level of data
181                     zup(jk) =  puvd(ji,jj,jpkm1,1)
182                     zvp(jk) =  puvd(ji,jj,jpkm1,2)
183                  ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1
184                     DO jkk = 1, jpkm1                      ! when  gdept(jkk) < zl < gdept(jkk+1)
185                        IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN
186                           zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk))
187                           zup(jk) = puvd(ji,jj,jkk,1) + ( puvd(ji,jj,jkk+1,1 ) - puvd(ji,jj,jkk,1) ) * zi 
188                           zvp(jk) = puvd(ji,jj,jkk,2) + ( puvd(ji,jj,jkk+1,2 ) - puvd(ji,jj,jkk,2) ) * zi
189                        ENDIF
190                     END DO
191                  ENDIF
192               END DO
193               DO jk = 1, jpkm1           ! apply mask
194                  puvd(ji,jj,jk,1) = zup(jk) * umask(ji,jj,jk)
195                  puvd(ji,jj,jk,2) = zvp(jk) * vmask(ji,jj,jk)
196               END DO
197               puvd(ji,jj,jpk,1) = 0._wp
198               puvd(ji,jj,jpk,2) = 0._wp
199            END DO
200         END DO
201         !
202         CALL wrk_dealloc( jpk, zup, zvp )
203         !
204      ELSE                                !==   z- or zps- coordinate   ==!
205         !                             
206         puvd(:,:,:,1) = puvd(:,:,:,1) * umask(:,:,:)       ! apply mask
207         puvd(:,:,:,2) = puvd(:,:,:,2) * vmask(:,:,:)
208         !
209         IF( ln_zps ) THEN                ! zps-coordinate (partial steps) interpolation at the last ocean level
210            DO jj = 1, jpj
211               DO ji = 1, jpi
212                  ik = mbkt(ji,jj) 
213                  IF( ik > 1 ) THEN
214                     zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) )
215                     puvd(ji,jj,ik,1) = (1.-zl) * puvd(ji,jj,ik,1) + zl * puvd(ji,jj,ik-1,1)
216                     puvd(ji,jj,ik,2) = (1.-zl) * puvd(ji,jj,ik,2) + zl * puvd(ji,jj,ik-1,2)
217                  ENDIF
218               END DO
219            END DO
220         ENDIF
221         !
222      ENDIF
223      !
224      IF( lwp .AND. kt == nit000 ) THEN   ! control print
225         WRITE(numout,*) ' U current '
226         WRITE(numout,*)
227         WRITE(numout,*)'  level = 1'
228         CALL prihre( puvd(:,:,1    ,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
229         WRITE(numout,*)'  level = ', jpk/2
230         CALL prihre( puvd(:,:,jpk/2,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
231         WRITE(numout,*)'  level = ', jpkm1
232         CALL prihre( puvd(:,:,jpkm1,1), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
233         WRITE(numout,*)
234         WRITE(numout,*) ' V current '
235         WRITE(numout,*)
236         WRITE(numout,*)'  level = 1'
237         CALL prihre( puvd(:,:,1    ,2), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
238         WRITE(numout,*)'  level = ', jpk/2
239         CALL prihre( puvd(:,:,jpk/2,2), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
240         WRITE(numout,*)'  level = ', jpkm1
241         CALL prihre( puvd(:,:,jpkm1,2), jpi, jpj, 1, jpi, 20, 1, jpj, 20, 1., numout )
242         WRITE(numout,*)
243      ENDIF
244      !
245      IF( .NOT. ln_uvd_dyndmp    ) THEN   !==   deallocate U & V current structure   ==!
246         !                                !==   (data used only for initialization)  ==!
247         IF(lwp) WRITE(numout,*) 'dta_uvd: deallocate U & V current arrays as they are only used to initialize the run'
248                                   DEALLOCATE( sf_uvd(1)%fnow )     ! U current arrays in the structure
249         IF( sf_uvd(1)%ln_tint )   DEALLOCATE( sf_uvd(1)%fdta )
250                                   DEALLOCATE( sf_uvd(2)%fnow )     ! V current arrays in the structure
251         IF( sf_uvd(2)%ln_tint )   DEALLOCATE( sf_uvd(2)%fdta )
252                                   DEALLOCATE( sf_uvd         )     ! the structure itself
253      ENDIF
254      !
255      IF( nn_timing == 1 )  CALL timing_stop('dta_uvd')
256      !
257   END SUBROUTINE dta_uvd
258
259   !!======================================================================
260END MODULE dtauvd
Note: See TracBrowser for help on using the repository browser.