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/trunk/src/OCE/SBC – NEMO

source: NEMO/trunk/src/OCE/SBC/sbcwave.F90 @ 12489

Last change on this file since 12489 was 12377, checked in by acc, 4 years ago

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

​svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
​svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.
  • Property svn:keywords set to Id
File size: 24.4 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 "do_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( Kmm )
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, INTENT(in) :: Kmm ! ocean time level index
95      INTEGER  ::   jj, ji, jk   ! dummy loop argument
96      INTEGER  ::   ik           ! local integer
97      REAL(wp) ::  ztransp, zfac, zsp0
98      REAL(wp) ::  zdepth, zsqrt_depth,  zexp_depth, z_two_thirds, zsqrtpi !sqrt of pi
99      REAL(wp) ::  zbot_u, zbot_v, zkb_u, zkb_v, zke3_u, zke3_v, zda_u, zda_v
100      REAL(wp) ::  zstokes_psi_u_bot, zstokes_psi_v_bot
101      REAL(wp) ::  zdep_u, zdep_v, zkh_u, zkh_v
102      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE ::   zk_t, zk_u, zk_v, zu0_sd, zv0_sd     ! 2D workspace
103      REAL(wp), DIMENSION(:,:)  , ALLOCATABLE ::   zstokes_psi_u_top, zstokes_psi_v_top ! 2D workspace
104      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ze3divh                              ! 3D workspace
105      !!---------------------------------------------------------------------
106      !
107      ALLOCATE( ze3divh(jpi,jpj,jpk) )
108      ALLOCATE( zk_t(jpi,jpj), zk_u(jpi,jpj), zk_v(jpi,jpj), zu0_sd(jpi,jpj), zv0_sd(jpi,jpj) )
109      !
110      ! select parameterization for the calculation of vertical Stokes drift
111      ! exp. wave number at t-point
112      IF( ll_st_bv_li ) THEN   ! (Eq. (19) in Breivik et al. (2014) )
113         zfac = 2.0_wp * rpi / 16.0_wp
114         DO_2D_11_11
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_2D
122         DO_2D_10_10
123            zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) )
124            zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) )
125            !
126            zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )
127            zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) )
128         END_2D
129      ELSE IF( ll_st_peakfr ) THEN    ! peak wave number calculated from the peak frequency received by the wave model
130         DO_2D_11_11
131            zk_t(ji,jj) = ( 2.0_wp * rpi * wfreq(ji,jj) ) * ( 2.0_wp * rpi * wfreq(ji,jj) ) / grav
132         END_2D
133         DO_2D_10_10
134            zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) )
135            zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) )
136            !
137            zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) )
138            zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) )
139         END_2D
140      ENDIF
141      !
142      !                       !==  horizontal Stokes Drift 3D velocity  ==!
143      IF( ll_st_bv2014 ) THEN
144         DO_3D_00_00( 1, jpkm1 )
145            zdep_u = 0.5_wp * ( gdept(ji,jj,jk,Kmm) + gdept(ji+1,jj,jk,Kmm) )
146            zdep_v = 0.5_wp * ( gdept(ji,jj,jk,Kmm) + gdept(ji,jj+1,jk,Kmm) )
147            !                         
148            zkh_u = zk_u(ji,jj) * zdep_u     ! k * depth
149            zkh_v = zk_v(ji,jj) * zdep_v
150            !                                ! Depth attenuation
151            zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u )
152            zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v )
153            !
154            usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)
155            vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)
156         END_3D
157      ELSE IF( ll_st_li2017 .OR. ll_st_peakfr ) THEN
158         ALLOCATE( zstokes_psi_u_top(jpi,jpj), zstokes_psi_v_top(jpi,jpj) )
159         DO_2D_10_10
160            zstokes_psi_u_top(ji,jj) = 0._wp
161            zstokes_psi_v_top(ji,jj) = 0._wp
162         END_2D
163         zsqrtpi = SQRT(rpi)
164         z_two_thirds = 2.0_wp / 3.0_wp
165         DO_3D_00_00( 1, jpkm1 )
166            zbot_u = ( gdepw(ji,jj,jk+1,Kmm) + gdepw(ji+1,jj,jk+1,Kmm) )  ! 2 * bottom depth
167            zbot_v = ( gdepw(ji,jj,jk+1,Kmm) + gdepw(ji,jj+1,jk+1,Kmm) )  ! 2 * bottom depth
168            zkb_u  = zk_u(ji,jj) * zbot_u                             ! 2 * k * bottom depth
169            zkb_v  = zk_v(ji,jj) * zbot_v                             ! 2 * k * bottom depth
170            !
171            zke3_u = MAX(1.e-8_wp, 2.0_wp * zk_u(ji,jj) * e3u(ji,jj,jk,Kmm))     ! 2k * thickness
172            zke3_v = MAX(1.e-8_wp, 2.0_wp * zk_v(ji,jj) * e3v(ji,jj,jk,Kmm))     ! 2k * thickness
173
174            ! Depth attenuation .... do u component first..
175            zdepth      = zkb_u
176            zsqrt_depth = SQRT(zdepth)
177            zexp_depth  = EXP(-zdepth)
178            zstokes_psi_u_bot = 1.0_wp - zexp_depth  &
179                 &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) &
180                 &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth )
181            zda_u                    = ( zstokes_psi_u_bot - zstokes_psi_u_top(ji,jj) ) / zke3_u
182            zstokes_psi_u_top(ji,jj) =   zstokes_psi_u_bot
183
184            !         ... and then v component
185            zdepth      =zkb_v
186            zsqrt_depth = SQRT(zdepth)
187            zexp_depth  = EXP(-zdepth)
188            zstokes_psi_v_bot = 1.0_wp - zexp_depth  &
189                 &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) &
190                 &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth )
191            zda_v                    = ( zstokes_psi_v_bot - zstokes_psi_v_top(ji,jj) ) / zke3_v
192            zstokes_psi_v_top(ji,jj) =   zstokes_psi_v_bot
193            !
194            usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk)
195            vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk)
196         END_3D
197         DEALLOCATE( zstokes_psi_u_top, zstokes_psi_v_top )
198      ENDIF
199
200      CALL lbc_lnk_multi( 'sbcwave', usd, 'U', -1., vsd, 'V', -1. )
201
202      !
203      !                       !==  vertical Stokes Drift 3D velocity  ==!
204      !
205      DO_3D_01_01( 1, jpkm1 )
206         ze3divh(ji,jj,jk) = (  e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm) * usd(ji  ,jj,jk)    &
207            &                 - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * usd(ji-1,jj,jk)    &
208            &                 + e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * vsd(ji,jj  ,jk)    &
209            &                 - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vsd(ji,jj-1,jk)  ) * r1_e1e2t(ji,jj)
210      END_3D
211      !
212#if defined key_agrif
213      IF( .NOT. Agrif_Root() ) THEN
214         IF( nbondi == -1 .OR. nbondi == 2 )   ze3divh( 2:nbghostcells+1,:        ,:) = 0._wp      ! west
215         IF( nbondi ==  1 .OR. nbondi == 2 )   ze3divh( nlci-nbghostcells:nlci-1,:,:) = 0._wp      ! east
216         IF( nbondj == -1 .OR. nbondj == 2 )   ze3divh( :,2:nbghostcells+1        ,:) = 0._wp      ! south
217         IF( nbondj ==  1 .OR. nbondj == 2 )   ze3divh( :,nlcj-nbghostcells:nlcj-1,:) = 0._wp      ! north
218      ENDIF
219#endif
220      !
221      CALL lbc_lnk( 'sbcwave', ze3divh, 'T', 1. )
222      !
223      IF( ln_linssh ) THEN   ;   ik = 1   ! none zero velocity through the sea surface
224      ELSE                   ;   ik = 2   ! w=0 at the surface (set one for all in sbc_wave_init)
225      ENDIF
226      DO jk = jpkm1, ik, -1          ! integrate from the bottom the hor. divergence (NB: at k=jpk w is always zero)
227         wsd(:,:,jk) = wsd(:,:,jk+1) - ze3divh(:,:,jk)
228      END DO
229      !
230      IF( ln_bdy ) THEN
231         DO jk = 1, jpkm1
232            wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:)
233         END DO
234      ENDIF
235      !                       !==  Horizontal divergence of barotropic Stokes transport  ==!
236      div_sd(:,:) = 0._wp
237      DO jk = 1, jpkm1                                 !
238        div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk)
239      END DO
240      !
241      CALL iom_put( "ustokes",  usd  )
242      CALL iom_put( "vstokes",  vsd  )
243      CALL iom_put( "wstokes",  wsd  )
244      !
245      DEALLOCATE( ze3divh )
246      DEALLOCATE( zk_t, zk_u, zk_v, zu0_sd, zv0_sd )
247      !
248   END SUBROUTINE sbc_stokes
249
250
251   SUBROUTINE sbc_wstress( )
252      !!---------------------------------------------------------------------
253      !!                     ***  ROUTINE sbc_wstress  ***
254      !!
255      !! ** Purpose :   Updates the ocean momentum modified by waves
256      !!
257      !! ** Method  : - Calculate u,v components of stress depending on stress
258      !!                model
259      !!              - Calculate the stress module
260      !!              - The wind module is not modified by waves
261      !! ** action 
262      !!---------------------------------------------------------------------
263      INTEGER  ::   jj, ji   ! dummy loop argument
264      !
265      IF( ln_tauwoc ) THEN
266         utau(:,:) = utau(:,:)*tauoc_wave(:,:)
267         vtau(:,:) = vtau(:,:)*tauoc_wave(:,:)
268         taum(:,:) = taum(:,:)*tauoc_wave(:,:)
269      ENDIF
270      !
271      IF( ln_tauw ) THEN
272         DO_2D_10_10
273            ! Stress components at u- & v-points
274            utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) )
275            vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) )
276            !
277            ! Stress module at t points
278            taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) )
279         END_2D
280         CALL lbc_lnk_multi( 'sbcwave', utau(:,:), 'U', -1. , vtau(:,:), 'V', -1. , taum(:,:) , 'T', -1. )
281      ENDIF
282      !
283   END SUBROUTINE sbc_wstress
284
285
286   SUBROUTINE sbc_wave( kt, Kmm )
287      !!---------------------------------------------------------------------
288      !!                     ***  ROUTINE sbc_wave  ***
289      !!
290      !! ** Purpose :   read wave parameters from wave model  in netcdf files.
291      !!
292      !! ** Method  : - Read namelist namsbc_wave
293      !!              - Read Cd_n10 fields in netcdf files
294      !!              - Read stokes drift 2d in netcdf files
295      !!              - Read wave number in netcdf files
296      !!              - Compute 3d stokes drift using Breivik et al.,2014
297      !!                formulation
298      !! ** action 
299      !!---------------------------------------------------------------------
300      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
301      INTEGER, INTENT(in   ) ::   Kmm  ! ocean time index
302      !!---------------------------------------------------------------------
303      !
304      IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN     !==  Neutral drag coefficient  ==!
305         CALL fld_read( kt, nn_fsbc, sf_cd )             ! read from external forcing
306         cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) * tmask(:,:,1)
307      ENDIF
308
309      IF( ln_tauwoc .AND. .NOT. cpl_tauwoc ) THEN  !==  Wave induced stress  ==!
310         CALL fld_read( kt, nn_fsbc, sf_tauwoc )         ! read wave norm stress from external forcing
311         tauoc_wave(:,:) = sf_tauwoc(1)%fnow(:,:,1) * tmask(:,:,1)
312      ENDIF
313
314      IF( ln_tauw .AND. .NOT. cpl_tauw ) THEN      !==  Wave induced stress  ==!
315         CALL fld_read( kt, nn_fsbc, sf_tauw )           ! read ocean stress components from external forcing (T grid)
316         tauw_x(:,:) = sf_tauw(1)%fnow(:,:,1) * tmask(:,:,1)
317         tauw_y(:,:) = sf_tauw(2)%fnow(:,:,1) * tmask(:,:,1)
318      ENDIF
319
320      IF( ln_sdw )  THEN                           !==  Computation of the 3d Stokes Drift  ==!
321         !
322         IF( jpfld > 0 ) THEN                            ! Read from file only if the field is not coupled
323            CALL fld_read( kt, nn_fsbc, sf_sd )          ! read wave parameters from external forcing
324            IF( jp_hsw > 0 )   hsw  (:,:) = sf_sd(jp_hsw)%fnow(:,:,1) * tmask(:,:,1)  ! significant wave height
325            IF( jp_wmp > 0 )   wmp  (:,:) = sf_sd(jp_wmp)%fnow(:,:,1) * tmask(:,:,1)  ! wave mean period
326            IF( jp_wfr > 0 )   wfreq(:,:) = sf_sd(jp_wfr)%fnow(:,:,1) * tmask(:,:,1)  ! Peak wave frequency
327            IF( jp_usd > 0 )   ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1) * tmask(:,:,1)  ! 2D zonal Stokes Drift at T point
328            IF( jp_vsd > 0 )   vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) * tmask(:,:,1)  ! 2D meridional Stokes Drift at T point
329         ENDIF
330         !
331         ! Read also wave number if needed, so that it is available in coupling routines
332         IF( ln_zdfswm .AND. .NOT.cpl_wnum ) THEN
333            CALL fld_read( kt, nn_fsbc, sf_wn )          ! read wave parameters from external forcing
334            wnum(:,:) = sf_wn(1)%fnow(:,:,1) * tmask(:,:,1)
335         ENDIF
336           
337         ! Calculate only if required fields have been read
338         ! In coupled wave model-NEMO case the call is done after coupling
339         !
340         IF( ( ll_st_bv_li   .AND. jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0 ) .OR. &
341           & ( ll_st_peakfr  .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0                ) ) CALL sbc_stokes( Kmm )
342         !
343      ENDIF
344      !
345   END SUBROUTINE sbc_wave
346
347
348   SUBROUTINE sbc_wave_init
349      !!---------------------------------------------------------------------
350      !!                     ***  ROUTINE sbc_wave_init  ***
351      !!
352      !! ** Purpose :   read wave parameters from wave model  in netcdf files.
353      !!
354      !! ** Method  : - Read namelist namsbc_wave
355      !!              - Read Cd_n10 fields in netcdf files
356      !!              - Read stokes drift 2d in netcdf files
357      !!              - Read wave number in netcdf files
358      !!              - Compute 3d stokes drift using Breivik et al.,2014
359      !!                formulation
360      !! ** action 
361      !!---------------------------------------------------------------------
362      INTEGER ::   ierror, ios   ! local integer
363      INTEGER ::   ifpr
364      !!
365      CHARACTER(len=100)     ::  cn_dir                            ! Root directory for location of drag coefficient files
366      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i, slf_j     ! array of namelist informations on the fields to read
367      TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd,  &
368                             &   sn_hsw, sn_wmp, sn_wfr, sn_wnum, &
369                             &   sn_tauwoc, sn_tauwx, sn_tauwy     ! informations about the fields to be read
370      !
371      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, &
372                             sn_wnum, sn_tauwoc, sn_tauwx, sn_tauwy
373      !!---------------------------------------------------------------------
374      !
375      READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901)
376901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc_wave in reference namelist' )
377         
378      READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 )
379902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist' )
380      IF(lwm) WRITE ( numond, namsbc_wave )
381      !
382      IF( ln_cdgw ) THEN
383         IF( .NOT. cpl_wdrag ) THEN
384            ALLOCATE( sf_cd(1), STAT=ierror )               !* allocate and fill sf_wave with sn_cdg
385            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' )
386            !
387                                   ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)   )
388            IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) )
389            CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' )
390         ENDIF
391         ALLOCATE( cdn_wave(jpi,jpj) )
392      ENDIF
393
394      IF( ln_tauwoc ) THEN
395         IF( .NOT. cpl_tauwoc ) THEN
396            ALLOCATE( sf_tauwoc(1), STAT=ierror )           !* allocate and fill sf_wave with sn_tauwoc
397            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' )
398            !
399                                     ALLOCATE( sf_tauwoc(1)%fnow(jpi,jpj,1)   )
400            IF( sn_tauwoc%ln_tint )  ALLOCATE( sf_tauwoc(1)%fdta(jpi,jpj,1,2) )
401            CALL fld_fill( sf_tauwoc, (/ sn_tauwoc /), cn_dir, 'sbc_wave_init', 'Wave module', 'namsbc_wave' )
402         ENDIF
403         ALLOCATE( tauoc_wave(jpi,jpj) )
404      ENDIF
405
406      IF( ln_tauw ) THEN
407         IF( .NOT. cpl_tauw ) THEN
408            ALLOCATE( sf_tauw(2), STAT=ierror )           !* allocate and fill sf_wave with sn_tauwx/y
409            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauw structure' )
410            !
411            ALLOCATE( slf_j(2) )
412            slf_j(1) = sn_tauwx
413            slf_j(2) = sn_tauwy
414                                    ALLOCATE( sf_tauw(1)%fnow(jpi,jpj,1)   )
415                                    ALLOCATE( sf_tauw(2)%fnow(jpi,jpj,1)   )
416            IF( slf_j(1)%ln_tint )  ALLOCATE( sf_tauw(1)%fdta(jpi,jpj,1,2) )
417            IF( slf_j(2)%ln_tint )  ALLOCATE( sf_tauw(2)%fdta(jpi,jpj,1,2) )
418            CALL fld_fill( sf_tauw, (/ slf_j /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
419         ENDIF
420         ALLOCATE( tauw_x(jpi,jpj) )
421         ALLOCATE( tauw_y(jpi,jpj) )
422      ENDIF
423
424      IF( ln_sdw ) THEN   ! Find out how many fields have to be read from file if not coupled
425         jpfld=0
426         jp_usd=0   ;   jp_vsd=0   ;   jp_hsw=0   ;   jp_wmp=0   ;   jp_wfr=0
427         IF( .NOT. cpl_sdrftx ) THEN
428            jpfld  = jpfld + 1
429            jp_usd = jpfld
430         ENDIF
431         IF( .NOT. cpl_sdrfty ) THEN
432            jpfld  = jpfld + 1
433            jp_vsd = jpfld
434         ENDIF
435         IF( .NOT. cpl_hsig  .AND. ll_st_bv_li  ) THEN
436            jpfld  = jpfld + 1
437            jp_hsw = jpfld
438         ENDIF
439         IF( .NOT. cpl_wper  .AND. ll_st_bv_li  ) THEN
440            jpfld  = jpfld + 1
441            jp_wmp = jpfld
442         ENDIF
443         IF( .NOT. cpl_wfreq .AND. ll_st_peakfr ) THEN
444            jpfld  = jpfld + 1
445            jp_wfr = jpfld
446         ENDIF
447
448         ! Read from file only the non-coupled fields
449         IF( jpfld > 0 ) THEN
450            ALLOCATE( slf_i(jpfld) )
451            IF( jp_usd > 0 )   slf_i(jp_usd) = sn_usd
452            IF( jp_vsd > 0 )   slf_i(jp_vsd) = sn_vsd
453            IF( jp_hsw > 0 )   slf_i(jp_hsw) = sn_hsw
454            IF( jp_wmp > 0 )   slf_i(jp_wmp) = sn_wmp
455            IF( jp_wfr > 0 )   slf_i(jp_wfr) = sn_wfr
456
457            ALLOCATE( sf_sd(jpfld), STAT=ierror )   !* allocate and fill sf_sd with stokes drift
458            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' )
459            !
460            DO ifpr= 1, jpfld
461               ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) )
462               IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) )
463            END DO
464            !
465            CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' )
466         ENDIF
467         ALLOCATE( usd  (jpi,jpj,jpk), vsd  (jpi,jpj,jpk), wsd(jpi,jpj,jpk) )
468         ALLOCATE( hsw  (jpi,jpj)    , wmp  (jpi,jpj)     )
469         ALLOCATE( wfreq(jpi,jpj) )
470         ALLOCATE( ut0sd(jpi,jpj)    , vt0sd(jpi,jpj)     )
471         ALLOCATE( div_sd(jpi,jpj) )
472         ALLOCATE( tsd2d (jpi,jpj) )
473
474         ut0sd(:,:) = 0._wp
475         vt0sd(:,:) = 0._wp
476         hsw(:,:) = 0._wp
477         wmp(:,:) = 0._wp
478
479         usd(:,:,:) = 0._wp
480         vsd(:,:,:) = 0._wp
481         wsd(:,:,:) = 0._wp
482         ! Wave number needed only if ln_zdfswm=T
483         IF( .NOT. cpl_wnum ) THEN
484            ALLOCATE( sf_wn(1), STAT=ierror )           !* allocate and fill sf_wave with sn_wnum
485            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wave structure' )
486                                   ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1)   )
487            IF( sn_wnum%ln_tint )  ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) )
488            CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' )
489         ENDIF
490         ALLOCATE( wnum(jpi,jpj) )
491      ENDIF
492      !
493   END SUBROUTINE sbc_wave_init
494
495   !!======================================================================
496END MODULE sbcwave
Note: See TracBrowser for help on using the repository browser.