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.
sbcwave.F90 in NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/SBC – NEMO

source: NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/SBC/sbcwave.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: 25.3 KB
Line 
1MODULE sbcwave
2   !!======================================================================
3   !!                       ***  MODULE  sbcwave  ***
4   !! Wave module
5   !!======================================================================
6   !! History :  3.3  !  2011-09  (M. Adani)  Original code: Drag Coefficient
7   !!         :  3.4  !  2012-10  (M. Adani)  Stokes Drift
8   !!            3.6  !  2014-09  (E. Clementi,P. Oddo) New Stokes Drift Computation
9   !!             -   !  2016-12  (G. Madec, E. Clementi) update Stoke drift computation
10   !!                                                    + add sbc_wave_ini routine
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   sbc_stokes    : calculate 3D Stokes-drift velocities
15   !!   sbc_wave      : wave data from wave model in netcdf files
16   !!   sbc_wave_init : initialisation fo surface waves
17   !!----------------------------------------------------------------------
18   USE phycst         ! physical constants
19   USE oce            ! ocean variables
20   USE sbc_oce        ! Surface boundary condition: ocean fields
21   USE zdf_oce,  ONLY : ln_zdfswm
22   USE bdy_oce        ! open boundary condition variables
23   USE domvvl         ! domain: variable volume layers
24   !
25   USE iom            ! I/O manager library
26   USE in_out_manager ! I/O manager
27   USE lib_mpp        ! distribued memory computing library
28   USE fldread        ! read input fields
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   sbc_stokes      ! routine called in sbccpl
34   PUBLIC   sbc_wstress     ! routine called in sbcmod
35   PUBLIC   sbc_wave        ! routine called in sbcmod
36   PUBLIC   sbc_wave_init   ! routine called in sbcmod
37   
38   ! Variables checking if the wave parameters are coupled (if not, they are read from file)
39   LOGICAL, PUBLIC ::   cpl_hsig   = .FALSE.
40   LOGICAL, PUBLIC ::   cpl_phioc  = .FALSE.
41   LOGICAL, PUBLIC ::   cpl_sdrftx = .FALSE.
42   LOGICAL, PUBLIC ::   cpl_sdrfty = .FALSE.
43   LOGICAL, PUBLIC ::   cpl_wper   = .FALSE.
44   LOGICAL, PUBLIC ::   cpl_wfreq  = .FALSE.
45   LOGICAL, PUBLIC ::   cpl_wnum   = .FALSE.
46   LOGICAL, PUBLIC ::   cpl_tauwoc = .FALSE.
47   LOGICAL, PUBLIC ::   cpl_tauw   = .FALSE.
48   LOGICAL, PUBLIC ::   cpl_wdrag  = .FALSE.
49
50   INTEGER ::   jpfld    ! number of files to read for stokes drift
51   INTEGER ::   jp_usd   ! index of stokes drift  (i-component) (m/s)    at T-point
52   INTEGER ::   jp_vsd   ! index of stokes drift  (j-component) (m/s)    at T-point
53   INTEGER ::   jp_hsw   ! index of significant wave hight      (m)      at T-point
54   INTEGER ::   jp_wmp   ! index of mean wave period            (s)      at T-point
55   INTEGER ::   jp_wfr   ! index of wave peak frequency         (1/s)    at T-point
56
57   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_cd      ! structure of input fields (file informations, fields read) Drag Coefficient
58   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sd      ! structure of input fields (file informations, fields read) Stokes Drift
59   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_wn      ! structure of input fields (file informations, fields read) wave number for Qiao
60   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tauwoc  ! structure of input fields (file informations, fields read) normalized wave stress into the ocean
61   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tauw    ! structure of input fields (file informations, fields read) ocean stress components from wave model
62
63   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   cdn_wave            !:
64   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   hsw, wmp, wnum      !:
65   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   wfreq               !:
66   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tauoc_wave          !: 
67   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tauw_x, tauw_y      !: 
68   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tsd2d               !:
69   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   div_sd              !: barotropic stokes drift divergence
70   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   ut0sd, vt0sd        !: surface Stokes drift velocities at t-point
71   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) ::   usd  , vsd  , wsd   !: Stokes drift velocities at u-, v- & w-points, resp.
72
73   !! * Substitutions
74#  include "vectopt_loop_substitute.h90"
75   !!----------------------------------------------------------------------
76   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
77   !! $Id$
78   !! Software governed by the CeCILL license (see ./LICENSE)
79   !!----------------------------------------------------------------------
80CONTAINS
81
82   SUBROUTINE sbc_stokes( )
83      !!---------------------------------------------------------------------
84      !!                     ***  ROUTINE sbc_stokes  ***
85      !!
86      !! ** Purpose :   compute the 3d Stokes Drift according to Breivik et al.,
87      !!                2014 (DOI: 10.1175/JPO-D-14-0020.1)
88      !!
89      !! ** Method  : - Calculate Stokes transport speed
90      !!              - Calculate horizontal divergence
91      !!              - Integrate the horizontal divergenze from the bottom
92      !! ** action 
93      !!---------------------------------------------------------------------
94      INTEGER  ::   jj, ji, jk   ! dummy loop argument
95      INTEGER  ::   ik           ! local integer
96      REAL(wp) ::  ztransp, zfac, zsp0
97      REAL(wp) ::  zdepth, zsqrt_depth,  zexp_depth, z_two_thirds, zsqrtpi !sqrt of pi
98      REAL(wp) ::  zbot_u, zbot_v, zkb_u, zkb_v, zke3_u, zke3_v, zda_u, zda_v
99      REAL(wp) ::  zstokes_psi_u_bot, zstokes_psi_v_bot
100      REAL(wp) ::  zdep_u, zdep_v, zkh_u, zkh_v
101      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE ::   zk_t, zk_u, zk_v, zu0_sd, zv0_sd     ! 2D workspace
102      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE ::   zstokes_psi_u_top, zstokes_psi_v_top ! 2D workspace
103      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ze3divh                              ! 3D workspace
104      !!---------------------------------------------------------------------
105      !
106      ALLOCATE( ze3divh(jpi,jpj,jpk) )
107      ALLOCATE( zk_t(jpi,jpj), zk_u(jpi,jpj), zk_v(jpi,jpj), zu0_sd(jpi,jpj), zv0_sd(jpi,jpj) )
108      !
109      ! select parameterization for the calculation of vertical Stokes drift
110      ! exp. wave number at t-point
111      IF( ll_st_bv_li ) THEN   ! (Eq. (19) in Breivik et al. (2014) )
112         zfac = 2.0_wp * rpi / 16.0_wp
113         DO jj = 1, jpj
114            DO ji = 1, jpi
115               ! Stokes drift velocity estimated from Hs and Tmean
116               ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp )
117               ! Stokes surface speed
118               tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj))
119               ! Wavenumber scale
120               zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / MAX( ABS( 5.97_wp*ztransp ), 0.0000001_wp )
121            END DO
122         END DO
123         DO jj = 1, jpjm1              ! exp. wave number & Stokes drift velocity at u- & v-points
124            DO ji = 1, jpim1
125               zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) )
126               zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) )
127               !
128               zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )
129               zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) )
130            END DO
131         END DO
132      ELSE IF( ll_st_peakfr ) THEN    ! peak wave number calculated from the peak frequency received by the wave model
133         DO jj = 1, jpj
134            DO ji = 1, jpi
135               zk_t(ji,jj) = ( 2.0_wp * rpi * wfreq(ji,jj) ) * ( 2.0_wp * rpi * wfreq(ji,jj) ) / grav
136            END DO
137         END DO
138         DO jj = 1, jpjm1
139            DO ji = 1, jpim1
140               zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) )
141               zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) )
142               !
143               zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )
144               zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) )
145            END DO
146         END DO
147      ENDIF
148      !
149      !                       !==  horizontal Stokes Drift 3D velocity  ==!
150      IF( ll_st_bv2014 ) THEN
151         DO jk = 1, jpkm1
152            DO jj = 2, jpjm1
153               DO ji = 2, jpim1
154                  zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) )
155                  zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) )
156                  !                         
157                  zkh_u = zk_u(ji,jj) * zdep_u     ! k * depth
158                  zkh_v = zk_v(ji,jj) * zdep_v
159                  !                                ! Depth attenuation
160                  zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u )
161                  zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v )
162                  !
163                  usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)
164                  vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)
165               END DO
166            END DO
167         END DO
168      ELSE IF( ll_st_li2017 .OR. ll_st_peakfr ) THEN
169         ALLOCATE( zstokes_psi_u_top(jpi,jpj), zstokes_psi_v_top(jpi,jpj) )
170         DO jj = 1, jpjm1              ! exp. wave number & Stokes drift velocity at u- & v-points
171            DO ji = 1, jpim1
172               zstokes_psi_u_top(ji,jj) = 0._wp
173               zstokes_psi_v_top(ji,jj) = 0._wp
174            END DO
175         END DO
176         zsqrtpi = SQRT(rpi)
177         z_two_thirds = 2.0_wp / 3.0_wp
178         DO jk = 1, jpkm1
179            DO jj = 2, jpjm1
180               DO ji = 2, jpim1
181                  zbot_u = ( gdepw_n(ji,jj,jk+1) + gdepw_n(ji+1,jj,jk+1) )  ! 2 * bottom depth
182                  zbot_v = ( gdepw_n(ji,jj,jk+1) + gdepw_n(ji,jj+1,jk+1) )  ! 2 * bottom depth
183                  zkb_u  = zk_u(ji,jj) * zbot_u                             ! 2 * k * bottom depth
184                  zkb_v  = zk_v(ji,jj) * zbot_v                             ! 2 * k * bottom depth
185                  !
186                  zke3_u = MAX(1.e-8_wp, 2.0_wp * zk_u(ji,jj) * e3u_n(ji,jj,jk))     ! 2k * thickness
187                  zke3_v = MAX(1.e-8_wp, 2.0_wp * zk_v(ji,jj) * e3v_n(ji,jj,jk))     ! 2k * thickness
188
189                  ! Depth attenuation .... do u component first..
190                  zdepth      = zkb_u
191                  zsqrt_depth = SQRT(zdepth)
192                  zexp_depth  = EXP(-zdepth)
193                  zstokes_psi_u_bot = 1.0_wp - zexp_depth  &
194                       &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) &
195                       &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth )
196                  zda_u                    = ( zstokes_psi_u_bot - zstokes_psi_u_top(ji,jj) ) / zke3_u
197                  zstokes_psi_u_top(ji,jj) =   zstokes_psi_u_bot
198
199                  !         ... and then v component
200                  zdepth      =zkb_v
201                  zsqrt_depth = SQRT(zdepth)
202                  zexp_depth  = EXP(-zdepth)
203                  zstokes_psi_v_bot = 1.0_wp - zexp_depth  &
204                       &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) &
205                       &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth )
206                  zda_v                    = ( zstokes_psi_v_bot - zstokes_psi_v_top(ji,jj) ) / zke3_v
207                  zstokes_psi_v_top(ji,jj) =   zstokes_psi_v_bot
208                  !
209                  usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)
210                  vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)
211               END DO
212            END DO
213         END DO
214         DEALLOCATE( zstokes_psi_u_top, zstokes_psi_v_top )
215      ENDIF
216
217      CALL lbc_lnk_multi( 'sbcwave', usd, 'U', -1., vsd, 'V', -1. )
218
219      !
220      !                       !==  vertical Stokes Drift 3D velocity  ==!
221      !
222      DO jk = 1, jpkm1               ! Horizontal e3*divergence
223         DO jj = 2, jpj
224            DO ji = fs_2, jpi
225               ze3divh(ji,jj,jk) = (  e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) * usd(ji  ,jj,jk)    &
226                  &                 - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd(ji-1,jj,jk)    &
227                  &                 + e1v(ji,jj  ) * e3v_n(ji,jj  ,jk) * vsd(ji,jj  ,jk)    &
228                  &                 - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd(ji,jj-1,jk)  ) * r1_e1e2t(ji,jj)
229            END DO
230         END DO
231      END DO
232      !
233#if defined key_agrif
234      IF( .NOT. Agrif_Root() ) THEN
235         IF( nbondi == -1 .OR. nbondi == 2 )   ze3divh( 2:nbghostcells+1,:        ,:) = 0._wp      ! west
236         IF( nbondi ==  1 .OR. nbondi == 2 )   ze3divh( nlci-nbghostcells:nlci-1,:,:) = 0._wp      ! east
237         IF( nbondj == -1 .OR. nbondj == 2 )   ze3divh( :,2:nbghostcells+1        ,:) = 0._wp      ! south
238         IF( nbondj ==  1 .OR. nbondj == 2 )   ze3divh( :,nlcj-nbghostcells:nlcj-1,:) = 0._wp      ! north
239      ENDIF
240#endif
241      !
242      CALL lbc_lnk( 'sbcwave', ze3divh, 'T', 1. )
243      !
244      IF( ln_linssh ) THEN   ;   ik = 1   ! none zero velocity through the sea surface
245      ELSE                   ;   ik = 2   ! w=0 at the surface (set one for all in sbc_wave_init)
246      ENDIF
247      DO jk = jpkm1, ik, -1          ! integrate from the bottom the hor. divergence (NB: at k=jpk w is always zero)
248         wsd(:,:,jk) = wsd(:,:,jk+1) - ze3divh(:,:,jk)
249      END DO
250      !
251      IF( ln_bdy ) THEN
252         DO jk = 1, jpkm1
253            wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:)
254         END DO
255      ENDIF
256      !                       !==  Horizontal divergence of barotropic Stokes transport  ==!
257      div_sd(:,:) = 0._wp
258      DO jk = 1, jpkm1                                 !
259        div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk)
260      END DO
261      !
262      CALL iom_put( "ustokes",  usd  )
263      CALL iom_put( "vstokes",  vsd  )
264      CALL iom_put( "wstokes",  wsd  )
265      !
266      DEALLOCATE( ze3divh )
267      DEALLOCATE( zk_t, zk_u, zk_v, zu0_sd, zv0_sd )
268      !
269   END SUBROUTINE sbc_stokes
270
271
272   SUBROUTINE sbc_wstress( )
273      !!---------------------------------------------------------------------
274      !!                     ***  ROUTINE sbc_wstress  ***
275      !!
276      !! ** Purpose :   Updates the ocean momentum modified by waves
277      !!
278      !! ** Method  : - Calculate u,v components of stress depending on stress
279      !!                model
280      !!              - Calculate the stress module
281      !!              - The wind module is not modified by waves
282      !! ** action 
283      !!---------------------------------------------------------------------
284      INTEGER  ::   jj, ji   ! dummy loop argument
285      !
286      IF( ln_tauwoc ) THEN
287         utau(:,:) = utau(:,:)*tauoc_wave(:,:)
288         vtau(:,:) = vtau(:,:)*tauoc_wave(:,:)
289         taum(:,:) = taum(:,:)*tauoc_wave(:,:)
290      ENDIF
291      !
292      IF( ln_tauw ) THEN
293         DO jj = 1, jpjm1
294            DO ji = 1, jpim1
295               ! Stress components at u- & v-points
296               utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) )
297               vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) )
298               !
299               ! Stress module at t points
300               taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) )
301            END DO
302         END DO
303         CALL lbc_lnk_multi( 'sbcwave', utau(:,:), 'U', -1. , vtau(:,:), 'V', -1. , taum(:,:) , 'T', -1. )
304      ENDIF
305      !
306   END SUBROUTINE sbc_wstress
307
308
309   SUBROUTINE sbc_wave( kt )
310      !!---------------------------------------------------------------------
311      !!                     ***  ROUTINE sbc_wave  ***
312      !!
313      !! ** Purpose :   read wave parameters from wave model  in netcdf files.
314      !!
315      !! ** Method  : - Read namelist namsbc_wave
316      !!              - Read Cd_n10 fields in netcdf files
317      !!              - Read stokes drift 2d in netcdf files
318      !!              - Read wave number in netcdf files
319      !!              - Compute 3d stokes drift using Breivik et al.,2014
320      !!                formulation
321      !! ** action 
322      !!---------------------------------------------------------------------
323      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
324      !!---------------------------------------------------------------------
325      !
326      IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN     !==  Neutral drag coefficient  ==!
327         CALL fld_read( kt, nn_fsbc, sf_cd )             ! read from external forcing
328         cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) * tmask(:,:,1)
329      ENDIF
330
331      IF( ln_tauwoc .AND. .NOT. cpl_tauwoc ) THEN  !==  Wave induced stress  ==!
332         CALL fld_read( kt, nn_fsbc, sf_tauwoc )         ! read wave norm stress from external forcing
333         tauoc_wave(:,:) = sf_tauwoc(1)%fnow(:,:,1) * tmask(:,:,1)
334      ENDIF
335
336      IF( ln_tauw .AND. .NOT. cpl_tauw ) THEN      !==  Wave induced stress  ==!
337         CALL fld_read( kt, nn_fsbc, sf_tauw )           ! read ocean stress components from external forcing (T grid)
338         tauw_x(:,:) = sf_tauw(1)%fnow(:,:,1) * tmask(:,:,1)
339         tauw_y(:,:) = sf_tauw(2)%fnow(:,:,1) * tmask(:,:,1)
340      ENDIF
341
342      IF( ln_sdw )  THEN                           !==  Computation of the 3d Stokes Drift  ==!
343         !
344         IF( jpfld > 0 ) THEN                            ! Read from file only if the field is not coupled
345            CALL fld_read( kt, nn_fsbc, sf_sd )          ! read wave parameters from external forcing
346            IF( jp_hsw > 0 )   hsw  (:,:) = sf_sd(jp_hsw)%fnow(:,:,1) * tmask(:,:,1)  ! significant wave height
347            IF( jp_wmp > 0 )   wmp  (:,:) = sf_sd(jp_wmp)%fnow(:,:,1) * tmask(:,:,1)  ! wave mean period
348            IF( jp_wfr > 0 )   wfreq(:,:) = sf_sd(jp_wfr)%fnow(:,:,1) * tmask(:,:,1)  ! Peak wave frequency
349            IF( jp_usd > 0 )   ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1) * tmask(:,:,1)  ! 2D zonal Stokes Drift at T point
350            IF( jp_vsd > 0 )   vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) * tmask(:,:,1)  ! 2D meridional Stokes Drift at T point
351         ENDIF
352         !
353         ! Read also wave number if needed, so that it is available in coupling routines
354         IF( ln_zdfswm .AND. .NOT.cpl_wnum ) THEN
355            CALL fld_read( kt, nn_fsbc, sf_wn )          ! read wave parameters from external forcing
356            wnum(:,:) = sf_wn(1)%fnow(:,:,1) * tmask(:,:,1)
357         ENDIF
358           
359         ! Calculate only if required fields have been read
360         ! In coupled wave model-NEMO case the call is done after coupling
361         !
362         IF( ( ll_st_bv_li   .AND. jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0 ) .OR. &
363           & ( ll_st_peakfr  .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0                ) ) CALL sbc_stokes()
364         !
365      ENDIF
366      !
367   END SUBROUTINE sbc_wave
368
369
370   SUBROUTINE sbc_wave_init
371      !!---------------------------------------------------------------------
372      !!                     ***  ROUTINE sbc_wave_init  ***
373      !!
374      !! ** Purpose :   read wave parameters from wave model  in netcdf files.
375      !!
376      !! ** Method  : - Read namelist namsbc_wave
377      !!              - Read Cd_n10 fields in netcdf files
378      !!              - Read stokes drift 2d in netcdf files
379      !!              - Read wave number in netcdf files
380      !!              - Compute 3d stokes drift using Breivik et al.,2014
381      !!                formulation
382      !! ** action 
383      !!---------------------------------------------------------------------
384      INTEGER ::   ierror, ios   ! local integer
385      INTEGER ::   ifpr
386      !!
387      CHARACTER(len=100)     ::  cn_dir                            ! Root directory for location of drag coefficient files
388      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i, slf_j     ! array of namelist informations on the fields to read
389      TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd,  &
390                             &   sn_hsw, sn_wmp, sn_wfr, sn_wnum, &
391                             &   sn_tauwoc, sn_tauwx, sn_tauwy     ! informations about the fields to be read
392      !
393      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, &
394                             sn_wnum, sn_tauwoc, sn_tauwx, sn_tauwy
395      !!---------------------------------------------------------------------
396      !
397      READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901)
398901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc_wave in reference namelist' )
399         
400      READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 )
401902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist' )
402      IF(lwm) WRITE ( numond, namsbc_wave )
403      !
404      IF( ln_cdgw ) THEN
405         IF( .NOT. cpl_wdrag ) THEN
406            ALLOCATE( sf_cd(1), STAT=ierror )               !* allocate and fill sf_wave with sn_cdg
407            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' )
408            !
409                                   ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)   )
410            IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) )
411            CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' )
412         ENDIF
413         ALLOCATE( cdn_wave(jpi,jpj) )
414      ENDIF
415
416      IF( ln_tauwoc ) THEN
417         IF( .NOT. cpl_tauwoc ) THEN
418            ALLOCATE( sf_tauwoc(1), STAT=ierror )           !* allocate and fill sf_wave with sn_tauwoc
419            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' )
420            !
421                                     ALLOCATE( sf_tauwoc(1)%fnow(jpi,jpj,1)   )
422            IF( sn_tauwoc%ln_tint )  ALLOCATE( sf_tauwoc(1)%fdta(jpi,jpj,1,2) )
423            CALL fld_fill( sf_tauwoc, (/ sn_tauwoc /), cn_dir, 'sbc_wave_init', 'Wave module', 'namsbc_wave' )
424         ENDIF
425         ALLOCATE( tauoc_wave(jpi,jpj) )
426      ENDIF
427
428      IF( ln_tauw ) THEN
429         IF( .NOT. cpl_tauw ) THEN
430            ALLOCATE( sf_tauw(2), STAT=ierror )           !* allocate and fill sf_wave with sn_tauwx/y
431            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauw structure' )
432            !
433            ALLOCATE( slf_j(2) )
434            slf_j(1) = sn_tauwx
435            slf_j(2) = sn_tauwy
436                                    ALLOCATE( sf_tauw(1)%fnow(jpi,jpj,1)   )
437                                    ALLOCATE( sf_tauw(2)%fnow(jpi,jpj,1)   )
438            IF( slf_j(1)%ln_tint )  ALLOCATE( sf_tauw(1)%fdta(jpi,jpj,1,2) )
439            IF( slf_j(2)%ln_tint )  ALLOCATE( sf_tauw(2)%fdta(jpi,jpj,1,2) )
440            CALL fld_fill( sf_tauw, (/ slf_j /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
441         ENDIF
442         ALLOCATE( tauw_x(jpi,jpj) )
443         ALLOCATE( tauw_y(jpi,jpj) )
444      ENDIF
445
446      IF( ln_sdw ) THEN   ! Find out how many fields have to be read from file if not coupled
447         jpfld=0
448         jp_usd=0   ;   jp_vsd=0   ;   jp_hsw=0   ;   jp_wmp=0   ;   jp_wfr=0
449         IF( .NOT. cpl_sdrftx ) THEN
450            jpfld  = jpfld + 1
451            jp_usd = jpfld
452         ENDIF
453         IF( .NOT. cpl_sdrfty ) THEN
454            jpfld  = jpfld + 1
455            jp_vsd = jpfld
456         ENDIF
457         IF( .NOT. cpl_hsig  .AND. ll_st_bv_li  ) THEN
458            jpfld  = jpfld + 1
459            jp_hsw = jpfld
460         ENDIF
461         IF( .NOT. cpl_wper  .AND. ll_st_bv_li  ) THEN
462            jpfld  = jpfld + 1
463            jp_wmp = jpfld
464         ENDIF
465         IF( .NOT. cpl_wfreq .AND. ll_st_peakfr ) THEN
466            jpfld  = jpfld + 1
467            jp_wfr = jpfld
468         ENDIF
469
470         ! Read from file only the non-coupled fields
471         IF( jpfld > 0 ) THEN
472            ALLOCATE( slf_i(jpfld) )
473            IF( jp_usd > 0 )   slf_i(jp_usd) = sn_usd
474            IF( jp_vsd > 0 )   slf_i(jp_vsd) = sn_vsd
475            IF( jp_hsw > 0 )   slf_i(jp_hsw) = sn_hsw
476            IF( jp_wmp > 0 )   slf_i(jp_wmp) = sn_wmp
477            IF( jp_wfr > 0 )   slf_i(jp_wfr) = sn_wfr
478
479            ALLOCATE( sf_sd(jpfld), STAT=ierror )   !* allocate and fill sf_sd with stokes drift
480            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' )
481            !
482            DO ifpr= 1, jpfld
483               ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) )
484               IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) )
485            END DO
486            !
487            CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' )
488         ENDIF
489         ALLOCATE( usd  (jpi,jpj,jpk), vsd  (jpi,jpj,jpk), wsd(jpi,jpj,jpk) )
490         ALLOCATE( hsw  (jpi,jpj)    , wmp  (jpi,jpj)     )
491         ALLOCATE( wfreq(jpi,jpj) )
492         ALLOCATE( ut0sd(jpi,jpj)    , vt0sd(jpi,jpj)     )
493         ALLOCATE( div_sd(jpi,jpj) )
494         ALLOCATE( tsd2d (jpi,jpj) )
495
496         ut0sd(:,:) = 0._wp
497         vt0sd(:,:) = 0._wp
498         hsw(:,:) = 0._wp
499         wmp(:,:) = 0._wp
500
501         usd(:,:,:) = 0._wp
502         vsd(:,:,:) = 0._wp
503         wsd(:,:,:) = 0._wp
504         ! Wave number needed only if ln_zdfswm=T
505         IF( .NOT. cpl_wnum ) THEN
506            ALLOCATE( sf_wn(1), STAT=ierror )           !* allocate and fill sf_wave with sn_wnum
507            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wave structure' )
508                                   ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1)   )
509            IF( sn_wnum%ln_tint )  ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) )
510            CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' )
511         ENDIF
512         ALLOCATE( wnum(jpi,jpj) )
513      ENDIF
514      !
515   END SUBROUTINE sbc_wave_init
516
517   !!======================================================================
518END MODULE sbcwave
Note: See TracBrowser for help on using the repository browser.