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.
Changeset 7646 for trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90 – NEMO

Ignore:
Timestamp:
2017-02-06T10:25:03+01:00 (7 years ago)
Author:
timgraham
Message:

Merge of dev_merge_2016 into trunk. UPDATE TO ARCHFILES NEEDED for XIOS2.
LIM_SRC_s/limrhg.F90 to follow in next commit due to change of kind (I'm unable to do it in this commit).
Merged using the following steps:

1) svn merge --reintegrate svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk .
2) Resolve minor conflicts in sette.sh and namelist_cfg for ORCA2LIM3 (due to a change in trunk after branch was created)
3) svn commit
4) svn switch svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk
5) svn merge svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/2016/dev_merge_2016 .
6) At this stage I checked out a clean copy of the branch to compare against what is about to be committed to the trunk.
6) svn commit #Commit code to the trunk

In this commit I have also reverted a change to Fcheck_archfile.sh which was causing problems on the Paris machine.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90

    r6140 r7646  
    44   !! Wave module  
    55   !!====================================================================== 
    6    !! History :  3.3  !   2011-09  (Adani M)  Original code: Drag Coefficient  
    7    !!         :  3.4  !   2012-10  (Adani M)                 Stokes Drift  
    8    !!---------------------------------------------------------------------- 
    9  
    10    !!---------------------------------------------------------------------- 
    11    !!   sbc_wave      : read drag coefficient from wave model in netcdf files  
    12    !!---------------------------------------------------------------------- 
    13    USE oce            !  
     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 
    1420   USE sbc_oce        ! Surface boundary condition: ocean fields 
    15    USE bdy_oce        ! 
    16    USE domvvl         ! 
     21   USE zdf_oce,  ONLY : ln_zdfqiao 
     22   USE bdy_oce        ! open boundary condition variables 
     23   USE domvvl         ! domain: variable volume layers 
    1724   ! 
    1825   USE iom            ! I/O manager library 
     
    2532   PRIVATE 
    2633 
    27    PUBLIC   sbc_wave    ! routine called in sbc_blk_core or sbc_blk_mfs 
     34   PUBLIC   sbc_stokes      ! routine called in sbccpl 
     35   PUBLIC   sbc_wave        ! routine called in sbcmod 
     36   PUBLIC   sbc_wave_init   ! routine called in sbcmod 
    2837    
    29    INTEGER , PARAMETER ::   jpfld  = 3   ! maximum number of files to read for srokes drift 
    30    INTEGER , PARAMETER ::   jp_usd = 1   ! index of stokes drift  (i-component) (m/s)    at T-point 
    31    INTEGER , PARAMETER ::   jp_vsd = 2   ! index of stokes drift  (j-component) (m/s)    at T-point 
    32    INTEGER , PARAMETER ::   jp_wn  = 3   ! index of wave number                 (1/m)    at T-point 
    33  
    34    TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_cd    ! structure of input fields (file informations, fields read) Drag Coefficient 
    35    TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_sd    ! structure of input fields (file informations, fields read) Stokes Drift 
    36  
    37    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION (:,:)   :: cdn_wave  
    38    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION (:,:,:) :: usd3d, vsd3d, wsd3d  
    39    REAL(wp),         ALLOCATABLE, DIMENSION (:,:)   :: usd2d, vsd2d, uwavenum, vwavenum  
     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_wnum   = .FALSE. 
     45   LOGICAL, PUBLIC ::   cpl_wstrf  = .FALSE. 
     46   LOGICAL, PUBLIC ::   cpl_wdrag  = .FALSE. 
     47 
     48   INTEGER ::   jpfld    ! number of files to read for stokes drift 
     49   INTEGER ::   jp_usd   ! index of stokes drift  (i-component) (m/s)    at T-point 
     50   INTEGER ::   jp_vsd   ! index of stokes drift  (j-component) (m/s)    at T-point 
     51   INTEGER ::   jp_hsw   ! index of significant wave hight      (m)      at T-point 
     52   INTEGER ::   jp_wmp   ! index of mean wave period            (s)      at T-point 
     53 
     54   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_cd      ! structure of input fields (file informations, fields read) Drag Coefficient 
     55   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_sd      ! structure of input fields (file informations, fields read) Stokes Drift 
     56   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_wn      ! structure of input fields (file informations, fields read) wave number for Qiao 
     57   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_tauoc   ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 
     58   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   cdn_wave            !: 
     59   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   hsw, wmp, wnum      !:  
     60   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tauoc_wave          !:   
     61   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tsd2d               !:  
     62   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   div_sd              !: barotropic stokes drift divergence 
     63   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   ut0sd, vt0sd        !: surface Stokes drift velocities at t-point 
     64   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) ::   usd  , vsd  , wsd   !: Stokes drift velocities at u-, v- & w-points, resp. 
    4065 
    4166   !! * Substitutions 
     
    4873CONTAINS 
    4974 
     75   SUBROUTINE sbc_stokes( ) 
     76      !!--------------------------------------------------------------------- 
     77      !!                     ***  ROUTINE sbc_stokes  *** 
     78      !! 
     79      !! ** Purpose :   compute the 3d Stokes Drift according to Breivik et al., 
     80      !!                2014 (DOI: 10.1175/JPO-D-14-0020.1) 
     81      !! 
     82      !! ** Method  : - Calculate Stokes transport speed  
     83      !!              - Calculate horizontal divergence  
     84      !!              - Integrate the horizontal divergenze from the bottom  
     85      !! ** action   
     86      !!--------------------------------------------------------------------- 
     87      INTEGER  ::   jj, ji, jk   ! dummy loop argument 
     88      INTEGER  ::   ik           ! local integer  
     89      REAL(wp) ::  ztransp, zfac, ztemp, zsp0 
     90      REAL(wp) ::  zdep_u, zdep_v, zkh_u, zkh_v, zda_u, zda_v 
     91      REAL(wp), DIMENSION(:,:)  , POINTER ::   zk_t, zk_u, zk_v, zu0_sd, zv0_sd   ! 2D workspace 
     92      REAL(wp), DIMENSION(:,:,:), POINTER ::   ze3divh                            ! 3D workspace 
     93      !!--------------------------------------------------------------------- 
     94      ! 
     95      CALL wrk_alloc( jpi,jpj,jpk,   ze3divh ) 
     96      CALL wrk_alloc( jpi,jpj,       zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 
     97      ! 
     98      ! 
     99      zfac =  2.0_wp * rpi / 16.0_wp 
     100      DO jj = 1, jpj                ! exp. wave number at t-point    (Eq. (19) in Breivick et al. (2014) ) 
     101         DO ji = 1, jpi 
     102               ! Stokes drift velocity estimated from Hs and Tmean 
     103               ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj) , 0.0000001_wp ) 
     104               ! Stokes surface speed 
     105               zsp0 = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj) ) 
     106               tsd2d(ji,jj) = zsp0 
     107               ! Wavenumber scale 
     108               zk_t(ji,jj) = ABS( zsp0 ) / MAX( ABS( 5.97_wp*ztransp ) , 0.0000001_wp ) 
     109         END DO 
     110      END DO       
     111      DO jj = 1, jpjm1              ! exp. wave number & Stokes drift velocity at u- & v-points 
     112         DO ji = 1, jpim1 
     113            zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 
     114            zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 
     115            ! 
     116            zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 
     117            zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 
     118         END DO 
     119      END DO 
     120      ! 
     121      !                       !==  horizontal Stokes Drift 3D velocity  ==! 
     122      DO jk = 1, jpkm1 
     123         DO jj = 2, jpjm1 
     124            DO ji = 2, jpim1 
     125               zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) 
     126               zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) 
     127               !                           
     128               zkh_u = zk_u(ji,jj) * zdep_u     ! k * depth 
     129               zkh_v = zk_v(ji,jj) * zdep_v 
     130               !                                ! Depth attenuation 
     131               zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u ) 
     132               zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v ) 
     133               ! 
     134               usd(ji,jj,jk) = zda_u * zk_u(ji,jj) * umask(ji,jj,jk) 
     135               vsd(ji,jj,jk) = zda_v * zk_v(ji,jj) * vmask(ji,jj,jk) 
     136            END DO 
     137         END DO 
     138      END DO    
     139      CALL lbc_lnk( usd(:,:,:), 'U', vsd(:,:,:), 'V', -1. ) 
     140      ! 
     141      !                       !==  vertical Stokes Drift 3D velocity  ==! 
     142      ! 
     143      DO jk = 1, jpkm1               ! Horizontal e3*divergence 
     144         DO jj = 2, jpj 
     145            DO ji = fs_2, jpi 
     146               ze3divh(ji,jj,jk) = (  e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) * usd(ji  ,jj,jk)    & 
     147                  &                 - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd(ji-1,jj,jk)    & 
     148                  &                 + e1v(ji,jj  ) * e3v_n(ji,jj  ,jk) * vsd(ji,jj  ,jk)    & 
     149                  &                 - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd(ji,jj-1,jk)  ) * r1_e1e2t(ji,jj) 
     150            END DO 
     151         END DO 
     152      END DO 
     153      ! 
     154      IF( .NOT. AGRIF_Root() ) THEN 
     155         IF( nbondi ==  1 .OR. nbondi == 2 )   ze3divh(nlci-1,   :  ,:) = 0._wp      ! east 
     156         IF( nbondi == -1 .OR. nbondi == 2 )   ze3divh(  2   ,   :  ,:) = 0._wp      ! west 
     157         IF( nbondj ==  1 .OR. nbondj == 2 )   ze3divh(  :   ,nlcj-1,:) = 0._wp      ! north 
     158         IF( nbondj == -1 .OR. nbondj == 2 )   ze3divh(  :   ,  2   ,:) = 0._wp      ! south 
     159      ENDIF 
     160      ! 
     161      CALL lbc_lnk( ze3divh, 'T', 1. ) 
     162      ! 
     163      IF( ln_linssh ) THEN   ;   ik = 1   ! none zero velocity through the sea surface 
     164      ELSE                   ;   ik = 2   ! w=0 at the surface (set one for all in sbc_wave_init) 
     165      ENDIF 
     166      DO jk = jpkm1, ik, -1          ! integrate from the bottom the hor. divergence (NB: at k=jpk w is always zero) 
     167         wsd(:,:,jk) = wsd(:,:,jk+1) - ze3divh(:,:,jk) 
     168      END DO 
     169      ! 
     170      IF( ln_bdy ) THEN 
     171         DO jk = 1, jpkm1 
     172            wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:) 
     173         END DO 
     174      ENDIF 
     175      !                       !==  Horizontal divergence of barotropic Stokes transport  ==! 
     176      div_sd(:,:) = 0._wp 
     177      DO jk = 1, jpkm1                                 !  
     178        div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk) 
     179      END DO 
     180      ! 
     181      CALL iom_put( "ustokes",  usd  ) 
     182      CALL iom_put( "vstokes",  vsd  ) 
     183      CALL iom_put( "wstokes",  wsd  ) 
     184      ! 
     185      CALL wrk_dealloc( jpi,jpj,jpk,   ze3divh ) 
     186      CALL wrk_dealloc( jpi,jpj,       zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 
     187      ! 
     188   END SUBROUTINE sbc_stokes 
     189 
     190 
    50191   SUBROUTINE sbc_wave( kt ) 
    51192      !!--------------------------------------------------------------------- 
    52       !!                     ***  ROUTINE sbc_apr  *** 
    53       !! 
    54       !! ** Purpose :   read drag coefficient from wave model  in netcdf files. 
     193      !!                     ***  ROUTINE sbc_wave  *** 
     194      !! 
     195      !! ** Purpose :   read wave parameters from wave model  in netcdf files. 
    55196      !! 
    56197      !! ** Method  : - Read namelist namsbc_wave 
    57198      !!              - Read Cd_n10 fields in netcdf files  
    58199      !!              - Read stokes drift 2d in netcdf files  
    59       !!              - Read wave number      in netcdf files  
    60       !!              - Compute 3d stokes drift using monochromatic 
    61       !! ** action  :    
    62       !!--------------------------------------------------------------------- 
    63       INTEGER, INTENT( in  ) ::   kt       ! ocean time step 
    64       ! 
    65       INTEGER                ::   ierror   ! return error code 
    66       INTEGER                ::   ifpr, jj,ji,jk  
    67       INTEGER                ::   ios     ! Local integer output status for namelist read 
    68       TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i     ! array of namelist informations on the fields to read 
     200      !!              - Read wave number in netcdf files  
     201      !!              - Compute 3d stokes drift using Breivik et al.,2014 
     202      !!                formulation 
     203      !! ** action   
     204      !!--------------------------------------------------------------------- 
     205      INTEGER, INTENT(in   ) ::   kt   ! ocean time step 
     206      !!--------------------------------------------------------------------- 
     207      ! 
     208      IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN     !==  Neutral drag coefficient  ==! 
     209         CALL fld_read( kt, nn_fsbc, sf_cd )             ! read from external forcing 
     210         cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 
     211      ENDIF 
     212 
     213      IF( ln_tauoc .AND. .NOT. cpl_wstrf ) THEN    !==  Wave induced stress  ==! 
     214         CALL fld_read( kt, nn_fsbc, sf_tauoc )          ! read wave norm stress from external forcing 
     215         tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 
     216      ENDIF 
     217 
     218      IF( ln_sdw )  THEN                           !==  Computation of the 3d Stokes Drift  ==!  
     219         ! 
     220         IF( jpfld > 0 ) THEN                            ! Read from file only if the field is not coupled 
     221            CALL fld_read( kt, nn_fsbc, sf_sd )          ! read wave parameters from external forcing 
     222            IF( jp_hsw > 0 )   hsw  (:,:) = sf_sd(jp_hsw)%fnow(:,:,1)   ! significant wave height 
     223            IF( jp_wmp > 0 )   wmp  (:,:) = sf_sd(jp_wmp)%fnow(:,:,1)   ! wave mean period 
     224            IF( jp_usd > 0 )   ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1)   ! 2D zonal Stokes Drift at T point 
     225            IF( jp_vsd > 0 )   vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1)   ! 2D meridional Stokes Drift at T point 
     226         ENDIF 
     227         ! 
     228         ! Read also wave number if needed, so that it is available in coupling routines 
     229         IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 
     230            CALL fld_read( kt, nn_fsbc, sf_wn )          ! read wave parameters from external forcing 
     231            wnum(:,:) = sf_wn(1)%fnow(:,:,1) 
     232         ENDIF 
     233            
     234         !                                         !==  Computation of the 3d Stokes Drift  ==!  
     235         ! 
     236         IF( jpfld == 4 )   CALL sbc_stokes()            ! Calculate only if required fields are read 
     237         !                                               ! In coupled wave model-NEMO case the call is done after coupling 
     238         ! 
     239      ENDIF 
     240      ! 
     241   END SUBROUTINE sbc_wave 
     242 
     243 
     244   SUBROUTINE sbc_wave_init 
     245      !!--------------------------------------------------------------------- 
     246      !!                     ***  ROUTINE sbc_wave_init  *** 
     247      !! 
     248      !! ** Purpose :   read wave parameters from wave model  in netcdf files. 
     249      !! 
     250      !! ** Method  : - Read namelist namsbc_wave 
     251      !!              - Read Cd_n10 fields in netcdf files  
     252      !!              - Read stokes drift 2d in netcdf files  
     253      !!              - Read wave number in netcdf files  
     254      !!              - Compute 3d stokes drift using Breivik et al.,2014 
     255      !!                formulation 
     256      !! ** action   
     257      !!--------------------------------------------------------------------- 
     258      INTEGER ::   ierror, ios   ! local integer 
     259      INTEGER ::   ifpr 
     260      !! 
    69261      CHARACTER(len=100)     ::  cn_dir                          ! Root directory for location of drag coefficient files 
    70       TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd, sn_wn   ! informations about the fields to be read 
    71       REAL(wp), DIMENSION(:,:,:), POINTER ::   zusd_t, zvsd_t, ze3hdiv   ! 3D workspace 
    72       !! 
    73       NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_wn, ln_cdgw , ln_sdw 
    74       !!--------------------------------------------------------------------- 
    75       ! 
    76       !                                         ! -------------------- ! 
    77       IF( kt == nit000 ) THEN                   ! First call kt=nit000 ! 
    78          !                                      ! -------------------- ! 
    79          REWIND( numnam_ref )              ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 
    80          READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 
    81 901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 
    82          ! 
    83          REWIND( numnam_cfg )              ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 
    84          READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 
    85 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 
    86          IF(lwm) WRITE ( numond, namsbc_wave ) 
    87          ! 
    88          IF(lwp) THEN               ! Control print 
    89             WRITE(numout,*) '        Namelist namsbc_wave : surface wave setting'  
    90             WRITE(numout,*) '           wave drag coefficient                      ln_cdgw  = ', ln_cdgw   
    91             WRITE(numout,*) '           wave stokes drift                          ln_sdw   = ', ln_sdw 
    92          ENDIF 
    93          ! 
    94          IF( .NOT.( ln_cdgw .OR. ln_sdw ) )    & 
    95             &  CALL ctl_warn( 'ln_sbcwave=T but nor drag coefficient (ln_cdgw=F) neither stokes drift activated (ln_sdw=F)' ) 
    96          IF( ln_cdgw .AND. .NOT.(ln_blk_mfs .OR. ln_blk_core) )   &        
    97             &  CALL ctl_stop( 'drag coefficient read from wave model definable only with mfs bulk formulae and core') 
    98          ! 
    99          IF( ln_cdgw ) THEN 
     262      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i     ! array of namelist informations on the fields to read 
     263      TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd,  & 
     264                             &   sn_hsw, sn_wmp, sn_wnum, sn_tauoc      ! informations about the fields to be read 
     265      ! 
     266      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wnum, sn_tauoc 
     267      !!--------------------------------------------------------------------- 
     268      ! 
     269      REWIND( numnam_ref )              ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 
     270      READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 
     271901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 
     272          
     273      REWIND( numnam_cfg )              ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 
     274      READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 
     275902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 
     276      IF(lwm) WRITE ( numond, namsbc_wave ) 
     277      ! 
     278      IF( ln_cdgw ) THEN 
     279         IF( .NOT. cpl_wdrag ) THEN 
    100280            ALLOCATE( sf_cd(1), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
    101             IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     281            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 
    102282            ! 
    103283                                   ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)   ) 
    104284            IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 
    105             CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    106             ALLOCATE( cdn_wave(jpi,jpj) ) 
    107             cdn_wave(:,:) = 0.0 
    108          ENDIF 
    109          IF( ln_sdw ) THEN 
    110             slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; slf_i(jp_wn) = sn_wn 
    111             ALLOCATE( sf_sd(3), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
    112             IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     285            CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 
     286         ENDIF 
     287         ALLOCATE( cdn_wave(jpi,jpj) ) 
     288      ENDIF 
     289 
     290      IF( ln_tauoc ) THEN 
     291         IF( .NOT. cpl_wstrf ) THEN 
     292            ALLOCATE( sf_tauoc(1), STAT=ierror )           !* allocate and fill sf_wave with sn_tauoc 
     293            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 
     294            ! 
     295                                    ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1)   ) 
     296            IF( sn_tauoc%ln_tint )  ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 
     297            CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave_init', 'Wave module', 'namsbc_wave' ) 
     298         ENDIF 
     299         ALLOCATE( tauoc_wave(jpi,jpj) ) 
     300      ENDIF 
     301 
     302      IF( ln_sdw ) THEN   ! Find out how many fields have to be read from file if not coupled 
     303         jpfld=0 
     304         jp_usd=0   ;   jp_vsd=0   ;   jp_hsw=0   ;   jp_wmp=0 
     305         IF( .NOT. cpl_sdrftx ) THEN 
     306            jpfld  = jpfld + 1 
     307            jp_usd = jpfld 
     308         ENDIF 
     309         IF( .NOT. cpl_sdrfty ) THEN 
     310            jpfld  = jpfld + 1 
     311            jp_vsd = jpfld 
     312         ENDIF 
     313         IF( .NOT. cpl_hsig ) THEN 
     314            jpfld  = jpfld + 1 
     315            jp_hsw = jpfld 
     316         ENDIF 
     317         IF( .NOT. cpl_wper ) THEN 
     318            jpfld  = jpfld + 1 
     319            jp_wmp = jpfld 
     320         ENDIF 
     321 
     322         ! Read from file only the non-coupled fields  
     323         IF( jpfld > 0 ) THEN 
     324            ALLOCATE( slf_i(jpfld) ) 
     325            IF( jp_usd > 0 )   slf_i(jp_usd) = sn_usd 
     326            IF( jp_vsd > 0 )   slf_i(jp_vsd) = sn_vsd 
     327            IF( jp_hsw > 0 )   slf_i(jp_hsw) = sn_hsw 
     328            IF( jp_wmp > 0 )   slf_i(jp_wmp) = sn_wmp 
     329            ALLOCATE( sf_sd(jpfld), STAT=ierror )   !* allocate and fill sf_sd with stokes drift 
     330            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 
    113331            ! 
    114332            DO ifpr= 1, jpfld 
     
    116334               IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 
    117335            END DO 
    118             CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    119             ALLOCATE( usd2d(jpi,jpj) , vsd2d(jpi,jpj) , uwavenum(jpi,jpj) , vwavenum(jpi,jpj) ) 
    120             ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 
    121             usd3d(:,:,:) = 0._wp   ;   usd2d(:,:) = 0._wp   ;    uwavenum(:,:) = 0._wp 
    122             vsd3d(:,:,:) = 0._wp   ;   vsd2d(:,:) = 0._wp   ;    vwavenum(:,:) = 0._wp 
    123             wsd3d(:,:,:) = 0._wp 
    124          ENDIF 
    125       ENDIF 
    126       ! 
    127       IF( ln_cdgw ) THEN               !==  Neutral drag coefficient  ==! 
    128          CALL fld_read( kt, nn_fsbc, sf_cd )      ! read from external forcing 
    129          cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 
    130       ENDIF 
    131       ! 
    132       IF( ln_sdw )  THEN               !==  Computation of the 3d Stokes Drift  ==! 
    133          ! 
    134          CALL wrk_alloc( jpi,jpj,jpk,   zusd_t, zvsd_t, ze3hdiv ) 
    135          ! 
    136          CALL fld_read( kt, nn_fsbc, sf_sd )    !* read drag coefficient from external forcing 
    137          ! 
    138          DO jk = 1, jpkm1                       !* distribute it on the vertical 
    139             zusd_t(:,:,jk) = sf_sd(jp_usd)%fnow(:,:,1) * EXP( -2._wp * sf_sd(jp_wn)%fnow(:,:,1) * gdept_n(:,:,jk) ) 
    140             zvsd_t(:,:,jk) = sf_sd(jp_vsd)%fnow(:,:,1) * EXP( -2._wp * sf_sd(jp_wn)%fnow(:,:,1) * gdept_n(:,:,jk) ) 
    141          END DO 
    142          DO jk = 1, jpkm1                       !* interpolate the stokes drift from t-point to u- and v-points 
    143             DO jj = 1, jpjm1 
    144                DO ji = 1, jpim1 
    145                    usd3d(ji,jj,jk) = 0.5_wp * ( zusd_t(ji  ,jj,jk) + zusd_t(ji+1,jj,jk) ) * umask(ji,jj,jk) 
    146                    vsd3d(ji,jj,jk) = 0.5_wp * ( zvsd_t(ji  ,jj,jk) + zvsd_t(ji,jj+1,jk) ) * vmask(ji,jj,jk) 
    147                END DO 
    148             END DO 
    149          END DO 
    150          CALL lbc_lnk( usd3d(:,:,:), 'U', -1. ) 
    151          CALL lbc_lnk( vsd3d(:,:,:), 'V', -1. ) 
    152          ! 
    153          DO jk = 1, jpkm1                       !* e3t * Horizontal divergence  ==! 
    154             DO jj = 2, jpjm1 
    155                DO ji = fs_2, fs_jpim1   ! vector opt. 
    156                   ze3hdiv(ji,jj,jk) = (  e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) * usd3d(ji  ,jj,jk)     & 
    157                      &                 - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd3d(ji-1,jj,jk)     & 
    158                      &                 + e1v(ji,jj  ) * e3v_n(ji,jj  ,jk) * vsd3d(ji,jj  ,jk)     & 
    159                      &                 - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd3d(ji,jj-1,jk)   ) * r1_e1e2t(ji,jj) 
    160                END DO   
    161             END DO   
    162             IF( .NOT. AGRIF_Root() ) THEN 
    163                IF( nbondi ==  1 .OR. nbondi == 2 )   ze3hdiv(nlci-1,   :  ,jk) = 0._wp      ! east 
    164                IF( nbondi == -1 .OR. nbondi == 2 )   ze3hdiv(  2   ,   :  ,jk) = 0._wp      ! west 
    165                IF( nbondj ==  1 .OR. nbondj == 2 )   ze3hdiv(  :   ,nlcj-1,jk) = 0._wp      ! north 
    166                IF( nbondj == -1 .OR. nbondj == 2 )   ze3hdiv(  :   ,  2   ,jk) = 0._wp      ! south 
    167             ENDIF 
    168          END DO 
    169          CALL lbc_lnk( ze3hdiv, 'T', 1. )  
    170          ! 
    171          DO jk = jpkm1, 1, -1                   !* integrate from the bottom the e3t * hor. divergence 
    172             wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - ze3hdiv(:,:,jk) 
    173          END DO 
    174 #if defined key_bdy 
    175          IF( lk_bdy ) THEN 
    176             DO jk = 1, jpkm1 
    177                wsd3d(:,:,jk) = wsd3d(:,:,jk) * bdytmask(:,:) 
    178             END DO 
    179          ENDIF 
    180 #endif 
    181          CALL wrk_dealloc( jpi,jpj,jpk,   zusd_t, zvsd_t, ze3hdiv ) 
    182          !  
    183       ENDIF 
    184       ! 
    185    END SUBROUTINE sbc_wave 
    186        
     336            ! 
     337            CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 
     338         ENDIF 
     339         ALLOCATE( usd  (jpi,jpj,jpk), vsd  (jpi,jpj,jpk), wsd(jpi,jpj,jpk) ) 
     340         ALLOCATE( hsw  (jpi,jpj)    , wmp  (jpi,jpj)     ) 
     341         ALLOCATE( ut0sd(jpi,jpj)    , vt0sd(jpi,jpj)     ) 
     342         ALLOCATE( div_sd(jpi,jpj) ) 
     343         ALLOCATE( tsd2d (jpi,jpj) ) 
     344         usd(:,:,:) = 0._wp 
     345         vsd(:,:,:) = 0._wp 
     346         wsd(:,:,:) = 0._wp 
     347         ! Wave number needed only if ln_zdfqiao=T 
     348         IF( .NOT. cpl_wnum ) THEN 
     349            ALLOCATE( sf_wn(1), STAT=ierror )           !* allocate and fill sf_wave with sn_wnum 
     350            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wave structure' ) 
     351                                   ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1)   ) 
     352            IF( sn_wnum%ln_tint )  ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 
     353            CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 
     354         ENDIF 
     355         ALLOCATE( wnum(jpi,jpj) ) 
     356      ENDIF 
     357      ! 
     358   END SUBROUTINE sbc_wave_init 
     359 
    187360   !!====================================================================== 
    188361END MODULE sbcwave 
Note: See TracChangeset for help on using the changeset viewer.