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 14037 for NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/SBC/sbcmod.F90 – NEMO

Ignore:
Timestamp:
2020-12-03T12:20:38+01:00 (3 years ago)
Author:
ayoung
Message:

Updated to trunk at 14020. Sette tests passed with change of results for configurations with non-linear ssh. Ticket #2506.

Location:
NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG

    • Property svn:externals
      •  

        old new  
        88 
        99# SETTE 
        10 ^/utils/CI/sette@13292        sette 
         10^/utils/CI/sette_wave@13990         sette 
  • NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/SBC/sbcmod.F90

    r13286 r14037  
    1616   !!            4.0  ! 2016-06  (L. Brodeau) new general bulk formulation 
    1717   !!            4.0  ! 2019-03  (F. Lemarié & G. Samson)  add ABL compatibility (ln_abl=TRUE) 
     18   !!            4.2  ! 2020-12  (G. Madec, E. Clementi) modified wave forcing and coupling   
    1819   !!---------------------------------------------------------------------- 
    1920 
     
    5455   USE usrdef_sbc     ! user defined: surface boundary condition 
    5556   USE closea         ! closed sea 
     57   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
    5658   ! 
    5759   USE prtctl         ! Print control                    (prt_ctl routine) 
     
    7072 
    7173   INTEGER ::   nsbc   ! type of surface boundary condition (deduced from namsbc informations) 
    72  
     74   !! * Substitutions 
     75#  include "do_loop_substitute.h90" 
    7376   !!---------------------------------------------------------------------- 
    7477   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    99102         &             nn_ice   , ln_ice_embd,                                       & 
    100103         &             ln_traqsr, ln_dm2dc ,                                         & 
    101          &             ln_rnf   , nn_fwb     , ln_ssr   , ln_apr_dyn,              & 
    102          &             ln_wave  , ln_cdgw  , ln_sdw   , ln_tauwoc  , ln_stcor  ,     & 
    103          &             ln_tauw  , nn_lsm, nn_sdrift 
     104         &             ln_rnf   , nn_fwb     , ln_ssr   , ln_apr_dyn,                & 
     105         &             ln_wave  , nn_lsm 
    104106      !!---------------------------------------------------------------------- 
    105107      ! 
     
    119121#if defined key_mpp_mpi 
    120122      ncom_fsbc = nn_fsbc    ! make nn_fsbc available for lib_mpp 
     123#endif 
     124#if ! defined key_si3 
     125      IF( nn_ice == 2 )    nn_ice = 0  ! without key key_si3 you cannot use si3... 
    121126#endif 
    122127      ! 
     
    130135         WRITE(numout,*) '         bulk         formulation                   ln_blk        = ', ln_blk 
    131136         WRITE(numout,*) '         ABL          formulation                   ln_abl        = ', ln_abl 
     137         WRITE(numout,*) '         Surface wave (forced or coupled)           ln_wave       = ', ln_wave 
    132138         WRITE(numout,*) '      Type of coupling (Ocean/Ice/Atmosphere) : ' 
    133139         WRITE(numout,*) '         ocean-atmosphere coupled formulation       ln_cpl        = ', ln_cpl 
     
    147153         WRITE(numout,*) '         runoff / runoff mouths                     ln_rnf        = ', ln_rnf 
    148154         WRITE(numout,*) '         nb of iterations if land-sea-mask applied  nn_lsm        = ', nn_lsm 
    149          WRITE(numout,*) '         surface wave                               ln_wave       = ', ln_wave 
    150          WRITE(numout,*) '               Stokes drift corr. to vert. velocity ln_sdw        = ', ln_sdw 
    151          WRITE(numout,*) '                  vertical parametrization          nn_sdrift     = ', nn_sdrift 
    152          WRITE(numout,*) '               wave modified ocean stress           ln_tauwoc     = ', ln_tauwoc 
    153          WRITE(numout,*) '               wave modified ocean stress component ln_tauw       = ', ln_tauw 
    154          WRITE(numout,*) '               Stokes coriolis term                 ln_stcor      = ', ln_stcor 
    155          WRITE(numout,*) '               neutral drag coefficient (CORE,NCAR) ln_cdgw       = ', ln_cdgw 
    156       ENDIF 
    157       ! 
    158       IF( .NOT.ln_wave ) THEN 
    159          ln_sdw = .false. ; ln_cdgw = .false. ; ln_tauwoc = .false. ; ln_tauw = .false. ; ln_stcor = .false. 
    160       ENDIF  
    161       IF( ln_sdw ) THEN 
    162          IF( .NOT.(nn_sdrift==jp_breivik_2014 .OR. nn_sdrift==jp_li_2017 .OR. nn_sdrift==jp_peakfr) ) & 
    163             CALL ctl_stop( 'The chosen nn_sdrift for Stokes drift vertical velocity must be 0, 1, or 2' ) 
    164       ENDIF 
    165       ll_st_bv2014  = ( nn_sdrift==jp_breivik_2014 ) 
    166       ll_st_li2017  = ( nn_sdrift==jp_li_2017 ) 
    167       ll_st_bv_li   = ( ll_st_bv2014 .OR. ll_st_li2017 ) 
    168       ll_st_peakfr  = ( nn_sdrift==jp_peakfr ) 
    169       IF( ln_tauwoc .AND. ln_tauw ) & 
    170          CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 
    171                                   '(ln_tauwoc=.true. and ln_tauw=.true.)' ) 
    172       IF( ln_tauwoc ) & 
    173          CALL ctl_warn( 'You are subtracting the wave stress to the ocean (ln_tauwoc=.true.)' ) 
    174       IF( ln_tauw ) & 
    175          CALL ctl_warn( 'The wave modified ocean stress components are used (ln_tauw=.true.) ', & 
    176                               'This will override any other specification of the ocean stress' ) 
     155      ENDIF 
    177156      ! 
    178157      IF( .NOT.ln_usr ) THEN     ! the model calendar needs some specificities (except in user defined case) 
     
    226205      CASE DEFAULT                     !- not supported 
    227206      END SELECT 
    228       IF( ln_diurnal .AND. .NOT. ln_blk )   CALL ctl_stop( "sbc_init: diurnal flux processing only implemented for bulk forcing" ) 
     207      IF( ln_diurnal .AND. .NOT. (ln_blk.OR.ln_abl) )   CALL ctl_stop( "sbc_init: diurnal flux processing only implemented for bulk forcing" ) 
    229208      ! 
    230209      !                       !**  allocate and set required variables 
     
    243222      ENDIF 
    244223      ! 
    245  
    246224      IF( nn_ice == 0 ) THEN        !* No sea-ice in the domain : ice fraction is always zero 
    247225         IF( nn_components /= jp_iam_opa )   fr_i(:,:) = 0._wp    ! except for OPA in SAS-OPA coupled case 
     
    250228      sfx   (:,:) = 0._wp           !* salt flux due to freezing/melting 
    251229      fmmflx(:,:) = 0._wp           !* freezing minus melting flux 
     230      cloud_fra(:,:) = pp_cldf      !* cloud fraction over sea ice (used in si3) 
    252231 
    253232      taum(:,:) = 0._wp             !* wind stress module (needed in GLS in case of reduced restart) 
     
    334313      IF( l_sbc_clo   )   CALL sbc_clo_init              ! closed sea surface initialisation 
    335314      ! 
    336       IF( ln_blk      )   CALL sbc_blk_init            ! bulk formulae initialization 
    337  
    338       IF( ln_abl      )   CALL sbc_abl_init            ! Atmospheric Boundary Layer (ABL) 
    339  
    340       IF( ln_ssr      )   CALL sbc_ssr_init            ! Sea-Surface Restoring initialization 
     315      IF( ln_blk      )   CALL sbc_blk_init              ! bulk formulae initialization 
     316 
     317      IF( ln_abl      )   CALL sbc_abl_init              ! Atmospheric Boundary Layer (ABL) 
     318 
     319      IF( ln_ssr      )   CALL sbc_ssr_init              ! Sea-Surface Restoring initialization 
    341320      ! 
    342321      ! 
     
    354333      IF( nn_ice == 3 )   CALL cice_sbc_init( nsbc, Kbb, Kmm )   ! CICE initialization 
    355334      ! 
    356       IF( ln_wave     )   CALL sbc_wave_init                     ! surface wave initialisation 
    357       ! 
    358       IF( lwxios ) THEN 
    359          CALL iom_set_rstw_var_active('utau_b') 
    360          CALL iom_set_rstw_var_active('vtau_b') 
    361          CALL iom_set_rstw_var_active('qns_b') 
    362          ! The 3D heat content due to qsr forcing is treated in traqsr 
    363          ! CALL iom_set_rstw_var_active('qsr_b') 
    364          CALL iom_set_rstw_var_active('emp_b') 
    365          CALL iom_set_rstw_var_active('sfx_b') 
    366       ENDIF 
    367  
     335      IF( ln_wave     ) THEN 
     336                          CALL sbc_wave_init                     ! surface wave initialisation 
     337      ELSE 
     338                          IF(lwp) WRITE(numout,*) 
     339                          IF(lwp) WRITE(numout,*) '   No surface waves : all wave related logical set to false' 
     340                          ln_sdw       = .false. 
     341                          ln_stcor     = .false. 
     342                          ln_cdgw      = .false. 
     343                          ln_tauoc     = .false. 
     344                          ln_wave_test = .false. 
     345                          ln_charn     = .false. 
     346                          ln_taw       = .false. 
     347                          ln_phioc     = .false. 
     348                          ln_bern_srfc = .false. 
     349                          ln_breivikFV_2016 = .false. 
     350                          ln_vortex_force = .false. 
     351                          ln_stshear  = .false. 
     352      ENDIF 
     353      ! 
    368354   END SUBROUTINE sbc_init 
    369355 
     
    387373      INTEGER, INTENT(in) ::   kt   ! ocean time step 
    388374      INTEGER, INTENT(in) ::   Kbb, Kmm   ! ocean time level indices 
     375      INTEGER  ::   jj, ji          ! dummy loop argument 
    389376      ! 
    390377      LOGICAL ::   ll_sas, ll_opa   ! local logical 
     
    419406      ! 
    420407      IF( .NOT.ll_sas )   CALL sbc_ssm ( kt, Kbb, Kmm )  ! mean ocean sea surface variables (sst_m, sss_m, ssu_m, ssv_m) 
    421       IF( ln_wave     )   CALL sbc_wave( kt, Kmm )       ! surface waves 
    422  
    423408      ! 
    424409      !                                            !==  sbc formulation  ==! 
    425410      !                                                    
     411      ! 
    426412      SELECT CASE( nsbc )                                ! Compute ocean surface boundary condition 
    427413      !                                                  ! (i.e. utau,vtau, qns, qsr, emp, sfx) 
     
    430416      CASE( jp_blk     ) 
    431417         IF( ll_sas    )       CALL sbc_cpl_rcv   ( kt, nn_fsbc, nn_ice, Kbb, Kmm )   ! OPA-SAS coupling: SAS receiving fields from OPA 
     418!!!!!!!!!!! ATTENTION:ln_wave is not only used for oasis coupling !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
     419         IF( ln_wave )   THEN 
     420             IF ( lk_oasis )  CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice, Kbb, Kmm )   ! OPA-wave coupling 
     421             CALL sbc_wave ( kt, Kmm ) 
     422         ENDIF 
    432423                               CALL sbc_blk       ( kt )                    ! bulk formulation for the ocean 
    433424                               ! 
     
    443434      IF( ln_mixcpl )          CALL sbc_cpl_rcv   ( kt, nn_fsbc, nn_ice, Kbb, Kmm )  ! forced-coupled mixed formulation after forcing 
    444435      ! 
    445       IF ( ln_wave .AND. (ln_tauwoc .OR. ln_tauw) ) CALL sbc_wstress( )              ! Wind stress provided by waves  
     436      IF( ln_wave .AND. ln_tauoc )  THEN            ! Wave stress reduction 
     437         DO_2D( 0, 0, 0, 0) 
     438            utau(ji,jj) = utau(ji,jj) * ( tauoc_wave(ji,jj) + tauoc_wave(ji-1,jj) ) * 0.5_wp 
     439            vtau(ji,jj) = vtau(ji,jj) * ( tauoc_wave(ji,jj) + tauoc_wave(ji,jj-1) ) * 0.5_wp 
     440         END_2D 
     441         ! 
     442         CALL lbc_lnk( 'sbcwave', utau, 'U', -1. ) 
     443         CALL lbc_lnk( 'sbcwave', vtau, 'V', -1. ) 
     444         ! 
     445         taum(:,:) = taum(:,:)*tauoc_wave(:,:) 
     446         ! 
     447         IF( kt == nit000 )   CALL ctl_warn( 'sbc: You are subtracting the wave stress to the ocean.',   & 
     448            &                                'If not requested select ln_tauoc=.false.' ) 
     449         ! 
     450      ELSEIF( ln_wave .AND. ln_taw ) THEN                  ! Wave stress reduction 
     451         utau(:,:) = utau(:,:) - tawx(:,:) + twox(:,:) 
     452         vtau(:,:) = vtau(:,:) - tawy(:,:) + twoy(:,:) 
     453         CALL lbc_lnk( 'sbcwave', utau, 'U', -1. ) 
     454         CALL lbc_lnk( 'sbcwave', vtau, 'V', -1. ) 
     455         ! 
     456         DO_2D( 0, 0, 0, 0) 
     457             taum(ji,jj) = sqrt((.5*(utau(ji-1,jj)+utau(ji,jj)))**2 + (.5*(vtau(ji,jj-1)+vtau(ji,jj)))**2) 
     458         END_2D 
     459         ! 
     460         IF( kt == nit000 )   CALL ctl_warn( 'sbc: You are subtracting the wave stress to the ocean.',   & 
     461            &                                'If not requested select ln_taw=.false.' ) 
     462         ! 
     463      ENDIF 
     464      CALL lbc_lnk( 'sbcmod', taum(:,:), 'T', 1. ) 
    446465      ! 
    447466      !                                            !==  Misc. Options  ==! 
     
    507526            & iom_varid( numror, 'utau_b', ldstop = .FALSE. ) > 0 ) THEN 
    508527            IF(lwp) WRITE(numout,*) '          nit000-1 surface forcing fields red in the restart file' 
    509             CALL iom_get( numror, jpdom_auto, 'utau_b', utau_b, ldxios = lrxios, cd_type = 'U', psgn = -1._wp )   ! before i-stress  (U-point) 
    510             CALL iom_get( numror, jpdom_auto, 'vtau_b', vtau_b, ldxios = lrxios, cd_type = 'V', psgn = -1._wp )   ! before j-stress  (V-point) 
    511             CALL iom_get( numror, jpdom_auto,  'qns_b',  qns_b, ldxios = lrxios )   ! before non solar heat flux (T-point) 
     528            CALL iom_get( numror, jpdom_auto, 'utau_b', utau_b )   ! before i-stress  (U-point) 
     529            CALL iom_get( numror, jpdom_auto, 'vtau_b', vtau_b )   ! before j-stress  (V-point) 
     530            CALL iom_get( numror, jpdom_auto,  'qns_b',  qns_b )   ! before non solar heat flux (T-point) 
    512531            ! The 3D heat content due to qsr forcing is treated in traqsr 
    513             ! CALL iom_get( numror, jpdom_auto, 'qsr_b' , qsr_b, ldxios = lrxios  ) ! before     solar heat flux (T-point) 
    514             CALL iom_get( numror, jpdom_auto, 'emp_b', emp_b, ldxios = lrxios  )    ! before     freshwater flux (T-point) 
     532            ! CALL iom_get( numror, jpdom_auto, 'qsr_b' , qsr_b  ) ! before     solar heat flux (T-point) 
     533            CALL iom_get( numror, jpdom_auto, 'emp_b', emp_b  )    ! before     freshwater flux (T-point) 
    515534            ! To ensure restart capability with 3.3x/3.4 restart files    !! to be removed in v3.6 
    516535            IF( iom_varid( numror, 'sfx_b', ldstop = .FALSE. ) > 0 ) THEN 
    517                CALL iom_get( numror, jpdom_auto, 'sfx_b', sfx_b, ldxios = lrxios )  ! before salt flux (T-point) 
     536               CALL iom_get( numror, jpdom_auto, 'sfx_b', sfx_b )  ! before salt flux (T-point) 
    518537            ELSE 
    519538               sfx_b (:,:) = sfx(:,:) 
     
    535554            &                    'at it= ', kt,' date= ', ndastp 
    536555         IF(lwp) WRITE(numout,*) '~~~~' 
    537          IF( lwxios ) CALL iom_swap(      cwxios_context          ) 
    538          CALL iom_rstput( kt, nitrst, numrow, 'utau_b' , utau, ldxios = lwxios ) 
    539          CALL iom_rstput( kt, nitrst, numrow, 'vtau_b' , vtau, ldxios = lwxios ) 
    540          CALL iom_rstput( kt, nitrst, numrow, 'qns_b'  , qns, ldxios = lwxios  ) 
     556         CALL iom_rstput( kt, nitrst, numrow, 'utau_b' , utau ) 
     557         CALL iom_rstput( kt, nitrst, numrow, 'vtau_b' , vtau ) 
     558         CALL iom_rstput( kt, nitrst, numrow, 'qns_b'  , qns  ) 
    541559         ! The 3D heat content due to qsr forcing is treated in traqsr 
    542560         ! CALL iom_rstput( kt, nitrst, numrow, 'qsr_b'  , qsr  ) 
    543          CALL iom_rstput( kt, nitrst, numrow, 'emp_b'  , emp, ldxios = lwxios  ) 
    544          CALL iom_rstput( kt, nitrst, numrow, 'sfx_b'  , sfx, ldxios = lwxios  ) 
    545          IF( lwxios ) CALL iom_swap(      cxios_context          ) 
     561         CALL iom_rstput( kt, nitrst, numrow, 'emp_b'  , emp  ) 
     562         CALL iom_rstput( kt, nitrst, numrow, 'sfx_b'  , sfx  ) 
    546563      ENDIF 
    547564      !                                                ! ---------------------------------------- ! 
     
    563580      ENDIF 
    564581      ! 
    565       CALL iom_put( "utau", utau )   ! i-wind stress   (stress can be updated at each time step in sea-ice) 
    566       CALL iom_put( "vtau", vtau )   ! j-wind stress 
    567       ! 
    568582      IF(sn_cfctl%l_prtctl) THEN     ! print mean trends (used for debugging) 
    569583         CALL prt_ctl(tab2d_1=fr_i             , clinfo1=' fr_i     - : ', mask1=tmask ) 
Note: See TracChangeset for help on using the changeset viewer.