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 10473 for branches/UKMO – NEMO

Changeset 10473 for branches/UKMO


Ignore:
Timestamp:
2019-01-08T18:02:36+01:00 (5 years ago)
Author:
jcastill
Message:

Merged branch UKMO/r6232_INGV1_WAVE-coupling@7620

Location:
branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM
Files:
1 added
20 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/CONFIG/SHARED/field_def.xml

    r5517 r10473  
    374374         <field id="uocet"        long_name="ocean transport along i-axis times temperature (CRS)"                                               unit="degC*m/s"   grid_ref="grid_U_3D" /> 
    375375         <field id="uoces"        long_name="ocean transport along i-axis times salinity (CRS)"                                                  unit="1e-3*m/s"   grid_ref="grid_U_3D" /> 
     376    <field id="ustokes"      long_name="Stokes Drift Velocity i-axis"                           standard_name="StokesDrift_x_velocity"      unit="m/s"        grid_ref="grid_U_3D" />   
     377    <field id="ustokes_e3u"  long_name="Stokes Drift Velocity i-axis  (thickness weighted)"                                                 unit="m/s"        grid_ref="grid_U_3D"  > ustokes * e3u </field> 
    376378 
    377379         <!-- variables available with MLE --> 
     
    409411         <field id="vocet"        long_name="ocean transport along j-axis times temperature (CRS)"                                               unit="degC*m/s"   grid_ref="grid_V_3D" /> 
    410412         <field id="voces"        long_name="ocean transport along j-axis times salinity (CRS)"                                                  unit="1e-3*m/s"   grid_ref="grid_V_3D" /> 
     413    <field id="vstokes"      long_name="Stokes Drift Velocity j-axis"                           standard_name="StokesDrift_y_velocity"      unit="m/s"        grid_ref="grid_V_3D" />   
     414    <field id="vstokes_e3v"  long_name="Stokes Drift Velocity j-axis  (thickness weighted)"                                                 unit="m/s"        grid_ref="grid_V_3D"  > vstokes * e3v </field> 
    411415 
    412416         <!-- variables available with MLE --> 
     
    437441        <field id="woce"         long_name="ocean vertical velocity"              standard_name="upward_sea_water_velocity"   unit="m/s"  /> 
    438442        <field id="wocetr_eff"   long_name="effective ocean vertical transport"                                               unit="m3/s" /> 
     443        <field id="wstokes"      long_name="Stokes Drift vertical velocity"       standard_name="upward_StokesDrift_velocity" unit="m/s" /> 
    439444 
    440445        <!-- woce_eiv: available with key_traldf_eiv and key_diaeiv --> 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/CONFIG/SHARED/namelist_ref

    r10268 r10473  
    264264                           !     =1 global mean of e-p-r set to zero at each time step 
    265265                           !     =2 annual global mean of e-p-r set to zero 
    266    ln_wave = .false.       !  Activate coupling with wave (either Stokes Drift or Drag coefficient, or both)  (T => fill namsbc_wave) 
    267    ln_cdgw = .false.       !  Neutral drag coefficient read from wave model (T => fill namsbc_wave) 
    268    ln_sdw  = .false.       !  Computation of 3D stokes drift                (T => fill namsbc_wave) 
     266   ln_wave = .false.       !  Activate coupling with wave (T => fill namsbc_wave)    
     267   ln_cdgw = .false.       !  Neutral drag coefficient read from wave model (T => ln_wave=.true. & fill namsbc_wave)    
     268   ln_sdw  = .false.       !  Read 2D Surf Stokes Drift & Computation of 3D stokes drift (T => ln_wave=.true. & fill namsbc_wave)     
     269   ln_tauoc= .false.       !  Activate ocean stress modified by external wave induced stress (T => ln_wave=.true. & fill namsbc_wave)    
     270   ln_stcor= .false.       !  Activate Stokes Coriolis term (T => ln_wave=.true. & ln_sdw=.true. & fill namsbc_wave 
    269271   nn_lsm  = 0             !  =0 land/sea mask for input fields is not applied (keep empty land/sea mask filename field) , 
    270272                           !  =1:n number of iterations of land/sea mask application for input fields (fill land/sea mask filename field) 
     
    363365   sn_snd_crt    =       'none'                 ,    'no'    , 'spherical' , 'eastward-northward' ,  'T' 
    364366   sn_snd_co2    =       'coupled'              ,    'no'    ,     ''      ,         ''           ,   '' 
     367   sn_snd_crtw   =       'none'                 ,    'no'    ,     ''      ,         ''           , 'U,V'   
     368   sn_snd_ifrac  =       'none'                 ,    'no'    ,     ''      ,         ''           ,   ''   
     369   sn_snd_wlev   =       'none'                 ,    'no'    ,     ''      ,         ''           ,   '' 
    365370! receive 
    366371   sn_rcv_w10m   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     
    374379   sn_rcv_cal    =       'coupled'              ,    'no'    ,     ''      ,         ''          ,   '' 
    375380   sn_rcv_co2    =       'coupled'              ,    'no'    ,     ''      ,         ''          ,   '' 
     381   sn_rcv_hsig   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     382   sn_rcv_iceflx =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     383   sn_rcv_mslp   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     384   sn_rcv_phioc  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     385   sn_rcv_sdrfx  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     386   sn_rcv_sdrfy  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     387   sn_rcv_wper   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     388   sn_rcv_wnum   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     389   sn_rcv_wstrf  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''   
     390   sn_rcv_wdrag  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   ''  
    376391! 
    377392   nn_cplmodel   =     1     !  Maximum number of models to/from which NEMO is potentialy sending/receiving data 
     
    921936   ln_zdfexp   = .false.   !  time-stepping: split-explicit (T) or implicit (F) time stepping 
    922937   nn_zdfexp   =    3            !  number of sub-timestep for ln_zdfexp=T 
     938   ln_zdfqiao  = .false.   !  Enhanced wave vertical mixing Qiao (2010) (T => ln_wave=.true. & ln_sdw=.true. & fill namsbc_wave) 
    923939/ 
    924940!----------------------------------------------------------------------- 
     
    9881004   rn_hsro       =  0.02   !  Minimum surface roughness 
    9891005   rn_frac_hs    =   1.3   !  Fraction of wave height as roughness (if nn_z0_met=2) 
    990    nn_z0_met     =     2   !  Method for surface roughness computation (0/1/2) 
     1006   nn_z0_met     =     2   !  Method for surface roughness computation (0/1/2/3)  
     1007      !                             ! =3 requires ln_wave=T  
    9911008   nn_bc_surf    =     1   !  surface condition (0/1=Dir/Neum) 
    9921009   nn_bc_bot     =     1   !  bottom condition (0/1=Dir/Neum) 
     
    12731290/ 
    12741291!----------------------------------------------------------------------- 
    1275 &namsbc_wave   ! External fields from wave model 
     1292&namsbc_wave   ! External fields from wave model                        (ln_wave=T) 
    12761293!----------------------------------------------------------------------- 
    12771294!              !  file name  ! frequency (hours) ! variable     ! time interp. !  clim   ! 'yearly'/ ! weights  ! rotation ! land/sea mask ! 
    12781295!              !             !  (if <0  months)  !   name       !   (logical)  !  (T/F)  ! 'monthly' ! filename ! pairing  ! filename      ! 
    1279    sn_cdg      =  'cdg_wave' ,        1          , 'drag_coeff' ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
     1296   sn_cdg      =  'cdw_wave' ,        1          , 'drag_coeff' ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    12801297   sn_usd      =  'sdw_wave' ,        1          , 'u_sd2d'     ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    12811298   sn_vsd      =  'sdw_wave' ,        1          , 'v_sd2d'     ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    1282    sn_wn       =  'sdw_wave' ,        1          , 'wave_num'   ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    1283 ! 
    1284    cn_dir_cdg  = './'  !  root directory for the location of drag coefficient files 
     1299   sn_hsw      =  'sdw_wave' ,        1          , 'hs'         ,     .true.   , .false. , 'daily'   ,  ''      , ''       , ''   
     1300   sn_wmp      =  'sdw_wave' ,        1          , 'wmp'        ,     .true.   , .false. , 'daily'   ,  ''      , ''       , ''   
     1301   sn_wnum     =  'sdw_wave' ,        1          , 'wave_num'   ,     .true.   , .false. , 'daily'   ,  ''      , ''       , ''   
     1302   sn_tauoc    =  'sdw_wave' ,        1          , 'wave_stress',     .true.   , .false. , 'daily'   ,  ''      , ''       , ''   
     1303 
     1304   cn_dir  = './'  !  root directory for the location of drag coefficient files 
    12851305/ 
    12861306!----------------------------------------------------------------------- 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/CONFIG/cfg.txt

    r8058 r10473  
    1111ORCA2_LIM OPA_SRC LIM_SRC_2 NST_SRC 
    1212ORCA2_OFF_PISCES OPA_SRC OFF_SRC TOP_SRC 
     13GYRE_LONG OPA_SRC 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r8561 r10473  
    4141   USE zdfddm          ! vertical  physics: double diffusion 
    4242   USE diahth          ! thermocline diagnostics 
     43   USE sbcwave         ! wave parameters 
    4344   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    4445   USE in_out_manager  ! I/O manager 
     
    709710         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un 
    710711            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
     712         IF( ln_wave .AND. ln_sdw) THEN   
     713            CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current"         , "m/s"    ,   &  ! usd   
     714               &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout )   
     715         ENDIF  
    711716         IF( ln_traldf_gdia ) THEN 
    712717            CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv 
     
    727732         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn 
    728733            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
     734         IF( ln_wave .AND. ln_sdw) THEN   
     735            CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current"    , "m/s"    ,   &  ! vsd   
     736               &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout )   
     737         ENDIF   
    729738         IF( ln_traldf_gdia ) THEN 
    730739            CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv 
     
    763772               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
    764773         ENDIF 
     774            
     775         IF( ln_wave .AND. ln_sdw) THEN   
     776            CALL histdef( nid_W, "sdvecrtz", "Stokes Drift Vertical Current"   , "m/s"    ,   &  ! wsd   
     777               &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout )   
     778         ENDIF   
    765779         !                                                                                      !!! nid_W : 2D 
    766780#if defined key_traldf_c2d 
     
    943957#endif 
    944958 
     959      IF( ln_wave .AND. ln_sdw ) THEN   
     960         CALL histwrite( nid_U, "sdzocrtx", it, usd           , ndim_U , ndex_U )    ! i-StokesDrift-current   
     961         CALL histwrite( nid_V, "sdmecrty", it, vsd           , ndim_V , ndex_V )    ! j-StokesDrift-current   
     962         CALL histwrite( nid_W, "sdvecrtz", it, wsd           , ndim_T , ndex_T )    ! StokesDrift vert. current   
     963      ENDIF   
     964        
    945965      ! 3. Close all files 
    946966      ! --------------------------------------- 
     
    10481068            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
    10491069      END IF 
     1070      !   
     1071      IF( ln_wave .AND. ln_sdw ) THEN   
     1072         CALL histdef( id_i, "sdzocrtx", "Stokes Drift Zonal", "m/s"    , &   ! StokesDrift zonal current   
     1073            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )   
     1074         CALL histdef( id_i, "sdmecrty", "Stokes Drift Merid", "m/s"    , &   ! StokesDrift meridonal current   
     1075            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )   
     1076         CALL histdef( id_i, "sdvecrtz", "Stokes Drift Vert", "m/s"    , &   ! StokesDrift vertical current   
     1077            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )   
     1078      ENDIF   
    10501079 
    10511080#if defined key_lim2 
     
    10811110         CALL histwrite( id_i, "vovvle3t", kt, fse3t_n (:,:,:), jpi*jpj*jpk, idex )!  T-cell thickness   
    10821111      END IF 
     1112        
     1113      IF( ln_wave .AND. ln_sdw ) THEN   
     1114         CALL histwrite( id_i, "sdzocrtx", kt, usd           , jpi*jpj*jpk, idex)     ! now StokesDrift i-velocity   
     1115         CALL histwrite( id_i, "sdmecrty", kt, vsd           , jpi*jpj*jpk, idex)     ! now StokesDrift j-velocity   
     1116         CALL histwrite( id_i, "sdvecrtz", kt, wsd           , jpi*jpj*jpk, idex)     ! now StokesDrift k-velocity   
     1117      ENDIF   
    10831118 
    10841119      ! 3. Close the file 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r10268 r10473  
    1111   !!             3.5  ! 2013-07  (J. Chanut) Switch to Forward-backward time stepping 
    1212   !!             3.6  ! 2013-11  (A. Coward) Update for z-tilde compatibility 
     13   !!              -   ! 2016-12  (G. Madec, E. Clementi) update for Stoke-Drift divergence 
    1314   !!--------------------------------------------------------------------- 
    1415#if defined key_dynspg_ts   ||   defined key_esopa 
     
    3132   USE sbctide         ! tides 
    3233   USE updtide         ! tide potential 
     34   USE sbcwave         ! surface wave  
     35   !  
    3336   USE lib_mpp         ! distributed memory computing library 
    3437   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
     
    333336         END DO 
    334337      END DO 
     338               
     339      !!gm  Question here when removing the Vertically integrated trends, we remove the vertically integrated NL trends on momentum....   
     340      !!gm  Is it correct to do so ?   I think so...   
     341               
     342               
    335343      !                                   !* barotropic Coriolis trends (vorticity scheme dependent) 
    336344      !                                   ! -------------------------------------------------------- 
     
    462470                &                        + fwfisf(:,:) + fwfisf_b(:,:)                     ) 
    463471      ENDIF 
     472      !   
     473      IF( ln_sdw ) THEN                         ! Stokes drift divergence added if necessary   
     474         zssh_frc(:,:) = zssh_frc(:,:) + div_sd(:,:)   
     475      ENDIF   
     476      !   
    464477#if defined key_asminc 
    465478      !                                         ! Include the IAU weighted SSH increment 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/DYN/dynvor.F90

    r8058 r10473  
    1616   !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase 
    1717   !!            3.7  ! 2014-04  (G. Madec) trend simplification: suppress jpdyn_trd_dat vorticity  
     18   !!             -   ! 2016-12  (G. Madec, E. Clementi) add Stokes-Coriolis trends (ln_stcor=T) 
    1819   !!---------------------------------------------------------------------- 
    1920 
     
    3233   USE trd_oce        ! trends: ocean variables 
    3334   USE trddyn         ! trend manager: dynamics 
     35   USE sbcwave        ! Surface Waves (add Stokes-Coriolis force)  
     36   USE sbc_oce , ONLY : ln_stcor    ! use Stoke-Coriolis force  
     37   !  
    3438   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    3539   USE prtctl         ! Print control 
     
    9195      ! 
    9296      CASE ( -1 )                                      ! esopa: test all possibility with control print 
    93          CALL vor_ene( kt, ntot, ua, va ) 
     97         CALL vor_ene( kt, ntot, un, vn, ua, va ) 
    9498         CALL prt_ctl( tab3d_1=ua, clinfo1=' vor0 - Ua: ', mask1=umask, & 
    9599            &          tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    96          CALL vor_ens( kt, ntot, ua, va ) 
     100         CALL vor_ens( kt, ntot, un, vn, ua, va ) 
    97101         CALL prt_ctl( tab3d_1=ua, clinfo1=' vor1 - Ua: ', mask1=umask, & 
    98102            &          tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
     
    100104         CALL prt_ctl( tab3d_1=ua, clinfo1=' vor2 - Ua: ', mask1=umask, & 
    101105            &          tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    102          CALL vor_een( kt, ntot, ua, va ) 
     106         CALL vor_een( kt, ntot, un, vn, ua, va ) 
    103107         CALL prt_ctl( tab3d_1=ua, clinfo1=' vor3 - Ua: ', mask1=umask, & 
    104108            &          tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
     
    108112            ztrdu(:,:,:) = ua(:,:,:) 
    109113            ztrdv(:,:,:) = va(:,:,:) 
    110             CALL vor_ene( kt, nrvm, ua, va )                ! relative vorticity or metric trend 
     114            CALL vor_ene( kt, nrvm, un, vn, ua, va )        ! relative vorticity or metric trend 
    111115            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    112116            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
     
    114118            ztrdu(:,:,:) = ua(:,:,:) 
    115119            ztrdv(:,:,:) = va(:,:,:) 
    116             CALL vor_ene( kt, ncor, ua, va )                ! planetary vorticity trend 
     120            CALL vor_ene( kt, ncor, un, vn, ua, va )        ! planetary vorticity trend 
    117121            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    118122            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    119123            CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt ) 
    120124         ELSE 
    121             CALL vor_ene( kt, ntot, ua, va )                ! total vorticity 
     125                             CALL vor_ene( kt, ntot, un, vn, ua, va )    ! total vorticity trend  
     126            IF( ln_stcor )   CALL vor_ene( kt, ncor, usd, vsd, ua, va )  ! add the Stokes-Coriolis trend  
    122127         ENDIF 
    123128         ! 
     
    126131            ztrdu(:,:,:) = ua(:,:,:) 
    127132            ztrdv(:,:,:) = va(:,:,:) 
    128             CALL vor_ens( kt, nrvm, ua, va )                ! relative vorticity or metric trend 
     133            CALL vor_ens( kt, nrvm, un, vn, ua, va )        ! relative vorticity or metric trend 
    129134            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    130135            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
     
    132137            ztrdu(:,:,:) = ua(:,:,:) 
    133138            ztrdv(:,:,:) = va(:,:,:) 
    134             CALL vor_ens( kt, ncor, ua, va )                ! planetary vorticity trend 
     139            CALL vor_ens( kt, ncor, un, vn, ua, va )        ! planetary vorticity trend 
    135140            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    136141            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    137142            CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt ) 
    138143         ELSE 
    139             CALL vor_ens( kt, ntot, ua, va )                ! total vorticity 
     144                             CALL vor_ens( kt, ntot, un, vn, ua, va )    ! total vorticity  
     145            IF( ln_stcor )   CALL vor_ens( kt, ncor, usd, vsd, ua, va )  ! add the Stokes-Coriolis trend 
    140146         ENDIF 
    141147         ! 
     
    144150            ztrdu(:,:,:) = ua(:,:,:) 
    145151            ztrdv(:,:,:) = va(:,:,:) 
    146             CALL vor_ens( kt, nrvm, ua, va )                ! relative vorticity or metric trend (ens) 
     152            CALL vor_ens( kt, nrvm, un, vn, ua, va )        ! relative vorticity or metric trend (ens) 
    147153            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    148154            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
     
    150156            ztrdu(:,:,:) = ua(:,:,:) 
    151157            ztrdv(:,:,:) = va(:,:,:) 
    152             CALL vor_ene( kt, ncor, ua, va )                ! planetary vorticity trend (ene) 
     158            CALL vor_ene( kt, ncor, un, vn, ua, va )        ! planetary vorticity trend (ene) 
    153159            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    154160            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    155161            CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt ) 
    156162         ELSE 
    157             CALL vor_mix( kt )                               ! total vorticity (mix=ens-ene) 
     163                             CALL vor_ens( kt, nrvm, un , vn , ua, va )   ! relative vorticity or metric trend (ens)  
     164                             CALL vor_ene( kt, ncor, un , vn , ua, va )   ! planetary vorticity trend (ene)  
     165            IF( ln_stcor )   CALL vor_ene( kt, ncor, usd, vsd, ua, va )   ! add the Stokes-Coriolis trend 
    158166         ENDIF 
    159167         ! 
     
    162170            ztrdu(:,:,:) = ua(:,:,:) 
    163171            ztrdv(:,:,:) = va(:,:,:) 
    164             CALL vor_een( kt, nrvm, ua, va )                ! relative vorticity or metric trend 
     172            CALL vor_een( kt, nrvm, un, vn, ua, va )        ! relative vorticity or metric trend 
    165173            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    166174            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
     
    168176            ztrdu(:,:,:) = ua(:,:,:) 
    169177            ztrdv(:,:,:) = va(:,:,:) 
    170             CALL vor_een( kt, ncor, ua, va )                ! planetary vorticity trend 
     178            CALL vor_een( kt, ncor, un, vn, ua, va )        ! planetary vorticity trend 
    171179            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 
    172180            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 
    173181            CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt ) 
    174182         ELSE 
    175             CALL vor_een( kt, ntot, ua, va )                ! total vorticity 
     183                             CALL vor_een( kt, ntot, un, vn, ua, va )    ! total vorticity  
     184            IF( ln_stcor )   CALL vor_ene( kt, ncor, usd, vsd, ua, va )  ! add the Stokes-Coriolis trend 
    176185         ENDIF 
    177186         ! 
     
    189198 
    190199 
    191    SUBROUTINE vor_ene( kt, kvor, pua, pva ) 
     200   SUBROUTINE vor_ene( kt, kvor, pun, pvn, pua, pva ) 
    192201      !!---------------------------------------------------------------------- 
    193202      !!                  ***  ROUTINE vor_ene  *** 
     
    214223      !!---------------------------------------------------------------------- 
    215224      ! 
    216       INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    217       INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
    218       !                                                           ! =nrvm (relative vorticity or metric) 
    219       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua    ! total u-trend 
    220       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva    ! total v-trend 
     225      INTEGER , INTENT(in   )                         ::   kt          ! ocean time-step index  
     226      INTEGER , INTENT(in   )                         ::   kvor        ! =ncor (planetary) ; =ntot (total) ;  
     227      !                                                                ! =nrvm (relative vorticity or metric)  
     228      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua         ! total u-trend  
     229      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva         ! total v-trend  
     230      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun, pvn    ! now velocities 
    221231      ! 
    222232      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     
    250260            DO jj = 1, jpjm1 
    251261               DO ji = 1, fs_jpim1   ! vector opt. 
    252                   zwz(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    253                        &         - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
     262                  zwz(ji,jj) = (   ( pvn(ji+1,jj  ,jk) + pvn(ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       &  
     263                       &         - ( pun(ji  ,jj+1,jk) + pun(ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    254264                       &     * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) 
    255265               END DO 
     
    260270               DO ji = 1, fs_jpim1   ! vector opt. 
    261271                  zwz(ji,jj) = ( ff (ji,jj)                                                                       & 
    262                        &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    263                        &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
     272                       &       + (   ( pvn(ji+1,jj  ,jk) + pvn(ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       &  
     273                       &           - ( pun(ji  ,jj+1,jk) + pun(ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   &  
    264274                       &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                               & 
    265275                       &       ) 
     
    270280         IF( ln_sco ) THEN 
    271281            zwz(:,:) = zwz(:,:) / fse3f(:,:,jk) 
    272             zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk) 
    273             zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
     282            zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * pun(:,:,jk)  
     283            zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * pvn(:,:,jk) 
    274284         ELSE 
    275             zwx(:,:) = e2u(:,:) * un(:,:,jk) 
    276             zwy(:,:) = e1v(:,:) * vn(:,:,jk) 
     285            zwx(:,:) = e2u(:,:) * pun(:,:,jk) 
     286            zwy(:,:) = e1v(:,:) * pvn(:,:,jk) 
    277287         ENDIF 
    278288 
     
    418428 
    419429 
    420    SUBROUTINE vor_ens( kt, kvor, pua, pva ) 
     430   SUBROUTINE vor_ens( kt, kvor, pun, pvn, pua, pva ) 
    421431      !!---------------------------------------------------------------------- 
    422432      !!                ***  ROUTINE vor_ens  *** 
     
    443453      !!---------------------------------------------------------------------- 
    444454      ! 
    445       INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    446       INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
    447          !                                                        ! =nrvm (relative vorticity or metric) 
    448       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua    ! total u-trend 
    449       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva    ! total v-trend 
     455      INTEGER , INTENT(in   )                         ::   kt          ! ocean time-step index  
     456      INTEGER , INTENT(in   )                         ::   kvor        ! =ncor (planetary) ; =ntot (total) ;  
     457         !                                                             ! =nrvm (relative vorticity or metric)  
     458      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua         ! total u-trend  
     459      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva         ! total v-trend  
     460      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun, pvn    ! now velocities 
    450461      ! 
    451462      INTEGER  ::   ji, jj, jk           ! dummy loop indices 
     
    479490            DO jj = 1, jpjm1 
    480491               DO ji = 1, fs_jpim1   ! vector opt. 
    481                   zwz(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    482                        &         - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
     492                  zwz(ji,jj) = (   ( pvn(ji+1,jj  ,jk) + pvn(ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       &  
     493                       &         - ( pun(ji  ,jj+1,jk) + pun(ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    483494                       &     * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) 
    484495               END DO 
     
    489500               DO ji = 1, fs_jpim1   ! vector opt. 
    490501                  zwz(ji,jj) = ( ff (ji,jj)                                                                       & 
    491                        &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    492                        &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
     502                       &       + (   ( pvn(ji+1,jj  ,jk) + pvn(ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       &  
     503                       &           - ( pun(ji  ,jj+1,jk) + pun(ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    493504                       &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                                & 
    494505                       &       ) 
     
    501512               DO ji = 1, jpi                   ! it causes optimization problems on NEC in auto-tasking 
    502513                  zwz(ji,jj) = zwz(ji,jj) / fse3f(ji,jj,jk) 
    503                   zwx(ji,jj) = e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk) 
    504                   zwy(ji,jj) = e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk) 
     514                  zwx(ji,jj) = e2u(ji,jj) * fse3u(ji,jj,jk) * pun(ji,jj,jk) 
     515                  zwy(ji,jj) = e1v(ji,jj) * fse3v(ji,jj,jk) * pvn(ji,jj,jk) 
    505516               END DO 
    506517            END DO 
     
    508519            DO jj = 1, jpj                      ! caution: don't use (:,:) for this loop  
    509520               DO ji = 1, jpi                   ! it causes optimization problems on NEC in auto-tasking 
    510                   zwx(ji,jj) = e2u(ji,jj) * un(ji,jj,jk) 
    511                   zwy(ji,jj) = e1v(ji,jj) * vn(ji,jj,jk) 
     521                  zwx(ji,jj) = e2u(ji,jj) * pun(ji,jj,jk) 
     522                  zwy(ji,jj) = e1v(ji,jj) * pvn(ji,jj,jk) 
    512523               END DO 
    513524            END DO 
     
    536547 
    537548 
    538    SUBROUTINE vor_een( kt, kvor, pua, pva ) 
     549   SUBROUTINE vor_een( kt, kvor, pun, pvn, pua, pva ) 
    539550      !!---------------------------------------------------------------------- 
    540551      !!                ***  ROUTINE vor_een  *** 
     
    554565      !!---------------------------------------------------------------------- 
    555566      ! 
    556       INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    557       INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
    558       !                                                           ! =nrvm (relative vorticity or metric) 
    559       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua    ! total u-trend 
    560       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva    ! total v-trend 
     567      INTEGER , INTENT(in   )                         ::   kt          ! ocean time-step index  
     568      INTEGER , INTENT(in   )                         ::   kvor        ! =ncor (planetary) ; =ntot (total) ;  
     569      !                                                                ! =nrvm (relative vorticity or metric)  
     570      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua         ! total u-trend  
     571      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva         ! total v-trend  
     572      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun, pvn    ! now velocities 
    561573      !! 
    562574      INTEGER  ::   ji, jj, jk                                    ! dummy loop indices 
     
    604616                     ze3  = ( fse3t(ji,jj+1,jk)*tmask(ji,jj+1,jk) + fse3t(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk)   & 
    605617                        &   + fse3t(ji,jj  ,jk)*tmask(ji,jj  ,jk) + fse3t(ji+1,jj  ,jk)*tmask(ji+1,jj  ,jk) ) 
    606                      IF( ze3 /= 0._wp )   ze3f(ji,jj,jk) = 4.0_wp / ze3 
     618                     IF( ze3 /= 0._wp ) THEN  ;   ze3f(ji,jj,jk) = 4.0_wp / ze3  
     619                     ELSE                     ;   ze3f(ji,jj,jk) = 0._wp  
     620                     ENDIF 
    607621                  END DO 
    608622               END DO 
     
    616630                     zmsk = (                   tmask(ji,jj+1,jk) +                     tmask(ji+1,jj+1,jk)   & 
    617631                        &                     + tmask(ji,jj  ,jk) +                     tmask(ji+1,jj  ,jk) ) 
    618                      IF( ze3 /= 0._wp )   ze3f(ji,jj,jk) = zmsk / ze3 
     632                     IF( ze3 /= 0._wp ) THEN  ;  ze3f(ji,jj,jk) = zmsk / ze3  
     633                     ELSE                     ;  ze3f(ji,jj,jk) = 0._wp  
     634                     ENDIF 
    619635                  END DO 
    620636               END DO 
     
    643659            DO jj = 1, jpjm1 
    644660               DO ji = 1, fs_jpim1   ! vector opt. 
    645                   zwz(ji,jj) = (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    646                        &         - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
     661                  zwz(ji,jj) = (   ( pvn(ji+1,jj  ,jk) + pvn(ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       &  
     662                       &         - ( pun(ji  ,jj+1,jk) + pun(ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    647663                       &     * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) ) * ze3f(ji,jj,jk) 
    648664               END DO 
     
    655671               DO ji = 1, fs_jpim1   ! vector opt. 
    656672                  zwz(ji,jj) = ( ff (ji,jj)                                                                       & 
    657                        &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    658                        &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
     673                       &       + (   ( pvn(ji+1,jj  ,jk) + pvn(ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       &  
     674                       &           - ( pun(ji  ,jj+1,jk) + pun(ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    659675                       &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                                & 
    660676                       &       ) * ze3f(ji,jj,jk) 
     
    664680         END SELECT 
    665681 
    666          zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk) 
    667          zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
     682         zwx(:,:) = e2u(:,:) * fse3u(:,:,jk) * pun(:,:,jk) 
     683         zwy(:,:) = e1v(:,:) * fse3v(:,:,jk) * pvn(:,:,jk) 
    668684 
    669685         ! Compute and add the vorticity term trend 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/cpl_oasis3.F90

    r10392 r10473  
    8181   INTEGER                    ::   nsnd         ! total number of fields sent  
    8282   INTEGER                    ::   ncplmodel    ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
    83    INTEGER, PUBLIC, PARAMETER ::   nmaxfld=50   ! Maximum number of coupling fields 
     83   INTEGER, PUBLIC, PARAMETER ::   nmaxfld=55   ! Maximum number of coupling fields 
    8484   INTEGER, PUBLIC, PARAMETER ::   nmaxcat=5    ! Maximum number of coupling fields 
    8585   INTEGER, PUBLIC, PARAMETER ::   nmaxcpl=5    ! Maximum number of coupling fields 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90

    r10392 r10473  
    6565   LOGICAL , PUBLIC ::   ln_cdgw        !: true if neutral drag coefficient from wave model 
    6666   LOGICAL , PUBLIC ::   ln_sdw         !: true if 3d stokes drift from wave model 
     67   LOGICAL , PUBLIC ::   ln_tauoc       !: true if normalized stress from wave is used   
     68   LOGICAL , PUBLIC ::   ln_stcor       !: true if Stokes-Coriolis term is used 
    6769   ! 
    6870   LOGICAL , PUBLIC ::   ln_icebergs    !: Icebergs 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r8058 r10473  
    737737 
    738738      !! Neutral coefficients at 10m: 
    739       IF( ln_cdgw ) THEN      ! wave drag case 
     739      IF( ln_wave .AND. ln_cdgw ) THEN      ! wave drag case 
    740740         cdn_wave(:,:) = cdn_wave(:,:) + rsmall * ( 1._wp - tmask(:,:,1) ) 
    741741         ztmp0   (:,:) = cdn_wave(:,:) 
     
    783783         END IF 
    784784        
    785          IF( ln_cdgw ) THEN      ! surface wave case 
     785         IF( ln_wave .AND. ln_cdgw ) THEN      ! surface wave case 
    786786            sqrt_Cd = vkarmn / ( vkarmn / sqrt_Cd_n10 - zpsi_m_u )  
    787787            Cd      = sqrt_Cd * sqrt_Cd 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_mfs.F90

    r8058 r10473  
    2424   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2525   USE prtctl          ! Print control 
    26    USE sbcwave,ONLY : cdn_wave !wave module 
    2726 
    2827   IMPLICIT NONE 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r10396 r10473  
    2323   USE sbcapr 
    2424   USE sbcdcy          ! surface boundary condition: diurnal cycle 
     25   USE sbcwave         ! surface boundary condition: waves 
    2526   USE phycst          ! physical constants 
    2627#if defined key_lim3 
     
    108109   INTEGER, PARAMETER ::   jpr_e3t1st = 41            ! first T level thickness  
    109110   INTEGER, PARAMETER ::   jpr_fraqsr = 42            ! fraction of solar net radiation absorbed in the first ocean level 
    110    INTEGER, PARAMETER ::   jprcv      = 42            ! total number of fields received 
     111   INTEGER, PARAMETER ::   jpr_mslp   = 43            ! mean sea level pressure    
     112   INTEGER, PARAMETER ::   jpr_hsig   = 44            ! Hsig    
     113   INTEGER, PARAMETER ::   jpr_phioc  = 45            ! Wave=>ocean energy flux    
     114   INTEGER, PARAMETER ::   jpr_sdrftx = 46            ! Stokes drift on grid 1    
     115   INTEGER, PARAMETER ::   jpr_sdrfty = 47            ! Stokes drift on grid 2    
     116   INTEGER, PARAMETER ::   jpr_wper   = 48            ! Mean wave period   
     117   INTEGER, PARAMETER ::   jpr_wnum   = 49            ! Mean wavenumber   
     118   INTEGER, PARAMETER ::   jpr_wstrf  = 50            ! Stress fraction adsorbed by waves   
     119   INTEGER, PARAMETER ::   jpr_wdrag  = 51            ! Neutral surface drag coefficient   
     120   INTEGER, PARAMETER ::   jprcv      = 51            ! total number of fields received 
    111121 
    112122   INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere 
     
    138148   INTEGER, PARAMETER ::   jps_e3t1st = 27            ! first level depth (vvl) 
    139149   INTEGER, PARAMETER ::   jps_fraqsr = 28            ! fraction of solar net radiation absorbed in the first ocean level 
    140    INTEGER, PARAMETER ::   jpsnd      = 28            ! total number of fields sended 
     150   INTEGER, PARAMETER ::   jps_ficet  = 29            ! total ice fraction     
     151   INTEGER, PARAMETER ::   jps_ocxw   = 30            ! currents on grid 1     
     152   INTEGER, PARAMETER ::   jps_ocyw   = 31            ! currents on grid 2   
     153   INTEGER, PARAMETER ::   jps_wlev   = 32            ! water level    
     154   INTEGER, PARAMETER ::   jpsnd      = 32            ! total number of fields sent 
    141155 
    142156   !                                                         !!** namelist namsbc_cpl ** 
     
    152166   ! Received from the atmosphere                     ! 
    153167   TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf 
    154    TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2                         
     168   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp                              
     169   ! Send to waves    
     170   TYPE(FLD_C) ::   sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev    
     171   ! Received from waves    
     172   TYPE(FLD_C) ::   sn_rcv_hsig,sn_rcv_phioc,sn_rcv_sdrfx,sn_rcv_sdrfy,sn_rcv_wper,sn_rcv_wnum,sn_rcv_wstrf,sn_rcv_wdrag 
    155173   ! Other namelist parameters                        ! 
    156174   INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
     
    164182 
    165183   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   albedo_oce_mix     ! ocean albedo sent to atmosphere (mix clear/overcast sky) 
     184     
     185   REAL(wp) ::   rpref = 101000._wp   ! reference atmospheric pressure[N/m2]    
     186   REAL(wp) ::   r1_grau              ! = 1.e0 / (grav * rau0 
    166187 
    167188   INTEGER , ALLOCATABLE, SAVE, DIMENSION(    :) ::   nrcvinfo           ! OASIS info argument 
     
    182203      !!             ***  FUNCTION sbc_cpl_alloc  *** 
    183204      !!---------------------------------------------------------------------- 
    184       INTEGER :: ierr(3) 
     205      INTEGER :: ierr(4) 
    185206      !!---------------------------------------------------------------------- 
    186207      ierr(:) = 0 
     
    195216      ALLOCATE( xcplmask(jpi,jpj,0:3) , STAT=ierr(3) ) 
    196217      ! 
     218      IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(4) ) 
     219 
    197220      sbc_cpl_alloc = MAXVAL( ierr ) 
    198221      IF( lk_mpp            )   CALL mpp_sum ( sbc_cpl_alloc ) 
     
    221244      REAL(wp), POINTER, DIMENSION(:,:) ::   zacs, zaos 
    222245      !! 
    223       NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2,      & 
    224          &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr,      & 
    225          &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal   , sn_rcv_iceflx,   & 
    226          &                  sn_rcv_co2 , nn_cplmodel  , ln_usecplmask 
     246      NAMELIST/namsbc_cpl/  sn_snd_temp , sn_snd_alb  , sn_snd_thick , sn_snd_crt   , sn_snd_co2,      &    
     247         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau   , sn_rcv_dqnsdt, sn_rcv_qsr,      &    
     248         &                  sn_snd_ifrac, sn_snd_crtw , sn_snd_wlev  , sn_rcv_hsig  , sn_rcv_phioc ,   &    
     249         &                  sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper  , sn_rcv_wnum  , sn_rcv_wstrf ,   &   
     250         &                  sn_rcv_wdrag, sn_rcv_qns  , sn_rcv_emp   , sn_rcv_rnf   , sn_rcv_cal   ,   &   
     251         &                  sn_rcv_iceflx,sn_rcv_co2  , nn_cplmodel  , ln_usecplmask, sn_rcv_mslp 
    227252      !!--------------------------------------------------------------------- 
    228253      ! 
     
    265290         WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 
    266291         WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')' 
     292         WRITE(numout,*)'      significant wave heigth         = ', TRIM(sn_rcv_hsig%cldes  ), ' (', TRIM(sn_rcv_hsig%clcat  ), ')'    
     293         WRITE(numout,*)'      wave to oce energy flux         = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')'    
     294         WRITE(numout,*)'      Surface Stokes drift grid u     = ', TRIM(sn_rcv_sdrfx%cldes ), ' (', TRIM(sn_rcv_sdrfx%clcat ), ')'    
     295         WRITE(numout,*)'      Surface Stokes drift grid v     = ', TRIM(sn_rcv_sdrfy%cldes ), ' (', TRIM(sn_rcv_sdrfy%clcat ), ')'    
     296         WRITE(numout,*)'      Mean wave period                = ', TRIM(sn_rcv_wper%cldes  ), ' (', TRIM(sn_rcv_wper%clcat  ), ')'    
     297         WRITE(numout,*)'      Mean wave number                = ', TRIM(sn_rcv_wnum%cldes  ), ' (', TRIM(sn_rcv_wnum%clcat  ), ')'    
     298         WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')'    
     299         WRITE(numout,*)'      Neutral surf drag coefficient   = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' 
    267300         WRITE(numout,*)'  sent fields (multiple ice categories)' 
    268301         WRITE(numout,*)'      surface temperature             = ', TRIM(sn_snd_temp%cldes  ), ' (', TRIM(sn_snd_temp%clcat  ), ')' 
    269302         WRITE(numout,*)'      albedo                          = ', TRIM(sn_snd_alb%cldes   ), ' (', TRIM(sn_snd_alb%clcat   ), ')' 
    270303         WRITE(numout,*)'      ice/snow thickness              = ', TRIM(sn_snd_thick%cldes ), ' (', TRIM(sn_snd_thick%clcat ), ')' 
     304         WRITE(numout,*)'      total ice fraction              = ', TRIM(sn_snd_ifrac%cldes ), ' (', TRIM(sn_snd_ifrac%clcat ), ')' 
    271305         WRITE(numout,*)'      surface current                 = ', TRIM(sn_snd_crt%cldes   ), ' (', TRIM(sn_snd_crt%clcat   ), ')' 
    272306         WRITE(numout,*)'                      - referential   = ', sn_snd_crt%clvref  
     
    274308         WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd 
    275309         WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')' 
     310         WRITE(numout,*)'      water level                     = ', TRIM(sn_snd_wlev%cldes  ), ' (', TRIM(sn_snd_wlev%clcat  ), ')'    
     311         WRITE(numout,*)'      mean sea level pressure         = ', TRIM(sn_rcv_mslp%cldes  ), ' (', TRIM(sn_rcv_mslp%clcat  ), ')'    
     312         WRITE(numout,*)'      surface current to waves        = ', TRIM(sn_snd_crtw%cldes  ), ' (', TRIM(sn_snd_crtw%clcat  ), ')'    
     313         WRITE(numout,*)'                      - referential   = ', sn_snd_crtw%clvref    
     314         WRITE(numout,*)'                      - orientation   = ', sn_snd_crtw%clvor    
     315         WRITE(numout,*)'                      - mesh          = ', sn_snd_crtw%clvgrd 
    276316         WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel 
    277317         WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask 
     
    312352      !  
    313353      ! Vectors: change of sign at north fold ONLY if on the local grid 
     354      IF( TRIM( sn_rcv_tau%cldes ) == 'oce only' .OR. TRIM(sn_rcv_tau%cldes ) == 'oce and ice') THEN ! avoid working with the atmospheric fields if they are not coupled 
    314355      IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' )   srcv(jpr_otx1:jpr_itz2)%nsgn = -1. 
    315356       
     
    383424         srcv(jpr_ity1)%clgrid = 'V'                  ! i.e. it is always at U- & V-points for i- & j-comp. resp. 
    384425      ENDIF 
     426      ENDIF 
    385427        
    386428      !                                                      ! ------------------------- ! 
     
    479521      !                                                      ! ------------------------- ! 
    480522      srcv(jpr_co2 )%clname = 'O_AtmCO2'   ;   IF( TRIM(sn_rcv_co2%cldes   ) == 'coupled' )    srcv(jpr_co2 )%laction = .TRUE. 
     523 
     524      !                                                      ! ------------------------- !    
     525      !                                                      ! Mean Sea Level Pressure   !    
     526      !                                                      ! ------------------------- !    
     527      srcv(jpr_mslp)%clname = 'O_MSLP'     ;   IF( TRIM(sn_rcv_mslp%cldes  ) == 'coupled' )    srcv(jpr_mslp)%laction = .TRUE.    
     528 
    481529      !                                                      ! ------------------------- ! 
    482530      !                                                      !   topmelt and botmelt     !    
     
    492540         srcv(jpr_topm:jpr_botm)%laction = .TRUE. 
    493541      ENDIF 
     542      !                                                      ! ------------------------- !   
     543      !                                                      !      Wave breaking        !       
     544      !                                                      ! ------------------------- !    
     545      srcv(jpr_hsig)%clname  = 'O_Hsigwa'    ! significant wave height   
     546      IF( TRIM(sn_rcv_hsig%cldes  ) == 'coupled' )  THEN   
     547         srcv(jpr_hsig)%laction = .TRUE.   
     548         cpl_hsig = .TRUE.   
     549      ENDIF   
     550      srcv(jpr_phioc)%clname = 'O_PhiOce'    ! wave to ocean energy   
     551      IF( TRIM(sn_rcv_phioc%cldes ) == 'coupled' )  THEN   
     552         srcv(jpr_phioc)%laction = .TRUE.   
     553         cpl_phioc = .TRUE.   
     554      ENDIF   
     555      srcv(jpr_sdrftx)%clname = 'O_Sdrfx'    ! Stokes drift in the u direction   
     556      IF( TRIM(sn_rcv_sdrfx%cldes ) == 'coupled' )  THEN   
     557         srcv(jpr_sdrftx)%laction = .TRUE.   
     558         cpl_sdrftx = .TRUE.   
     559      ENDIF   
     560      srcv(jpr_sdrfty)%clname = 'O_Sdrfy'    ! Stokes drift in the v direction   
     561      IF( TRIM(sn_rcv_sdrfy%cldes ) == 'coupled' )  THEN   
     562         srcv(jpr_sdrfty)%laction = .TRUE.   
     563         cpl_sdrfty = .TRUE.   
     564      ENDIF   
     565      srcv(jpr_wper)%clname = 'O_WPer'       ! mean wave period   
     566      IF( TRIM(sn_rcv_wper%cldes  ) == 'coupled' )  THEN   
     567         srcv(jpr_wper)%laction = .TRUE.   
     568         cpl_wper = .TRUE.   
     569      ENDIF   
     570      srcv(jpr_wnum)%clname = 'O_WNum'       ! mean wave number   
     571      IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' )  THEN   
     572         srcv(jpr_wnum)%laction = .TRUE.   
     573         cpl_wnum = .TRUE.   
     574      ENDIF   
     575      srcv(jpr_wstrf)%clname = 'O_WStrf'     ! stress fraction adsorbed by the wave   
     576      IF( TRIM(sn_rcv_wstrf%cldes ) == 'coupled' )  THEN   
     577         srcv(jpr_wstrf)%laction = .TRUE.   
     578         cpl_wstrf = .TRUE.   
     579      ENDIF   
     580      srcv(jpr_wdrag)%clname = 'O_WDrag'     ! neutral surface drag coefficient   
     581      IF( TRIM(sn_rcv_wdrag%cldes ) == 'coupled' )  THEN   
     582         srcv(jpr_wdrag)%laction = .TRUE.   
     583         cpl_wdrag = .TRUE.   
     584      ENDIF   
     585      !    
    494586      !                                                      ! ------------------------------- ! 
    495587      !                                                      !   OPA-SAS coupling - rcv by opa !    
     
    646738      !                                                      ! ------------------------- ! 
    647739      ssnd(jps_fice)%clname = 'OIceFrc' 
     740      ssnd(jps_ficet)%clname = 'OIceFrcT' 
    648741      ssnd(jps_hice)%clname = 'OIceTck' 
    649742      ssnd(jps_hsnw)%clname = 'OSnwTck' 
     
    654747      ENDIF 
    655748       
     749      IF (TRIM( sn_snd_ifrac%cldes )  == 'coupled') ssnd(jps_ficet)%laction = .TRUE. 
     750 
    656751      SELECT CASE ( TRIM( sn_snd_thick%cldes ) ) 
    657752      CASE( 'none'         )       ! nothing to do 
     
    674769      ssnd(jps_ocy1)%clname = 'O_OCury1'   ;   ssnd(jps_ivy1)%clname = 'O_IVely1' 
    675770      ssnd(jps_ocz1)%clname = 'O_OCurz1'   ;   ssnd(jps_ivz1)%clname = 'O_IVelz1' 
     771      ssnd(jps_ocxw)%clname = 'O_OCurxw'    
     772      ssnd(jps_ocyw)%clname = 'O_OCuryw' 
    676773      ! 
    677774      ssnd(jps_ocx1:jps_ivz1)%nsgn = -1.   ! vectors: change of the sign at the north fold 
     
    694791      END SELECT 
    695792 
     793      ssnd(jps_ocxw:jps_ocyw)%nsgn = -1.   ! vectors: change of the sign at the north fold    
     794               
     795      IF( sn_snd_crtw%clvgrd == 'U,V' ) THEN    
     796         ssnd(jps_ocxw)%clgrid = 'U' ; ssnd(jps_ocyw)%clgrid = 'V'    
     797      ELSE IF( sn_snd_crtw%clvgrd /= 'T' ) THEN    
     798         CALL ctl_stop( 'sn_snd_crtw%clvgrd must be equal to T' )    
     799      ENDIF    
     800      IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) ssnd(jps_ocxw:jps_ocyw)%nsgn = 1.    
     801      SELECT CASE( TRIM( sn_snd_crtw%cldes ) )    
     802         CASE( 'none'                 )   ; ssnd(jps_ocxw:jps_ocyw)%laction = .FALSE.    
     803         CASE( 'oce only'             )   ; ssnd(jps_ocxw:jps_ocyw)%laction = .TRUE.    
     804         CASE( 'weighted oce and ice' )   !   nothing to do    
     805         CASE( 'mixed oce-ice'        )   ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.    
     806         CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crtw%cldes' )    
     807      END SELECT 
     808 
    696809      !                                                      ! ------------------------- ! 
    697810      !                                                      !          CO2 flux         ! 
    698811      !                                                      ! ------------------------- ! 
    699812      ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE. 
     813 
     814      !                                                      ! ------------------------- !    
     815      !                                                      !     Sea surface height    !    
     816      !                                                      ! ------------------------- !    
     817      ssnd(jps_wlev)%clname = 'O_Wlevel' ;  IF( TRIM(sn_snd_wlev%cldes) == 'coupled' )   ssnd(jps_wlev)%laction = .TRUE. 
    700818 
    701819      !                                                      ! ------------------------------- ! 
     
    792910      IF( ln_dm2dc .AND. ln_cpl .AND. ncpl_qsr_freq /= 86400 )   & 
    793911         &   CALL ctl_stop( 'sbc_cpl_init: diurnal cycle reconstruction (ln_dm2dc) needs daily couping for solar radiation' ) 
    794       ncpl_qsr_freq = 86400 / ncpl_qsr_freq 
     912      IF( ln_dm2dc .AND. ln_cpl ) ncpl_qsr_freq = 86400 / ncpl_qsr_freq 
    795913 
    796914      CALL wrk_dealloc( jpi,jpj, zacs, zaos ) 
     
    846964      !!                        emp          upward mass flux [evap. - precip. (- runoffs) (- calving)] (ocean only case) 
    847965      !!---------------------------------------------------------------------- 
     966      USE zdf_oce,  ONLY : ln_zdfqiao 
     967 
    848968      INTEGER, INTENT(in)           ::   kt          ! ocean model time step index 
    849969      INTEGER, INTENT(in)           ::   k_fsbc      ! frequency of sbc (-> ice model) computation  
     
    10281148      IF( srcv(jpr_co2)%laction )   atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1) 
    10291149#endif 
     1150      !    
     1151      !                                                      ! ========================= !    
     1152      !                                                      ! Mean Sea Level Pressure   !   (taum)    
     1153      !                                                      ! ========================= !    
     1154      !    
     1155      IF( srcv(jpr_mslp)%laction ) THEN                    ! UKMO SHELF effect of atmospheric pressure on SSH    
     1156          IF( kt /= nit000 )   ssh_ibb(:,:) = ssh_ib(:,:)    !* Swap of ssh_ib fields    
     1157        
     1158          r1_grau = 1.e0 / (grav * rau0)               !* constant for optimization    
     1159          ssh_ib(:,:) = - ( frcv(jpr_mslp)%z3(:,:,1) - rpref ) * r1_grau    ! equivalent ssh (inverse barometer)    
     1160          apr   (:,:) =     frcv(jpr_mslp)%z3(:,:,1) !atmospheric pressure    
     1161        
     1162          IF( kt == nit000 ) ssh_ibb(:,:) = ssh_ib(:,:)  ! correct this later (read from restart if possible)    
     1163      END IF    
     1164      !   
     1165      IF( ln_sdw ) THEN  ! Stokes Drift correction activated   
     1166      !                                                      ! ========================= !    
     1167      !                                                      !       Stokes drift u      !   
     1168      !                                                      ! ========================= !    
     1169         IF( srcv(jpr_sdrftx)%laction ) ut0sd(:,:) = frcv(jpr_sdrftx)%z3(:,:,1)   
     1170      !   
     1171      !                                                      ! ========================= !    
     1172      !                                                      !       Stokes drift v      !   
     1173      !                                                      ! ========================= !    
     1174         IF( srcv(jpr_sdrfty)%laction ) vt0sd(:,:) = frcv(jpr_sdrfty)%z3(:,:,1)   
     1175      !   
     1176      !                                                      ! ========================= !    
     1177      !                                                      !      Wave mean period     !   
     1178      !                                                      ! ========================= !    
     1179         IF( srcv(jpr_wper)%laction ) wmp(:,:) = frcv(jpr_wper)%z3(:,:,1)   
     1180      !   
     1181      !                                                      ! ========================= !    
     1182      !                                                      !  Significant wave height  !   
     1183      !                                                      ! ========================= !    
     1184         IF( srcv(jpr_hsig)%laction ) hsw(:,:) = frcv(jpr_hsig)%z3(:,:,1)   
     1185      !   
     1186      !                                                      ! ========================= !    
     1187      !                                                      !    Vertical mixing Qiao   !   
     1188      !                                                      ! ========================= !    
     1189         IF( srcv(jpr_wnum)%laction .AND. ln_zdfqiao ) wnum(:,:) = frcv(jpr_wnum)%z3(:,:,1)   
     1190        
     1191         ! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode   
     1192         IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. srcv(jpr_wper)%laction &   
     1193                                                                    .OR. srcv(jpr_hsig)%laction ) &   
     1194            CALL sbc_stokes()   
     1195      ENDIF   
     1196      !                                                      ! ========================= !    
     1197      !                                                      ! Stress adsorbed by waves  !   
     1198      !                                                      ! ========================= !    
     1199      IF( srcv(jpr_wstrf)%laction .AND. ln_tauoc ) tauoc_wave(:,:) = frcv(jpr_wstrf)%z3(:,:,1)   
     1200        
     1201      !                                                      ! ========================= !    
     1202      !                                                      ! Wave drag coefficient   !   
     1203      !                                                      ! ========================= !    
     1204      IF( srcv(jpr_wdrag)%laction .AND. ln_cdgw ) cdn_wave(:,:) = frcv(jpr_wdrag)%z3(:,:,1) 
    10301205 
    10311206      !  Fields received by SAS when OASIS coupling 
     
    21012276      ENDIF 
    21022277      ! 
     2278      !                                                      ! ------------------------- !    
     2279      !                                                      !  Surface current to waves !    
     2280      !                                                      ! ------------------------- !    
     2281      IF( ssnd(jps_ocxw)%laction .OR. ssnd(jps_ocyw)%laction ) THEN    
     2282          !        
     2283          !                                                  j+1  j     -----V---F    
     2284          ! surface velocity always sent from T point                    !       |    
     2285          !                                                       j      |   T   U    
     2286          !                                                              |       |    
     2287          !                                                   j   j-1   -I-------|    
     2288          !                                               (for I)        |       |    
     2289          !                                                             i-1  i   i    
     2290          !                                                              i      i+1 (for I)    
     2291          SELECT CASE( TRIM( sn_snd_crtw%cldes ) )    
     2292          CASE( 'oce only'             )      ! C-grid ==> T    
     2293             DO jj = 2, jpjm1    
     2294                DO ji = fs_2, fs_jpim1   ! vector opt.    
     2295                   zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj  ,1) )    
     2296                   zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji , jj-1,1) )     
     2297                END DO    
     2298             END DO    
     2299          CASE( 'weighted oce and ice' )       
     2300             SELECT CASE ( cp_ice_msh )    
     2301             CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T    
     2302                DO jj = 2, jpjm1    
     2303                   DO ji = fs_2, fs_jpim1   ! vector opt.    
     2304                      zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un (ji-1,jj  ,1) ) * zfr_l(ji,jj)      
     2305                      zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn (ji  ,jj-1,1) ) * zfr_l(ji,jj)    
     2306                      zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)    
     2307                      zity1(ji,jj) = 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)    
     2308                   END DO    
     2309                END DO    
     2310             CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T    
     2311                DO jj = 2, jpjm1    
     2312                   DO ji = 2, jpim1   ! NO vector opt.    
     2313                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)      
     2314                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)      
     2315                      zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &    
     2316                         &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2317                      zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &    
     2318                         &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2319                   END DO    
     2320                END DO    
     2321             CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T    
     2322                DO jj = 2, jpjm1    
     2323                   DO ji = 2, jpim1   ! NO vector opt.    
     2324                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)      
     2325                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)      
     2326                      zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &    
     2327                         &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2328                      zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &    
     2329                         &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2330                   END DO    
     2331                END DO    
     2332             END SELECT    
     2333             CALL lbc_lnk( zitx1, 'T', -1. )   ;   CALL lbc_lnk( zity1, 'T', -1. )    
     2334          CASE( 'mixed oce-ice'        )    
     2335             SELECT CASE ( cp_ice_msh )    
     2336             CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T    
     2337                DO jj = 2, jpjm1    
     2338                   DO ji = fs_2, fs_jpim1   ! vector opt.    
     2339                      zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &    
     2340                         &         + 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)    
     2341                      zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &    
     2342                         &         + 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)    
     2343                   END DO    
     2344                END DO    
     2345             CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T    
     2346                DO jj = 2, jpjm1    
     2347                   DO ji = 2, jpim1   ! NO vector opt.    
     2348                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &       
     2349                         &         + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &    
     2350                         &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2351                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &     
     2352                         &         + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &    
     2353                         &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2354                   END DO    
     2355                END DO    
     2356             CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T    
     2357                DO jj = 2, jpjm1    
     2358                   DO ji = 2, jpim1   ! NO vector opt.    
     2359                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &       
     2360                         &         + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &    
     2361                         &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2362                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &     
     2363                         &         + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &    
     2364                         &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2365                   END DO    
     2366                END DO    
     2367             END SELECT    
     2368          END SELECT    
     2369         CALL lbc_lnk( zotx1, ssnd(jps_ocxw)%clgrid, -1. )   ; CALL lbc_lnk( zoty1, ssnd(jps_ocyw)%clgrid, -1. )    
     2370         !    
     2371         !    
     2372         IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) THEN             ! Rotation of the components    
     2373         !                                                                        ! Ocean component    
     2374            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->e', ztmp1 )       ! 1st component     
     2375            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->n', ztmp2 )       ! 2nd component     
     2376            zotx1(:,:) = ztmp1(:,:)                                                   ! overwrite the components     
     2377            zoty1(:,:) = ztmp2(:,:)     
     2378            IF( ssnd(jps_ivx1)%laction ) THEN                                     ! Ice component    
     2379               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 )    ! 1st component     
     2380               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 )    ! 2nd component     
     2381               zitx1(:,:) = ztmp1(:,:)                                                ! overwrite the components     
     2382               zity1(:,:) = ztmp2(:,:)    
     2383            ENDIF    
     2384         ENDIF    
     2385         !    
     2386!         ! spherical coordinates to cartesian -> 2 components to 3 components    
     2387!         IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN    
     2388!            ztmp1(:,:) = zotx1(:,:)                     ! ocean currents    
     2389!            ztmp2(:,:) = zoty1(:,:)    
     2390!            CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 )    
     2391!            !    
     2392!            IF( ssnd(jps_ivx1)%laction ) THEN           ! ice velocities    
     2393!               ztmp1(:,:) = zitx1(:,:)    
     2394!               ztmp1(:,:) = zity1(:,:)    
     2395!               CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 )    
     2396!            ENDIF    
     2397!         ENDIF    
     2398         !    
     2399         IF( ssnd(jps_ocxw)%laction )   CALL cpl_snd( jps_ocxw, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info )   ! ocean x current 1st grid    
     2400         IF( ssnd(jps_ocyw)%laction )   CALL cpl_snd( jps_ocyw, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info )   ! ocean y current 1st grid    
     2401         !     
     2402      ENDIF    
     2403      !    
     2404      IF( ssnd(jps_ficet)%laction ) THEN    
     2405         CALL cpl_snd( jps_ficet, isec, RESHAPE ( fr_i, (/jpi,jpj,1/) ), info )    
     2406      END IF    
     2407      !                                                      ! ------------------------- !    
     2408      !                                                      !   Water levels to waves   !    
     2409      !                                                      ! ------------------------- !    
     2410      IF( ssnd(jps_wlev)%laction ) THEN    
     2411         IF( ln_apr_dyn ) THEN     
     2412            IF( kt /= nit000 ) THEN     
     2413               ztmp1(:,:) = sshb(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) )     
     2414            ELSE     
     2415               ztmp1(:,:) = sshb(:,:)     
     2416            ENDIF     
     2417         ELSE     
     2418            ztmp1(:,:) = sshn(:,:)     
     2419         ENDIF     
     2420         CALL cpl_snd( jps_wlev  , isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )    
     2421      END IF 
    21032422      ! 
    21042423      !  Fields sent by OPA to SAS when doing OPA<->SAS coupling 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r10394 r10473  
    8989         &             ln_blk_mfs, ln_apr_dyn, nn_ice, nn_ice_embd, ln_dm2dc   , ln_rnf   ,   & 
    9090         &             ln_ssr    , nn_isf    , nn_fwb, ln_cdgw    , ln_wave    , ln_sdw   ,   & 
    91          &             nn_lsm    , nn_limflx , nn_components, ln_cpl 
     91         &             ln_tauoc  , ln_stcor  , nn_lsm, nn_limflx , nn_components, ln_cpl 
    9292      INTEGER  ::   ios 
    9393      INTEGER  ::   ierr, ierr0, ierr1, ierr2, ierr3, jpm 
     
    132132         WRITE(numout,*) '              ocean-atmosphere coupled formulation       ln_cpl      = ', ln_cpl 
    133133         WRITE(numout,*) '              forced-coupled mixed formulation           ln_mixcpl   = ', ln_mixcpl 
     134         WRITE(numout,*) '              wave physics                               ln_wave     = ', ln_wave   
     135         WRITE(numout,*) '                 Stokes drift corr. to vert. velocity    ln_sdw      = ', ln_sdw   
     136         WRITE(numout,*) '                 wave modified ocean stress              ln_tauoc    = ', ln_tauoc   
     137         WRITE(numout,*) '                 Stokes coriolis term                    ln_stcor    = ', ln_stcor   
     138         WRITE(numout,*) '                 neutral drag coefficient (CORE, MFS)    ln_cdgw     = ', ln_cdgw 
    134139         WRITE(numout,*) '              OASIS coupling (with atm or sas)           lk_oasis    = ', lk_oasis 
    135140         WRITE(numout,*) '              components of your executable              nn_components = ', nn_components 
     
    216221       
    217222      IF ( ln_wave ) THEN 
    218       !Activated wave module but neither drag nor stokes drift activated 
    219          IF ( .NOT.(ln_cdgw .OR. ln_sdw) )   THEN 
    220             CALL ctl_warn( 'Ask for wave coupling but nor drag coefficient (ln_cdgw=F) neither stokes drift activated (ln_sdw=F)' ) 
    221       !drag coefficient read from wave model definable only with mfs bulk formulae and core  
     223         !Activated wave module but neither drag nor stokes drift activated   
     224         IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor ) )   THEN    
     225             CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_stcor=F')   
     226         !drag coefficient read from wave model definable only with mfs bulk formulae and core 
    222227         ELSEIF (ln_cdgw .AND. .NOT.(ln_blk_mfs .OR. ln_blk_core) )       THEN        
    223228             CALL ctl_stop( 'drag coefficient read from wave model definable only with mfs bulk formulae and core') 
     229         ELSEIF (ln_stcor .AND. .NOT. ln_sdw) THEN    
     230             CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T') 
    224231         ENDIF 
    225232      ELSE 
    226       IF ( ln_cdgw .OR. ln_sdw  )                                         &  
    227          &   CALL ctl_stop('Not Activated Wave Module (ln_wave=F) but     & 
    228          & asked coupling with drag coefficient (ln_cdgw =T) or Stokes drift (ln_sdw=T) ') 
     233         IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor ) &    
     234            &   CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ',    &   
     235            &                  'with drag coefficient (ln_cdgw =T) '  ,                        &   
     236            &                  'or Stokes Drift (ln_sdw=T) ' ,                                 &   
     237            &                  'or ocean stress modification due to waves (ln_tauoc=T) ',      &   
     238            &                  'or Stokes-Coriolis term (ln_stcori=T)'  ) 
    229239      ENDIF  
    230240      !                          ! Choice of the Surface Boudary Condition (set nsbc) 
     
    307317 
    308318      IF( nn_ice == 4      )   CALL cice_sbc_init( nsbc )      ! CICE initialisation 
     319      !   
     320      IF( ln_wave          )   CALL sbc_wave_init              ! surface wave initialisation   
     321      ! 
    309322       
    310323   END SUBROUTINE sbc_init 
     
    386399      IF( ln_mixcpl )        CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! forced-coupled mixed formulation after forcing 
    387400 
     401      IF ( ln_wave .AND. ln_tauoc) THEN                 ! Wave stress subctracted   
     402            utau(:,:) = utau(:,:)*tauoc_wave(:,:)   
     403            vtau(:,:) = vtau(:,:)*tauoc_wave(:,:)   
     404            taum(:,:) = taum(:,:)*tauoc_wave(:,:)   
     405      !   
     406            SELECT CASE( nsbc )   
     407            CASE(  0,1,2,3,5,-1 )  ;   
     408                IF(lwp .AND. kt == nit000 ) WRITE(numout,*) 'WARNING: You are subtracting the wave stress to the ocean. &   
     409                        & If not requested select ln_tauoc=.false'   
     410            END SELECT   
     411      !   
     412      END IF 
    388413 
    389414      !                                            !==  Misc. Options  ==! 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90

    r8058 r10473  
    44   !! Wave module  
    55   !!====================================================================== 
    6    !! History :  3.3.1  !   2011-09  (Adani M)  Original code: Drag Coefficient  
    7    !!         :  3.4    !   2012-10  (Adani M)                 Stokes Drift  
    8    !!---------------------------------------------------------------------- 
    9    USE iom             ! I/O manager library 
    10    USE in_out_manager  ! I/O manager 
    11    USE lib_mpp         ! distribued memory computing library 
    12    USE fldread        ! read input fields 
    13    USE oce 
    14    USE sbc_oce        ! Surface boundary condition: ocean fields 
    15    USE domvvl 
    16  
    17     
    18    !!---------------------------------------------------------------------- 
    19    !!   sbc_wave       : read drag coefficient from wave model in netcdf files  
    20    !!---------------------------------------------------------------------- 
     6   !! History :  3.3  !   2011-09  (Adani M)  Original code: Drag Coefficient  
     7   !!         :  3.4  !   2012-10  (Adani M)                 Stokes Drift  
     8   !!            3.6  !   2014-09  (Clementi E, Oddo P)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 oce            ! ocean variables 
     19   USE sbc_oce        ! Surface boundary condition: ocean fields 
     20   USE zdf_oce,  ONLY : ln_zdfqiao 
     21   USE bdy_oce        ! open boundary condition variables 
     22   USE domvvl         ! domain: variable volume layers 
     23   ! 
     24   USE iom            ! I/O manager library 
     25   USE in_out_manager ! I/O manager 
     26   USE lib_mpp        ! distribued memory computing library 
     27   USE fldread        ! read input fields 
     28   USE wrk_nemo       ! 
     29   USE phycst         ! physical constants  
    2130 
    2231   IMPLICIT NONE 
    2332   PRIVATE 
    2433 
    25    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 
    2637    
    27    INTEGER , PARAMETER ::   jpfld  = 3           ! maximum number of files to read for srokes drift 
    28    INTEGER , PARAMETER ::   jp_usd = 1           ! index of stokes drift  (i-component) (m/s)    at T-point 
    29    INTEGER , PARAMETER ::   jp_vsd = 2           ! index of stokes drift  (j-component) (m/s)    at T-point 
    30    INTEGER , PARAMETER ::   jp_wn  = 3           ! index of wave number                 (1/m)    at T-point 
    31    TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_cd    ! structure of input fields (file informations, fields read) Drag Coefficient 
    32    TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_sd    ! structure of input fields (file informations, fields read) Stokes Drift 
    33    REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:)       :: cdn_wave  
    34    REAL(wp),ALLOCATABLE,DIMENSION (:,:)              :: usd2d,vsd2d,uwavenum,vwavenum  
    35    REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:,:)     :: usd3d,vsd3d,wsd3d  
     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. 
    3665 
    3766   !! * Substitutions 
    38 #  include "domzgr_substitute.h90" 
     67#  include "vectopt_loop_substitute.h90" 
    3968   !!---------------------------------------------------------------------- 
    4069   !! NEMO/OPA 4.0 , NEMO Consortium (2011)  
     
    4473CONTAINS 
    4574 
     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 
     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      ! 
     96 
     97      CALL wrk_alloc( jpi,jpj,jpk,   ze3divh ) 
     98      CALL wrk_alloc( jpi,jpj,       zk_t, zk_u, zk_v, zu0_sd, zv0_sd )  
     99      ! 
     100      zfac = 2.0_wp * rpi / 16.0_wp 
     101      DO jj = 1, jpj               ! exp. wave number at t-point    (Eq. (19) in Breivick et al. (2014) ) 
     102         DO ji = 1, jpi 
     103            ! Stokes drift velocity estimated from Hs and Tmean 
     104            ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp ) 
     105            ! Stokes surface speed 
     106            tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj)) 
     107            ! Wavenumber scale 
     108            zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / 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_e12t(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( .NOT. lk_vvl ) 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#if defined key_bdy 
     170      IF( lk_bdy ) THEN 
     171         DO jk = 1, jpkm1 
     172            wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:) 
     173         END DO 
     174      ENDIF 
     175#endif 
     176      !                       !==  Horizontal divergence of barotropic Stokes transport  ==! 
     177      div_sd(:,:) = 0._wp 
     178      DO jk = 1, jpkm1                                 !  
     179        div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk) 
     180      END DO 
     181      ! 
     182      CALL iom_put( "ustokes",  usd  ) 
     183      CALL iom_put( "vstokes",  vsd  ) 
     184      CALL iom_put( "wstokes",  wsd  ) 
     185      ! 
     186      CALL wrk_dealloc( jpi,jpj,jpk,   ze3divh ) 
     187      CALL wrk_dealloc( jpi,jpj,       zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 
     188      ! 
     189   END SUBROUTINE sbc_stokes 
     190 
     191 
    46192   SUBROUTINE sbc_wave( kt ) 
    47193      !!--------------------------------------------------------------------- 
    48       !!                     ***  ROUTINE sbc_apr  *** 
    49       !! 
    50       !! ** Purpose :   read drag coefficient from wave model  in netcdf files. 
     194      !!                     ***  ROUTINE sbc_wave  *** 
     195      !! 
     196      !! ** Purpose :   read wave parameters from wave model  in netcdf files. 
    51197      !! 
    52198      !! ** Method  : - Read namelist namsbc_wave 
    53199      !!              - Read Cd_n10 fields in netcdf files  
    54200      !!              - Read stokes drift 2d in netcdf files  
    55       !!              - Read wave number      in netcdf files  
    56       !!              - Compute 3d stokes drift using monochromatic 
    57       !! ** action  :    
    58       !!                
    59       !!--------------------------------------------------------------------- 
    60       USE oce,  ONLY : un,vn,hdivn,rotn 
    61       USE divcur 
    62       USE wrk_nemo 
    63 #if defined key_bdy 
    64       USE bdy_oce, ONLY : bdytmask 
    65 #endif 
    66       INTEGER, INTENT( in  ) ::  kt       ! ocean time step 
    67       INTEGER                ::  ierror   ! return error code 
    68       INTEGER                ::  ifpr, jj,ji,jk  
    69       INTEGER                ::   ios     ! Local integer output status for namelist read 
    70       REAL(wp),DIMENSION(:,:,:),POINTER             ::  udummy,vdummy,hdivdummy,rotdummy 
    71       REAL                                          ::  z2dt,z1_2dt 
    72       TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i     ! array of namelist informations on the fields to read 
     201      !!              - Read wave number in netcdf files  
     202      !!              - Compute 3d stokes drift using Breivik et al.,2014 
     203      !!                formulation 
     204      !! ** action   
     205      !!--------------------------------------------------------------------- 
     206      INTEGER, INTENT(in   ) ::   kt   ! ocean time step 
     207      !!--------------------------------------------------------------------- 
     208      ! 
     209      IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN     !==  Neutral drag coefficient  ==! 
     210         CALL fld_read( kt, nn_fsbc, sf_cd )             ! read from external forcing 
     211         cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 
     212      ENDIF 
     213 
     214      IF( ln_tauoc .AND. .NOT. cpl_wstrf ) THEN    !==  Wave induced stress  ==! 
     215         CALL fld_read( kt, nn_fsbc, sf_tauoc )          ! read wave norm stress from external forcing 
     216         tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 
     217      ENDIF 
     218 
     219      IF( ln_sdw )  THEN                           !==  Computation of the 3d Stokes Drift  ==!  
     220         ! 
     221         IF( jpfld > 0 ) THEN                            ! Read from file only if the field is not coupled 
     222            CALL fld_read( kt, nn_fsbc, sf_sd )          ! read wave parameters from external forcing 
     223            IF( jp_hsw > 0 )   hsw  (:,:) = sf_sd(jp_hsw)%fnow(:,:,1)   ! significant wave height 
     224            IF( jp_wmp > 0 )   wmp  (:,:) = sf_sd(jp_wmp)%fnow(:,:,1)   ! wave mean period 
     225            IF( jp_usd > 0 )   ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1)   ! 2D zonal Stokes Drift at T point 
     226            IF( jp_vsd > 0 )   vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1)   ! 2D meridional Stokes Drift at T point 
     227         ENDIF 
     228         ! 
     229         ! Read also wave number if needed, so that it is available in coupling routines 
     230         IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 
     231            CALL fld_read( kt, nn_fsbc, sf_wn )          ! read wave parameters from external forcing 
     232            wnum(:,:) = sf_wn(1)%fnow(:,:,1) 
     233         ENDIF 
     234            
     235         !                                         !==  Computation of the 3d Stokes Drift  ==!  
     236         ! 
     237         IF( jpfld == 4 )   CALL sbc_stokes()            ! Calculate only if required fields are read 
     238         !                                               ! In coupled wave model-NEMO case the call is done after coupling 
     239         ! 
     240      ENDIF 
     241      ! 
     242   END SUBROUTINE sbc_wave 
     243 
     244 
     245   SUBROUTINE sbc_wave_init 
     246      !!--------------------------------------------------------------------- 
     247      !!                     ***  ROUTINE sbc_wave_init  *** 
     248      !! 
     249      !! ** Purpose :   read wave parameters from wave model  in netcdf files. 
     250      !! 
     251      !! ** Method  : - Read namelist namsbc_wave 
     252      !!              - Read Cd_n10 fields in netcdf files  
     253      !!              - Read stokes drift 2d in netcdf files  
     254      !!              - Read wave number in netcdf files  
     255      !!              - Compute 3d stokes drift using Breivik et al.,2014 
     256      !!                formulation 
     257      !! ** action   
     258      !!--------------------------------------------------------------------- 
     259      INTEGER ::   ierror, ios   ! local integer 
     260      INTEGER ::   ifpr 
     261      !! 
    73262      CHARACTER(len=100)     ::  cn_dir                          ! Root directory for location of drag coefficient files 
    74       TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd, sn_wn   ! informations about the fields to be read 
    75       !!--------------------------------------------------------------------- 
    76       NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_wn 
    77       !!--------------------------------------------------------------------- 
    78  
    79       !!---------------------------------------------------------------------- 
    80       ! 
    81       ! 
    82       !                                         ! -------------------- ! 
    83       IF( kt == nit000 ) THEN                   ! First call kt=nit000 ! 
    84          !                                      ! -------------------- ! 
    85          REWIND( numnam_ref )              ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 
    86          READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 
    87 901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 
    88  
    89          REWIND( numnam_cfg )              ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 
    90          READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 
    91 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 
    92          IF(lwm) WRITE ( numond, namsbc_wave ) 
    93          ! 
    94  
    95          IF ( ln_cdgw ) THEN 
     263      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i     ! array of namelist informations on the fields to read 
     264      TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd,  & 
     265                             &   sn_hsw, sn_wmp, sn_wnum, sn_tauoc      ! informations about the fields to be read 
     266      ! 
     267      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wnum, sn_tauoc 
     268      !!--------------------------------------------------------------------- 
     269      ! 
     270      REWIND( numnam_ref )              ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 
     271      READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 
     272901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 
     273          
     274      REWIND( numnam_cfg )              ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 
     275      READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 
     276902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 
     277      IF(lwm) WRITE ( numond, namsbc_wave ) 
     278      ! 
     279      IF( ln_cdgw ) THEN 
     280         IF( .NOT. cpl_wdrag ) THEN 
    96281            ALLOCATE( sf_cd(1), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
    97             IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     282            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 
    98283            ! 
    99284                                   ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)   ) 
    100285            IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 
    101             CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    102             ALLOCATE( cdn_wave(jpi,jpj) ) 
    103             cdn_wave(:,:) = 0.0 
    104         ENDIF 
    105          IF ( ln_sdw ) THEN 
    106             slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; slf_i(jp_wn) = sn_wn 
    107             ALLOCATE( sf_sd(3), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
    108             IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     286            CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 
     287         ENDIF 
     288         ALLOCATE( cdn_wave(jpi,jpj) ) 
     289      ENDIF 
     290 
     291      IF( ln_tauoc ) THEN 
     292         IF( .NOT. cpl_wstrf ) THEN 
     293            ALLOCATE( sf_tauoc(1), STAT=ierror )           !* allocate and fill sf_wave with sn_tauoc 
     294            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 
     295            ! 
     296                                    ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1)   ) 
     297            IF( sn_tauoc%ln_tint )  ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 
     298            CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave_init', 'Wave module', 'namsbc_wave' ) 
     299         ENDIF 
     300         ALLOCATE( tauoc_wave(jpi,jpj) ) 
     301      ENDIF 
     302 
     303      IF( ln_sdw ) THEN   ! Find out how many fields have to be read from file if not coupled 
     304         jpfld=0 
     305         jp_usd=0   ;   jp_vsd=0   ;   jp_hsw=0   ;   jp_wmp=0 
     306         IF( .NOT. cpl_sdrftx ) THEN 
     307            jpfld  = jpfld + 1 
     308            jp_usd = jpfld 
     309         ENDIF 
     310         IF( .NOT. cpl_sdrfty ) THEN 
     311            jpfld  = jpfld + 1 
     312            jp_vsd = jpfld 
     313         ENDIF 
     314         IF( .NOT. cpl_hsig ) THEN 
     315            jpfld  = jpfld + 1 
     316            jp_hsw = jpfld 
     317         ENDIF 
     318         IF( .NOT. cpl_wper ) THEN 
     319            jpfld  = jpfld + 1 
     320            jp_wmp = jpfld 
     321         ENDIF 
     322 
     323         ! Read from file only the non-coupled fields  
     324         IF( jpfld > 0 ) THEN 
     325            ALLOCATE( slf_i(jpfld) ) 
     326            IF( jp_usd > 0 )   slf_i(jp_usd) = sn_usd 
     327            IF( jp_vsd > 0 )   slf_i(jp_vsd) = sn_vsd 
     328            IF( jp_hsw > 0 )   slf_i(jp_hsw) = sn_hsw 
     329            IF( jp_wmp > 0 )   slf_i(jp_wmp) = sn_wmp 
     330            ALLOCATE( sf_sd(jpfld), STAT=ierror )   !* allocate and fill sf_sd with stokes drift 
     331            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_wave structure' ) 
    109332            ! 
    110333            DO ifpr= 1, jpfld 
     
    112335               IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 
    113336            END DO 
    114             CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    115             ALLOCATE( usd2d(jpi,jpj),vsd2d(jpi,jpj),uwavenum(jpi,jpj),vwavenum(jpi,jpj) ) 
    116             ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 
    117             usd2d(:,:) = 0.0 ;  vsd2d(:,:) = 0.0 ; uwavenum(:,:) = 0.0 ; vwavenum(:,:) = 0.0 
    118             usd3d(:,:,:) = 0.0 ;vsd3d(:,:,:) = 0.0 ; wsd3d(:,:,:) = 0.0 
    119          ENDIF 
    120       ENDIF 
    121          ! 
    122          ! 
    123       IF ( ln_cdgw ) THEN 
    124          CALL fld_read( kt, nn_fsbc, sf_cd )      !* read drag coefficient from external forcing 
    125          cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 
    126       ENDIF 
    127       IF ( ln_sdw )  THEN 
    128           CALL fld_read( kt, nn_fsbc, sf_sd )      !* read drag coefficient from external forcing 
    129  
    130          ! Interpolate wavenumber, stokes drift into the grid_V and grid_V 
    131          !------------------------------------------------- 
    132  
    133          DO jj = 1, jpjm1 
    134             DO ji = 1, jpim1 
    135                uwavenum(ji,jj)=0.5 * ( 2. - umask(ji,jj,1) ) * ( sf_sd(3)%fnow(ji,jj,1) * tmask(ji,jj,1) & 
    136                &                                + sf_sd(3)%fnow(ji+1,jj,1) * tmask(ji+1,jj,1) ) 
    137  
    138                vwavenum(ji,jj)=0.5 * ( 2. - vmask(ji,jj,1) ) * ( sf_sd(3)%fnow(ji,jj,1) * tmask(ji,jj,1) & 
    139                &                                + sf_sd(3)%fnow(ji,jj+1,1) * tmask(ji,jj+1,1) ) 
    140  
    141                usd2d(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * ( sf_sd(1)%fnow(ji,jj,1) * tmask(ji,jj,1) & 
    142                &                                + sf_sd(1)%fnow(ji+1,jj,1) * tmask(ji+1,jj,1) ) 
    143  
    144                vsd2d(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( sf_sd(2)%fnow(ji,jj,1) * tmask(ji,jj,1) & 
    145                &                                + sf_sd(2)%fnow(ji,jj+1,1) * tmask(ji,jj+1,1) ) 
    146             END DO 
    147          END DO 
    148  
    149           !Computation of the 3d Stokes Drift 
    150           DO jk = 1, jpk 
    151              DO jj = 1, jpj-1 
    152                 DO ji = 1, jpi-1 
    153                    usd3d(ji,jj,jk) = usd2d(ji,jj)*exp(2.0*uwavenum(ji,jj)*(-MIN( gdept_0(ji,jj,jk) , gdept_0(ji+1,jj  ,jk)))) 
    154                    vsd3d(ji,jj,jk) = vsd2d(ji,jj)*exp(2.0*vwavenum(ji,jj)*(-MIN( gdept_0(ji,jj,jk) , gdept_0(ji  ,jj+1,jk)))) 
    155                 END DO 
    156              END DO 
    157              usd3d(jpi,:,jk) = usd2d(jpi,:)*exp( 2.0*uwavenum(jpi,:)*(-gdept_0(jpi,:,jk)) ) 
    158              vsd3d(:,jpj,jk) = vsd2d(:,jpj)*exp( 2.0*vwavenum(:,jpj)*(-gdept_0(:,jpj,jk)) ) 
    159           END DO 
    160  
    161           CALL wrk_alloc( jpi,jpj,jpk,udummy,vdummy,hdivdummy,rotdummy) 
    162            
    163           udummy(:,:,:)=un(:,:,:) 
    164           vdummy(:,:,:)=vn(:,:,:) 
    165           hdivdummy(:,:,:)=hdivn(:,:,:) 
    166           rotdummy(:,:,:)=rotn(:,:,:) 
    167           un(:,:,:)=usd3d(:,:,:) 
    168           vn(:,:,:)=vsd3d(:,:,:) 
    169           CALL div_cur(kt) 
    170       !                                           !------------------------------! 
    171       !                                           !     Now Vertical Velocity    ! 
    172       !                                           !------------------------------! 
    173           z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog) 
    174  
    175           z1_2dt = 1.e0 / z2dt 
    176           DO jk = jpkm1, 1, -1                             ! integrate from the bottom the hor. divergence 
    177              ! - ML - need 3 lines here because replacement of fse3t by its expression yields too long lines otherwise 
    178              wsd3d(:,:,jk) = wsd3d(:,:,jk+1) -   fse3t_n(:,:,jk) * hdivn(:,:,jk)        & 
    179                 &                      - ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) )    & 
    180                 &                         * tmask(:,:,jk) * z1_2dt 
    181 #if defined key_bdy 
    182              wsd3d(:,:,jk) = wsd3d(:,:,jk) * bdytmask(:,:) 
    183 #endif 
    184           END DO 
    185           hdivn(:,:,:)=hdivdummy(:,:,:) 
    186           rotn(:,:,:)=rotdummy(:,:,:) 
    187           vn(:,:,:)=vdummy(:,:,:) 
    188           un(:,:,:)=udummy(:,:,:) 
    189           CALL wrk_dealloc( jpi,jpj,jpk,udummy,vdummy,hdivdummy,rotdummy) 
    190       ENDIF 
    191    END SUBROUTINE sbc_wave 
    192        
     337            ! 
     338            CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'Wave module ', 'namsbc_wave' ) 
     339         ENDIF 
     340         ALLOCATE( usd  (jpi,jpj,jpk), vsd  (jpi,jpj,jpk), wsd(jpi,jpj,jpk) ) 
     341         ALLOCATE( hsw  (jpi,jpj)    , wmp  (jpi,jpj)     ) 
     342         ALLOCATE( ut0sd(jpi,jpj)    , vt0sd(jpi,jpj)     ) 
     343         ALLOCATE( div_sd(jpi,jpj) ) 
     344         ALLOCATE( tsd2d (jpi,jpj) ) 
     345         usd(:,:,:) = 0._wp 
     346         vsd(:,:,:) = 0._wp 
     347         wsd(:,:,:) = 0._wp 
     348         ! Wave number needed only if ln_zdfqiao=T 
     349         IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 
     350            ALLOCATE( sf_wn(1), STAT=ierror )           !* allocate and fill sf_wave with sn_wnum 
     351            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wave structure' ) 
     352                                   ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1)   ) 
     353            IF( sn_wnum%ln_tint )  ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 
     354            CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 
     355         ENDIF 
     356         ALLOCATE( wnum(jpi,jpj) ) 
     357      ENDIF 
     358      ! 
     359   END SUBROUTINE sbc_wave_init 
     360 
    193361   !!====================================================================== 
    194362END MODULE sbcwave 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90

    r8058 r10473  
    66   !! History :  2.0  !  2005-11  (G. Madec)  Original code 
    77   !!            3.3  !  2010-09  (C. Ethe, G. Madec)  merge TRC-TRA + switch from velocity to transport 
     8   !!            3.6  !  2015-06  (E. Clementi) Addition of Stokes drift in case of wave coupling 
    89   !!            4.0  !  2011-06  (G. Madec)  Addition of Mixed Layer Eddy parameterisation 
    910   !!---------------------------------------------------------------------- 
     
    3435   USE timing          ! Timing 
    3536   USE sbc_oce 
    36    USE diaptr          ! Poleward heat transport  
    37  
     37   USE sbcwave        ! wave module   
     38   USE sbc_oce        ! surface boundary condition: ocean   
     39   USE diaptr         ! Poleward heat transport 
    3840 
    3941   IMPLICIT NONE 
     
    8587      CALL wrk_alloc( jpi, jpj, jpk, zun, zvn, zwn ) 
    8688      !                                          ! set time step 
     89      zun(:,:,:) = 0.0   
     90      zvn(:,:,:) = 0.0   
     91      zwn(:,:,:) = 0.0   
     92      ! 
    8793      IF( neuler == 0 .AND. kt == nit000 ) THEN     ! at nit000 
    8894         r2dtra(:) =  rdttra(:)                          ! = rdtra (restarting with Euler time stepping) 
     
    94100      ! 
    95101      !                                               !==  effective transport  ==! 
    96       DO jk = 1, jpkm1 
    97          zun(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)                  ! eulerian transport only 
    98          zvn(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
    99          zwn(:,:,jk) = e1t(:,:) * e2t(:,:)      * wn(:,:,jk) 
    100       END DO 
     102      IF(ln_wave .AND. ln_sdw) THEN   
     103         DO jk = 1, jpkm1                                                                      ! eulerian transport + Stokes Drift   
     104            zun(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) * ( un(:,:,jk) + usd(:,:,jk) )  
     105            zvn(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) * ( vn(:,:,jk) + vsd(:,:,jk) )   
     106            zwn(:,:,jk) = e1e2t(:,:)               * ( wn(:,:,jk) + wsd(:,:,jk) )   
     107         END DO   
     108      ELSE   
     109         DO jk = 1, jpkm1   
     110            zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)             ! eulerian transport only   
     111            zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk)   
     112            zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk)   
     113         END DO   
     114      ENDIF 
    101115      ! 
    102116      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90

    r8058 r10473  
    3535   INTEGER , PUBLIC ::   nn_npc      !: non penetrative convective scheme call  frequency 
    3636   INTEGER , PUBLIC ::   nn_npcp     !: non penetrative convective scheme print frequency 
     37   LOGICAL , PUBLIC ::   ln_zdfqiao  !: Enhanced wave vertical mixing Qiao(2010) formulation flag 
    3738 
    3839 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r8058 r10473  
    2424   USE phycst         ! physical constants 
    2525   USE zdfmxl         ! mixed layer 
     26   USE sbcwave ,  ONLY: hsw   ! significant wave height   
     27   ! 
    2628   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    2729   USE lib_mpp        ! MPP manager 
     
    197199         zdep(:,:)  = 30.*TANH(2.*0.3/(28.*SQRT(MAX(ustars2(:,:),rsmall))))             ! Wave age (eq. 10) 
    198200         zhsro(:,:) = MAX(rsbc_zs2 * ustars2(:,:) * zdep(:,:)**1.5, rn_hsro) ! zhsro = rn_frac_hs * Hsw (eq. 11) 
    199       ! 
     201      CASE ( 3 )             ! Roughness given by the wave model (coupled or read in file)   
     202         zhsro(:,:) = hsw(:,:) 
    200203      END SELECT 
    201204 
     
    909912 
    910913      !                                !* Check of some namelist values 
    911       IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'bad flag: nn_bc_surf is 0 or 1' ) 
    912       IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'bad flag: nn_bc_surf is 0 or 1' ) 
    913       IF( nn_z0_met < 0 .OR. nn_z0_met > 2 ) CALL ctl_stop( 'bad flag: nn_z0_met is 0, 1 or 2' ) 
    914       IF( nn_stab_func  < 0 .OR. nn_stab_func  > 3 ) CALL ctl_stop( 'bad flag: nn_stab_func is 0, 1, 2 and 3' ) 
    915       IF( nn_clos       < 0 .OR. nn_clos       > 3 ) CALL ctl_stop( 'bad flag: nn_clos is 0, 1, 2 or 3' ) 
     914      IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_bc_surf is 0 or 1' )   
     915      IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_bc_surf is 0 or 1' )   
     916      IF( nn_z0_met < 0 .OR. nn_z0_met > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_z0_met is 0, 1, 2 or 3' )   
     917      IF( nn_z0_met == 3 .AND. .NOT.ln_wave ) CALL ctl_stop( 'zdf_gls_init: nn_z0_met=3 requires ln_wave=T' )   
     918      IF( nn_stab_func  < 0 .OR. nn_stab_func  > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_stab_func is 0, 1, 2 and 3' )   
     919      IF( nn_clos       < 0 .OR. nn_clos       > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_clos is 0, 1, 2 or 3' ) 
    916920 
    917921      SELECT CASE ( nn_clos )          !* set the parameters for the chosen closure 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfini.F90

    r8058 r10473  
    5353      INTEGER ::   ios 
    5454      !! 
    55       NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp,   & 
    56          &              ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp 
     55      NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp,  &   
     56         &        ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp,       &   
     57         &        ln_zdfqiao 
    5758      !!---------------------------------------------------------------------- 
    5859 
     
    8384         WRITE(numout,*) '      npc call  frequency                 nn_npc    = ', nn_npc 
    8485         WRITE(numout,*) '      npc print frequency                 nn_npcp   = ', nn_npcp 
     86         WRITE(numout,*) '      Qiao formulation flag               ln_zdfqiao=', ln_zdfqiao 
    8587      ENDIF 
    8688 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/step.F90

    r10394 r10473  
    2525   !!            3.4  !  2011-04  (G. Madec, C. Ethe) Merge of dtatem and dtasal 
    2626   !!                 !  2012-07  (J. Simeon, G. Madec, C. Ethe) Online coarsening of outputs 
     27   !!            3.6  !  2014-10  (E. Clementi, P. Oddo) Add Qiao vertical mixing in case of waves 
    2728   !!            3.7  !  2014-04  (F. Roquet, G. Madec) New equations of state 
    2829   !!---------------------------------------------------------------------- 
     
    7273      !!              -8- Outputs and diagnostics 
    7374      !!---------------------------------------------------------------------- 
    74       INTEGER ::   jk      ! dummy loop indice 
     75      INTEGER ::   ji,jj,jk ! dummy loop indice 
    7576      INTEGER ::   indic    ! error indicator if < 0 
    7677      INTEGER ::   kcall    ! optional integer argument (dom_vvl_sf_nxt) 
     
    132133      IF( lk_zdftke  )   CALL zdf_tke( kstp )            ! TKE closure scheme for Kz 
    133134      IF( lk_zdfgls  )   CALL zdf_gls( kstp )            ! GLS closure scheme for Kz 
     135      IF( ln_zdfqiao )   CALL zdf_qiao( kstp )           ! Qiao vertical mixing    
     136      ! 
    134137      IF( lk_zdfkpp  )   CALL zdf_kpp( kstp )            ! KPP closure scheme for Kz 
    135138      IF( lk_zdfcst  ) THEN                              ! Constant Kz (reset avt, avm[uv] to the background value) 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/step_oce.F90

    r8561 r10473  
    2929   USE sbctide          ! Tide initialisation 
    3030   USE sbcapr           ! surface boundary condition: ssh_ib required by bdydta  
     31   USE sbcwave          ! Wave intialisation 
    3132 
    3233   USE traqsr           ! solar radiation penetration      (tra_qsr routine) 
     
    8687   USE zdfric           ! Richardson vertical mixing       (zdf_ric routine) 
    8788   USE zdfmxl           ! Mixed-layer depth                (zdf_mxl routine) 
     89   USE zdfqiao          !Qiao module wave induced mixing   (zdf_qiao routine 
    8890 
    8991   USE zpshde           ! partial step: hor. derivative     (zps_hde routine) 
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/SETTE/sette.sh

    r5588 r10473  
    8888# 
    8989# Compiler among those in NEMOGCM/ARCH 
    90 COMPILER=X64_ADA 
    91 export BATCH_COMMAND_PAR="llsubmit" 
     90module load cray-netcdf-hdf5parallel   
     91COMPILER=XC40_METO   
     92export BATCH_COMMAND_PAR="qsub 
    9293export BATCH_COMMAND_SEQ=$BATCH_COMMAND_PAR 
    9394export INTERACT_FLAG="no" 
Note: See TracChangeset for help on using the changeset viewer.