Changeset 12583


Ignore:
Timestamp:
2020-03-21T15:40:52+01:00 (3 weeks ago)
Author:
techene
Message:

OCE/DOM/domqe.F90: add gdep at time level Kbb in dom_qe_sf_update, OCE/DOM/domzgr_substitute.h90: create the substitute module, OCE/DYN/dynatfLF.F90, OCE/TRA/traatfLF.F90: move boundary condition management and agrif management from atf modules to OCE/steplf.F90, OCE/SBC/sbcice_cice.F90, ICE/iceistate.F90 : remove dom_vvl_interpol and replace by dom_vvl_zgr ?

Location:
NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src
Files:
1 added
6 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/ICE/iceistate.F90

    r12377 r12583  
    1818   USE oce            ! dynamics and tracers variables 
    1919   USE dom_oce        ! ocean domain 
    20    USE sbc_oce , ONLY : sst_m, sss_m, ln_ice_embd  
     20   USE sbc_oce , ONLY : sst_m, sss_m, ln_ice_embd 
    2121   USE sbc_ice , ONLY : tn_ice, snwice_mass, snwice_mass_b 
    2222   USE eosbn2         ! equation of state 
     
    6060   INTEGER , PARAMETER ::   jp_hpd = 9           ! index of pnd depth        (m) 
    6161   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   si  ! structure of input fields (file informations, fields read) 
    62    !    
     62   ! 
    6363   !! * Substitutions 
    6464#  include "do_loop_substitute.h90" 
     
    7777      !! 
    7878      !! ** Method  :   This routine will put some ice where ocean 
    79       !!                is at the freezing point, then fill in ice  
    80       !!                state variables using prescribed initial  
    81       !!                values in the namelist             
     79      !!                is at the freezing point, then fill in ice 
     80      !!                state variables using prescribed initial 
     81      !!                values in the namelist 
    8282      !! 
    8383      !! ** Steps   :   1) Set initial surface and basal temperatures 
     
    9191      !!              where there is no ice 
    9292      !!-------------------------------------------------------------------- 
    93       INTEGER, INTENT(in) :: kt            ! time step  
     93      INTEGER, INTENT(in) :: kt            ! time step 
    9494      INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 
    9595      ! 
     
    117117      ! basal temperature (considered at freezing point)   [Kelvin] 
    118118      CALL eos_fzp( sss_m(:,:), t_bo(:,:) ) 
    119       t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1)  
     119      t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1) 
    120120      ! 
    121121      ! surface temperature and conductivity 
     
    142142      e_i (:,:,:,:) = 0._wp 
    143143      e_s (:,:,:,:) = 0._wp 
    144        
     144 
    145145      ! general fields 
    146146      a_i (:,:,:) = 0._wp 
     
    215215            IF( TRIM(si(jp_apd)%clrootname) == 'NOT USED' ) & 
    216216               &     si(jp_apd)%fnow(:,:,1) = ( rn_apd_ini_n * zswitch + rn_apd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) & ! rn_apd = pond fraction => rn_apnd * a_i = pond conc. 
    217                &                              * si(jp_ati)%fnow(:,:,1)  
     217               &                              * si(jp_ati)%fnow(:,:,1) 
    218218            zapnd_ini(:,:) = si(jp_apd)%fnow(:,:,1) 
    219219            ! 
     
    224224            ! 
    225225            ! change the switch for the following 
    226             WHERE( zat_i_ini(:,:) > 0._wp )   ;   zswitch(:,:) = tmask(:,:,1)  
     226            WHERE( zat_i_ini(:,:) > 0._wp )   ;   zswitch(:,:) = tmask(:,:,1) 
    227227            ELSEWHERE                         ;   zswitch(:,:) = 0._wp 
    228228            END WHERE 
     
    231231            !                          !---------------! 
    232232            ! no ice if (sst - Tfreez) >= thresold 
    233             WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst )   ;   zswitch(:,:) = 0._wp  
     233            WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst )   ;   zswitch(:,:) = 0._wp 
    234234            ELSEWHERE                                                                    ;   zswitch(:,:) = tmask(:,:,1) 
    235235            END WHERE 
     
    244244               zt_su_ini(:,:) = rn_tsu_ini_n * zswitch(:,:) 
    245245               ztm_s_ini(:,:) = rn_tms_ini_n * zswitch(:,:) 
    246                zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc.  
     246               zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 
    247247               zhpnd_ini(:,:) = rn_hpd_ini_n * zswitch(:,:) 
    248248            ELSEWHERE 
     
    265265            zhpnd_ini(:,:) = 0._wp 
    266266         ENDIF 
    267           
     267 
    268268         !-------------! 
    269269         ! fill fields ! 
     
    292292         ALLOCATE( zhi_2d(npti,jpl), zhs_2d(npti,jpl), zai_2d (npti,jpl), & 
    293293            &      zti_2d(npti,jpl), zts_2d(npti,jpl), ztsu_2d(npti,jpl), zsi_2d(npti,jpl), zaip_2d(npti,jpl), zhip_2d(npti,jpl) ) 
    294           
     294 
    295295         ! distribute 1-cat into jpl-cat: (jpi*jpj) -> (jpi*jpj,jpl) 
    296296         CALL ice_var_itd( h_i_1d(1:npti)  , h_s_1d(1:npti)  , at_i_1d(1:npti),                                                   & 
     
    338338         DO jl = 1, jpl 
    339339            DO_3D_11_11( 1, nlay_i ) 
    340                t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl)  
     340               t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl) 
    341341               ztmelts          = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 ! melting temperature in K 
    342342               e_i(ji,jj,jk,jl) = zswitch(ji,jj) * v_i(ji,jj,jl) * r1_nlay_i * & 
     
    354354         END WHERE 
    355355         v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:) 
    356            
     356 
    357357         ! specific temperatures for coupled runs 
    358358         tn_ice(:,:,:) = t_su(:,:,:) 
     
    374374         ssh(:,:,Kbb) = ssh(:,:,Kbb) - snwice_mass(:,:) * r1_rau0 
    375375         ! 
    376          IF( .NOT.ln_linssh ) THEN 
    377             ! 
    378             WHERE( ht_0(:,:) > 0 )   ;   z2d(:,:) = 1._wp + ssh(:,:,Kmm)*tmask(:,:,1) / ht_0(:,:) 
    379             ELSEWHERE                ;   z2d(:,:) = 1._wp   ;   END WHERE 
    380             ! 
    381             DO jk = 1,jpkm1                     ! adjust initial vertical scale factors                 
    382                e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * z2d(:,:) 
    383                e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 
    384                e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kmm) 
    385             END DO 
    386             ! 
    387             ! Reconstruction of all vertical scale factors at now and before time-steps 
    388             ! ========================================================================= 
    389             ! Horizontal scale factor interpolations 
    390             ! -------------------------------------- 
    391             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 
    392             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 
    393             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
    394             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
    395             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 
    396             ! Vertical scale factor interpolations 
    397             ! ------------------------------------ 
    398             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  ) 
    399             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
    400             CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
    401             CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
    402             CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    403             ! t- and w- points depth 
    404             ! ---------------------- 
    405             !!gm not sure of that.... 
    406             gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
    407             gdepw(:,:,1,Kmm) = 0.0_wp 
    408             gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    409             DO jk = 2, jpk 
    410                gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk  ,Kmm) 
    411                gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 
    412                gde3w(:,:,jk) = gdept(:,:,jk  ,Kmm) - ssh (:,:,Kmm) 
    413             END DO 
    414          ENDIF 
     376         IF( .NOT.ln_linssh )   CALL dom_vvl_zgr( Kbb, Kmm, Kaa )   ! interpolation scale factor, depth and water column 
     377! !!st 
     378!          IF( .NOT.ln_linssh ) THEN 
     379!             ! 
     380!             WHERE( ht_0(:,:) > 0 )   ;   z2d(:,:) = 1._wp + ssh(:,:,Kmm)*tmask(:,:,1) / ht_0(:,:) 
     381!             ELSEWHERE                ;   z2d(:,:) = 1._wp   ;   END WHERE 
     382!             ! 
     383!             DO jk = 1,jpkm1                     ! adjust initial vertical scale factors 
     384!                e3t(:,:,jk,Kmm) = e3t_0(:,:,jk) * z2d(:,:) 
     385!                e3t(:,:,jk,Kbb) = e3t(:,:,jk,Kmm) 
     386!                e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kmm) 
     387!             END DO 
     388!             ! 
     389!             ! Reconstruction of all vertical scale factors at now and before time-steps 
     390!             ! ========================================================================= 
     391!             ! Horizontal scale factor interpolations 
     392!             ! -------------------------------------- 
     393!             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 
     394!             CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 
     395!             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
     396!             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
     397!             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 
     398!             ! Vertical scale factor interpolations 
     399!             ! ------------------------------------ 
     400!             CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  ) 
     401!             CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
     402!             CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
     403!             CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     404!             CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
     405!             ! t- and w- points depth 
     406!             ! ---------------------- 
     407!             !!gm not sure of that.... 
     408!             gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
     409!             gdepw(:,:,1,Kmm) = 0.0_wp 
     410!             gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
     411!             DO jk = 2, jpk 
     412!                gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk  ,Kmm) 
     413!                gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 
     414!                gde3w(:,:,jk) = gdept(:,:,jk  ,Kmm) - ssh (:,:,Kmm) 
     415!             END DO 
     416!          ENDIF 
    415417      ENDIF 
    416        
     418 
    417419      !------------------------------------ 
    418420      ! 4) store fields at before time-step 
     
    429431      v_ice_b(:,:)     = v_ice(:,:) 
    430432      ! total concentration is needed for Lupkes parameterizations 
    431       at_i_b (:,:)     = at_i (:,:)  
     433      at_i_b (:,:)     = at_i (:,:) 
    432434 
    433435!!clem: output of initial state should be written here but it is impossible because 
     
    441443      !!------------------------------------------------------------------- 
    442444      !!                   ***  ROUTINE ice_istate_init  *** 
    443       !!         
    444       !! ** Purpose :   Definition of initial state of the ice  
    445       !! 
    446       !! ** Method  :   Read the namini namelist and check the parameter  
     445      !! 
     446      !! ** Purpose :   Definition of initial state of the ice 
     447      !! 
     448      !! ** Method  :   Read the namini namelist and check the parameter 
    447449      !!              values called at the first timestep (nit000) 
    448450      !! 
     
    485487         WRITE(numout,*) '      max ocean temp. above Tfreeze with initial ice   rn_thres_sst   = ', rn_thres_sst 
    486488         IF( ln_iceini .AND. .NOT.ln_iceini_file ) THEN 
    487             WRITE(numout,*) '      initial snw thickness in the north-south         rn_hts_ini     = ', rn_hts_ini_n,rn_hts_ini_s  
     489            WRITE(numout,*) '      initial snw thickness in the north-south         rn_hts_ini     = ', rn_hts_ini_n,rn_hts_ini_s 
    488490            WRITE(numout,*) '      initial ice thickness in the north-south         rn_hti_ini     = ', rn_hti_ini_n,rn_hti_ini_s 
    489491            WRITE(numout,*) '      initial ice concentr  in the north-south         rn_ati_ini     = ', rn_ati_ini_n,rn_ati_ini_s 
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/DOM/domqe.F90

    r12581 r12583  
    1313 
    1414   !!---------------------------------------------------------------------- 
    15    !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
    16    !!   dom_vvl_sf_nxt   : Compute next vertical scale factors 
    17    !!   dom_vvl_sf_update   : Swap vertical scale factors and update the vertical grid 
    18    !!   dom_vvl_interpol : Interpolate vertical scale factors from one grid point to another 
    19    !!   dom_vvl_r3c      : Compute ssh/h_0 ratioat t-, u-, v-, and optionally f-points 
    20    !!   dom_vvl_rst      : read/write restart file 
    21    !!   dom_vvl_ctl      : Check the vvl options 
     15   !!   dom_qe_init     : define initial vertical scale factors, depths and column thickness 
     16   !!   dom_qe_sf_nxt   : Compute next vertical scale factors 
     17   !!   dom_qe_sf_update   : Swap vertical scale factors and update the vertical grid 
     18   !!   dom_qe_interpol : Interpolate vertical scale factors from one grid point to another 
     19   !!   dom_qe_r3c      : Compute ssh/h_0 ratioat t-, u-, v-, and optionally f-points 
     20   !!   dom_qe_rst      : read/write restart file 
     21   !!   dom_qe_ctl      : Check the vvl options 
    2222   !!---------------------------------------------------------------------- 
    2323   USE oce             ! ocean dynamics and tracers 
     
    374374                              + MAX(   0._wp    , gdepw_0(:,:,jk)-risfdep(:,:) ) * ( 1._wp + r3t(:,:,Kmm) ) 
    375375            gde3w(:,:,jk)     = gdept(:,:,jk,Kmm) - ssh(:,:,Kmm) 
     376            gdept(:,:,jk,Kbb) = MIN( risfdep(:,:) , gdept_0(:,:,jk) )            & 
     377                              + MAX(   0._wp    , gdept_0(:,:,jk)-risfdep(:,:) ) * ( 1._wp + r3t(:,:,Kbb) ) 
     378            gdepw(:,:,jk,Kbb) = MIN( risfdep(:,:) , gdepw_0(:,:,jk) )    & 
     379                              + MAX(   0._wp    , gdepw_0(:,:,jk)-risfdep(:,:) ) * ( 1._wp + r3t(:,:,Kbb) ) 
    376380         END DO 
    377381         ! 
     
    383387            gdepw(:,:,jk,Kmm) = gdepw_0(:,:,jk) * ( 1._wp + r3t(:,:,Kmm) ) 
    384388            gde3w(:,:,jk)     = gdept  (:,:,jk,Kmm)       - ssh(:,:,Kmm) 
     389            gdept(:,:,jk,Kbb) = gdept_0(:,:,jk) * ( 1._wp + r3t(:,:,Kbb) ) 
     390            gdepw(:,:,jk,Kbb) = gdepw_0(:,:,jk) * ( 1._wp + r3t(:,:,Kbb) ) 
    385391         END DO 
    386392         ! 
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/DYN/dynatfLF.F90

    r12581 r12583  
    116116      ENDIF 
    117117 
    118       IF ( ln_dynspg_ts ) THEN 
    119          ! Ensure below that barotropic velocities match time splitting estimate 
    120          ! Compute actual transport and replace it with ts estimate at "after" time step 
    121          zue(:,:) = pe3u(:,:,1,Kaa) * puu(:,:,1,Kaa) * umask(:,:,1) 
    122          zve(:,:) = pe3v(:,:,1,Kaa) * pvv(:,:,1,Kaa) * vmask(:,:,1) 
    123          DO jk = 2, jpkm1 
    124             zue(:,:) = zue(:,:) + pe3u(:,:,jk,Kaa) * puu(:,:,jk,Kaa) * umask(:,:,jk) 
    125             zve(:,:) = zve(:,:) + pe3v(:,:,jk,Kaa) * pvv(:,:,jk,Kaa) * vmask(:,:,jk) 
    126          END DO 
    127          DO jk = 1, jpkm1 
    128             puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kaa) - zue(:,:) * r1_hu(:,:,Kaa) + uu_b(:,:,Kaa) ) * umask(:,:,jk) 
    129             pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kaa) - zve(:,:) * r1_hv(:,:,Kaa) + vv_b(:,:,Kaa) ) * vmask(:,:,jk) 
    130          END DO 
    131          ! 
    132          IF( .NOT.ln_bt_fw ) THEN 
    133             ! Remove advective velocity from "now velocities" 
    134             ! prior to asselin filtering 
    135             ! In the forward case, this is done below after asselin filtering 
    136             ! so that asselin contribution is removed at the same time 
    137             DO jk = 1, jpkm1 
    138                puu(:,:,jk,Kmm) = ( puu(:,:,jk,Kmm) - un_adv(:,:)*r1_hu(:,:,Kmm) + uu_b(:,:,Kmm) )*umask(:,:,jk) 
    139                pvv(:,:,jk,Kmm) = ( pvv(:,:,jk,Kmm) - vn_adv(:,:)*r1_hv(:,:,Kmm) + vv_b(:,:,Kmm) )*vmask(:,:,jk) 
    140             END DO 
    141          ENDIF 
    142       ENDIF 
    143  
    144       ! Update after velocity on domain lateral boundaries 
    145       ! -------------------------------------------------- 
    146 # if defined key_agrif 
    147       CALL Agrif_dyn( kt )             !* AGRIF zoom boundaries 
    148 # endif 
    149       ! 
    150       CALL lbc_lnk_multi( 'dynatfLF', puu(:,:,:,Kaa), 'U', -1., pvv(:,:,:,Kaa), 'V', -1. )     !* local domain boundaries 
    151       ! 
    152       !                                !* BDY open boundaries 
    153       IF( ln_bdy .AND. ln_dynspg_exp )   CALL bdy_dyn( kt, Kbb, puu, pvv, Kaa ) 
    154       IF( ln_bdy .AND. ln_dynspg_ts  )   CALL bdy_dyn( kt, Kbb, puu, pvv, Kaa, dyn3d_only=.true. ) 
     118!       IF ( ln_dynspg_ts ) THEN 
     119!          ! Ensure below that barotropic velocities match time splitting estimate 
     120!          ! Compute actual transport and replace it with ts estimate at "after" time step 
     121!          zue(:,:) = pe3u(:,:,1,Kaa) * puu(:,:,1,Kaa) * umask(:,:,1) 
     122!          zve(:,:) = pe3v(:,:,1,Kaa) * pvv(:,:,1,Kaa) * vmask(:,:,1) 
     123!          DO jk = 2, jpkm1 
     124!             zue(:,:) = zue(:,:) + pe3u(:,:,jk,Kaa) * puu(:,:,jk,Kaa) * umask(:,:,jk) 
     125!             zve(:,:) = zve(:,:) + pe3v(:,:,jk,Kaa) * pvv(:,:,jk,Kaa) * vmask(:,:,jk) 
     126!          END DO 
     127!          DO jk = 1, jpkm1 
     128!             puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kaa) - zue(:,:) * r1_hu(:,:,Kaa) + uu_b(:,:,Kaa) ) * umask(:,:,jk) 
     129!             pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kaa) - zve(:,:) * r1_hv(:,:,Kaa) + vv_b(:,:,Kaa) ) * vmask(:,:,jk) 
     130!          END DO 
     131!          ! 
     132!          IF( .NOT.ln_bt_fw ) THEN 
     133!             ! Remove advective velocity from "now velocities" 
     134!             ! prior to asselin filtering 
     135!             ! In the forward case, this is done below after asselin filtering 
     136!             ! so that asselin contribution is removed at the same time 
     137!             DO jk = 1, jpkm1 
     138!                puu(:,:,jk,Kmm) = ( puu(:,:,jk,Kmm) - un_adv(:,:)*r1_hu(:,:,Kmm) + uu_b(:,:,Kmm) )*umask(:,:,jk) 
     139!                pvv(:,:,jk,Kmm) = ( pvv(:,:,jk,Kmm) - vn_adv(:,:)*r1_hv(:,:,Kmm) + vv_b(:,:,Kmm) )*vmask(:,:,jk) 
     140!             END DO 
     141!          ENDIF 
     142!       ENDIF 
     143! 
     144!       ! Update after velocity on domain lateral boundaries 
     145!       ! -------------------------------------------------- 
     146! # if defined key_agrif 
     147!       CALL Agrif_dyn( kt )             !* AGRIF zoom boundaries 
     148! # endif 
     149!       ! 
     150!       CALL lbc_lnk_multi( 'dynatfLF', puu(:,:,:,Kaa), 'U', -1., pvv(:,:,:,Kaa), 'V', -1. )     !* local domain boundaries 
     151!       ! 
     152!       !                                !* BDY open boundaries 
     153!       IF( ln_bdy .AND. ln_dynspg_exp )   CALL bdy_dyn( kt, Kbb, puu, pvv, Kaa ) 
     154!       IF( ln_bdy .AND. ln_dynspg_ts  )   CALL bdy_dyn( kt, Kbb, puu, pvv, Kaa, dyn3d_only=.true. ) 
    155155 
    156156!!$   Do we need a call to bdy_vol here?? 
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/SBC/sbcice_cice.F90

    r12377 r12583  
    3636# if defined key_cice4 
    3737   USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  & 
    38                 strocnxT,strocnyT,                               &  
     38                strocnxT,strocnyT,                               & 
    3939                sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_gbm,     & 
    4040                fresh_gbm,fhocn_gbm,fswthru_gbm,frzmlt,          & 
     
    4545#else 
    4646   USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  & 
    47                 strocnxT,strocnyT,                               &  
     47                strocnxT,strocnyT,                               & 
    4848                sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_ai,     & 
    4949                fresh_ai,fhocn_ai,fswthru_ai,frzmlt,          & 
     
    7070   INTEGER             ::   jj_off 
    7171 
    72    INTEGER , PARAMETER ::   jpfld   = 13   ! maximum number of files to read  
     72   INTEGER , PARAMETER ::   jpfld   = 13   ! maximum number of files to read 
    7373   INTEGER , PARAMETER ::   jp_snow = 1    ! index of snow file 
    7474   INTEGER , PARAMETER ::   jp_rain = 2    ! index of rain file 
     
    109109      !!--------------------------------------------------------------------- 
    110110      !!                  ***  ROUTINE sbc_ice_cice  *** 
    111       !!                    
    112       !! ** Purpose :   update the ocean surface boundary condition via the  
    113       !!                CICE Sea Ice Model time stepping  
    114       !! 
    115       !! ** Method  : - Get any extra forcing fields for CICE   
     111      !! 
     112      !! ** Purpose :   update the ocean surface boundary condition via the 
     113      !!                CICE Sea Ice Model time stepping 
     114      !! 
     115      !! ** Method  : - Get any extra forcing fields for CICE 
    116116      !!              - Prepare forcing fields 
    117117      !!              - CICE model time stepping 
    118       !!              - call the routine that computes mass and  
     118      !!              - call the routine that computes mass and 
    119119      !!                heat fluxes at the ice/ocean interface 
    120120      !! 
     
    171171      ! there is no restart file. 
    172172      ! Values from a CICE restart file would overwrite this 
    173       IF( .NOT. ln_rstart ) THEN     
    174          CALL nemo2cice( ts(:,:,1,jp_tem,Kmm) , sst , 'T' , 1.)  
    175       ENDIF   
     173      IF( .NOT. ln_rstart ) THEN 
     174         CALL nemo2cice( ts(:,:,1,jp_tem,Kmm) , sst , 'T' , 1.) 
     175      ENDIF 
    176176#endif 
    177177 
     
    233233!!gm This should be put elsewhere....   (same remark for limsbc) 
    234234!!gm especially here it is assumed zstar coordinate, but it can be ztilde.... 
    235             IF( .NOT.ln_linssh ) THEN 
    236                ! 
    237                DO jk = 1,jpkm1                     ! adjust initial vertical scale factors 
    238                   e3t(:,:,jk,Kmm) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kmm)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
    239                   e3t(:,:,jk,Kbb) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kbb)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
    240                ENDDO 
    241                e3t(:,:,:,Krhs) = e3t(:,:,:,Kbb) 
    242                ! Reconstruction of all vertical scale factors at now and before time-steps 
    243                ! ============================================================================= 
    244                ! Horizontal scale factor interpolations 
    245                ! -------------------------------------- 
    246                CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 
    247                CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 
    248                CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
    249                CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
    250                CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 
    251                ! Vertical scale factor interpolations 
    252                ! ------------------------------------ 
    253                CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  ) 
    254                CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
    255                CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
    256                CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
    257                CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    258                ! t- and w- points depth 
    259                ! ---------------------- 
    260                gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
    261                gdepw(:,:,1,Kmm) = 0.0_wp 
    262                gde3w(:,:,1)     = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    263                DO jk = 2, jpk 
    264                   gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk,Kmm) 
    265                   gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 
    266                   gde3w(:,:,jk)     = gdept(:,:,jk  ,Kmm) - sshn   (:,:) 
    267                END DO 
    268             ENDIF 
     235            IF( .NOT.ln_linssh )   CALL dom_vvl_zgr( Kbb, Kmm, Kaa )   ! interpolation scale factor, depth and water column 
     236            ! IF( .NOT.ln_linssh ) THEN 
     237            !    ! 
     238            !    DO jk = 1,jpkm1                     ! adjust initial vertical scale factors 
     239            !       e3t(:,:,jk,Kmm) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kmm)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
     240            !       e3t(:,:,jk,Kbb) = e3t_0(:,:,jk)*( 1._wp + ssh(:,:,Kbb)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) ) 
     241            !    ENDDO 
     242            !    e3t(:,:,:,Krhs) = e3t(:,:,:,Kbb) 
     243            !    ! Reconstruction of all vertical scale factors at now and before time-steps 
     244            !    ! ============================================================================= 
     245            !    ! Horizontal scale factor interpolations 
     246            !    ! -------------------------------------- 
     247            !    CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' ) 
     248            !    CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' ) 
     249            !    CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
     250            !    CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
     251            !    CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' ) 
     252            !    ! Vertical scale factor interpolations 
     253            !    ! ------------------------------------ 
     254            !    CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  ) 
     255            !    CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
     256            !    CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
     257            !    CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     258            !    CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
     259            !    ! t- and w- points depth 
     260            !    ! ---------------------- 
     261            !    gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
     262            !    gdepw(:,:,1,Kmm) = 0.0_wp 
     263            !    gde3w(:,:,1)     = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
     264            !    DO jk = 2, jpk 
     265            !       gdept(:,:,jk,Kmm) = gdept(:,:,jk-1,Kmm) + e3w(:,:,jk,Kmm) 
     266            !       gdepw(:,:,jk,Kmm) = gdepw(:,:,jk-1,Kmm) + e3t(:,:,jk-1,Kmm) 
     267            !       gde3w(:,:,jk)     = gdept(:,:,jk  ,Kmm) - sshn   (:,:) 
     268            !    END DO 
     269            ! ENDIF 
    269270         ENDIF 
    270271      ENDIF 
     
    272273   END SUBROUTINE cice_sbc_init 
    273274 
    274     
     275 
    275276   SUBROUTINE cice_sbc_in( kt, ksbc ) 
    276277      !!--------------------------------------------------------------------- 
     
    281282      INTEGER, INTENT(in   ) ::   ksbc ! surface forcing type 
    282283      ! 
    283       INTEGER  ::   ji, jj, jl                   ! dummy loop indices       
     284      INTEGER  ::   ji, jj, jl                   ! dummy loop indices 
    284285      REAL(wp), DIMENSION(jpi,jpj) :: ztmp, zpice 
    285286      REAL(wp), DIMENSION(jpi,jpj,ncat) :: ztmpn 
     
    293294      ztmp(:,:)=0.0 
    294295 
    295 ! Aggregate ice concentration already set in cice_sbc_out (or cice_sbc_init on  
     296! Aggregate ice concentration already set in cice_sbc_out (or cice_sbc_init on 
    296297! the first time-step) 
    297298 
    298 ! forced and coupled case  
     299! forced and coupled case 
    299300 
    300301      IF( (ksbc == jp_flx).OR.(ksbc == jp_purecpl) ) THEN 
     
    356357!  Convert to GBM 
    357358            IF(ksbc == jp_flx) THEN 
    358                ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl))*a_i(:,:,jl)  
     359               ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl))*a_i(:,:,jl) 
    359360            ELSE 
    360361               ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl)) 
     
    380381         CALL nemo2cice(ztmp,Tair,'T', 1. )    ! Air temperature (K) 
    381382         CALL nemo2cice(ztmp,potT,'T', 1. )    ! Potential temp (K) 
    382 ! Following line uses MAX(....) to avoid problems if tatm_ice has unset halo rows   
    383          ztmp(:,:) = 101000. / ( 287.04 * MAX(1.0,tatm_ice(:,:)) )     
     383! Following line uses MAX(....) to avoid problems if tatm_ice has unset halo rows 
     384         ztmp(:,:) = 101000. / ( 287.04 * MAX(1.0,tatm_ice(:,:)) ) 
    384385                                               ! Constant (101000.) atm pressure assumed 
    385386         CALL nemo2cice(ztmp,rhoa,'T', 1. )    ! Air density (kg/m^3) 
     
    389390         CALL nemo2cice(ztmp,zlvl,'T', 1. )    ! Atmos level height (m) 
    390391 
    391 ! May want to check all values are physically realistic (as in CICE routine  
     392! May want to check all values are physically realistic (as in CICE routine 
    392393! prepare_forcing)? 
    393394 
    394395! Divide shortwave into spectral bands (as in prepare_forcing) 
    395396         ztmp(:,:)=qsr_ice(:,:,1)*frcvdr       ! visible direct 
    396          CALL nemo2cice(ztmp,swvdr,'T', 1. )              
     397         CALL nemo2cice(ztmp,swvdr,'T', 1. ) 
    397398         ztmp(:,:)=qsr_ice(:,:,1)*frcvdf       ! visible diffuse 
    398          CALL nemo2cice(ztmp,swvdf,'T', 1. )               
     399         CALL nemo2cice(ztmp,swvdf,'T', 1. ) 
    399400         ztmp(:,:)=qsr_ice(:,:,1)*frcidr       ! near IR direct 
    400401         CALL nemo2cice(ztmp,swidr,'T', 1. ) 
     
    406407! Snowfall 
    407408! Ensure fsnow is positive (as in CICE routine prepare_forcing) 
    408       IF( iom_use('snowpre') )   CALL iom_put('snowpre',MAX( (1.0-fr_i(:,:))*sprecip(:,:) ,0.0)) !!Joakim edit   
    409       ztmp(:,:)=MAX(fr_i(:,:)*sprecip(:,:),0.0)   
    410       CALL nemo2cice(ztmp,fsnow,'T', 1. )  
     409      IF( iom_use('snowpre') )   CALL iom_put('snowpre',MAX( (1.0-fr_i(:,:))*sprecip(:,:) ,0.0)) !!Joakim edit 
     410      ztmp(:,:)=MAX(fr_i(:,:)*sprecip(:,:),0.0) 
     411      CALL nemo2cice(ztmp,fsnow,'T', 1. ) 
    411412 
    412413! Rainfall 
    413414      IF( iom_use('precip') )   CALL iom_put('precip', (1.0-fr_i(:,:))*(tprecip(:,:)-sprecip(:,:)) ) !!Joakim edit 
    414415      ztmp(:,:)=fr_i(:,:)*(tprecip(:,:)-sprecip(:,:)) 
    415       CALL nemo2cice(ztmp,frain,'T', 1. )  
     416      CALL nemo2cice(ztmp,frain,'T', 1. ) 
    416417 
    417418! Freezing/melting potential 
     
    482483      INTEGER, INTENT( in  ) ::   kt   ! ocean time step 
    483484      INTEGER, INTENT( in  ) ::   ksbc ! surface forcing type 
    484        
     485 
    485486      INTEGER  ::   ji, jj, jl                 ! dummy loop indices 
    486487      REAL(wp), DIMENSION(jpi,jpj) :: ztmp1, ztmp2 
     
    490491         IF(lwp) WRITE(numout,*)'cice_sbc_out' 
    491492      ENDIF 
    492        
    493 ! x comp of ocean-ice stress  
     493 
     494! x comp of ocean-ice stress 
    494495      CALL cice2nemo(strocnx,ztmp1,'F', -1. ) 
    495496      ss_iou(:,:)=0.0 
     
    500501      CALL lbc_lnk( 'sbcice_cice', ss_iou , 'U', -1. ) 
    501502 
    502 ! y comp of ocean-ice stress  
     503! y comp of ocean-ice stress 
    503504      CALL cice2nemo(strocny,ztmp1,'F', -1. ) 
    504505      ss_iov(:,:)=0.0 
     
    513514! Combine wind stress and ocean-ice stress 
    514515! [Note that fr_iu hasn't yet been updated, so still from start of CICE timestep] 
    515 ! strocnx and strocny already weighted by ice fraction in CICE so not done here  
     516! strocnx and strocny already weighted by ice fraction in CICE so not done here 
    516517 
    517518      utau(:,:)=(1.0-fr_iu(:,:))*utau(:,:)-ss_iou(:,:) 
    518       vtau(:,:)=(1.0-fr_iv(:,:))*vtau(:,:)-ss_iov(:,:)      
    519   
    520 ! Also need ice/ocean stress on T points so that taum can be updated  
    521 ! This interpolation is already done in CICE so best to use those values  
    522       CALL cice2nemo(strocnxT,ztmp1,'T',-1.)  
    523       CALL cice2nemo(strocnyT,ztmp2,'T',-1.)  
    524   
    525 ! Update taum with modulus of ice-ocean stress  
    526 ! strocnxT and strocnyT are not weighted by ice fraction in CICE so must be done here  
    527 taum(:,:)=(1.0-fr_i(:,:))*taum(:,:)+fr_i(:,:)*SQRT(ztmp1*ztmp1 + ztmp2*ztmp2)  
    528  
    529 ! Freshwater fluxes  
     519      vtau(:,:)=(1.0-fr_iv(:,:))*vtau(:,:)-ss_iov(:,:) 
     520 
     521! Also need ice/ocean stress on T points so that taum can be updated 
     522! This interpolation is already done in CICE so best to use those values 
     523      CALL cice2nemo(strocnxT,ztmp1,'T',-1.) 
     524      CALL cice2nemo(strocnyT,ztmp2,'T',-1.) 
     525 
     526! Update taum with modulus of ice-ocean stress 
     527! strocnxT and strocnyT are not weighted by ice fraction in CICE so must be done here 
     528taum(:,:)=(1.0-fr_i(:,:))*taum(:,:)+fr_i(:,:)*SQRT(ztmp1*ztmp1 + ztmp2*ztmp2) 
     529 
     530! Freshwater fluxes 
    530531 
    531532      IF(ksbc == jp_flx) THEN 
    532533! Note that emp from the forcing files is evap*(1-aice)-(tprecip-aice*sprecip) 
    533534! What we want here is evap*(1-aice)-tprecip*(1-aice) hence manipulation below 
    534 ! Not ideal since aice won't be the same as in the atmosphere.   
     535! Not ideal since aice won't be the same as in the atmosphere. 
    535536! Better to use evap and tprecip? (but for now don't read in evap in this case) 
    536537         emp(:,:)  = emp(:,:)+fr_i(:,:)*(tprecip(:,:)-sprecip(:,:)) 
    537538      ELSE IF(ksbc == jp_blk) THEN 
    538          emp(:,:)  = (1.0-fr_i(:,:))*emp(:,:)         
     539         emp(:,:)  = (1.0-fr_i(:,:))*emp(:,:) 
    539540      ELSE IF(ksbc == jp_purecpl) THEN 
    540 ! emp_tot is set in sbc_cpl_ice_flx (called from cice_sbc_in above)  
     541! emp_tot is set in sbc_cpl_ice_flx (called from cice_sbc_in above) 
    541542! This is currently as required with the coupling fields from the UM atmosphere 
    542          emp(:,:) = emp_tot(:,:)+tprecip(:,:)*fr_i(:,:)  
     543         emp(:,:) = emp_tot(:,:)+tprecip(:,:)*fr_i(:,:) 
    543544      ENDIF 
    544545 
     
    560561      emp(:,:)=emp(:,:)-ztmp1(:,:) 
    561562      fmmflx(:,:) = ztmp1(:,:) !!Joakim edit 
    562        
     563 
    563564      CALL lbc_lnk_multi( 'sbcice_cice', emp , 'T', 1., sfx , 'T', 1. ) 
    564565 
     
    634635      !!                    ***  ROUTINE cice_sbc_hadgam  *** 
    635636      !! ** Purpose: Prepare fields needed to pass to HadGAM3 atmosphere 
    636       !!  
     637      !! 
    637638      !! 
    638639      !!--------------------------------------------------------------------- 
     
    657658      CALL cice2nemo(vvel,v_ice,'F', -1. ) 
    658659      ! 
    659       ! Ice concentration (CO_1) = a_i calculated at end of cice_sbc_out   
     660      ! Ice concentration (CO_1) = a_i calculated at end of cice_sbc_out 
    660661      ! 
    661662      ! Snow and ice thicknesses (CO_2 and CO_3) 
     
    689690      !!--------------------------------------------------------------------- 
    690691      !! ** Method  :   READ monthly flux file in NetCDF files 
    691       !!       
    692       !!  snowfall     
    693       !!  rainfall     
    694       !!  sublimation rate     
     692      !! 
     693      !!  snowfall 
     694      !!  rainfall 
     695      !!  sublimation rate 
    695696      !!  topmelt (category) 
    696697      !!  botmelt (category) 
     
    709710      TYPE(FLD_N) ::   sn_snow, sn_rain, sn_sblm               ! informations about the fields to be read 
    710711      TYPE(FLD_N) ::   sn_top1, sn_top2, sn_top3, sn_top4, sn_top5 
    711       TYPE(FLD_N) ::   sn_bot1, sn_bot2, sn_bot3, sn_bot4, sn_bot5  
     712      TYPE(FLD_N) ::   sn_bot1, sn_bot2, sn_bot3, sn_bot4, sn_bot5 
    712713      !! 
    713714      NAMELIST/namsbc_cice/ cn_dir, sn_snow, sn_rain, sn_sblm,   & 
     
    727728         !            !    file          ! frequency !  variable    ! time intep !  clim   ! 'yearly' or ! weights  ! rotation   ! landmask 
    728729         !            !    name          !  (hours)  !   name       !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs      ! file 
    729          sn_snow = FLD_N( 'snowfall_1m'  ,    -1.    ,  'snowfall'  ,  .true.    , .true.  ,  ' yearly'  , ''       , ''         ,  ''    )  
    730          sn_rain = FLD_N( 'rainfall_1m'  ,    -1.    ,  'rainfall'  ,  .true.    , .true.  ,  ' yearly'  , ''       , ''         ,  ''    )  
     730         sn_snow = FLD_N( 'snowfall_1m'  ,    -1.    ,  'snowfall'  ,  .true.    , .true.  ,  ' yearly'  , ''       , ''         ,  ''    ) 
     731         sn_rain = FLD_N( 'rainfall_1m'  ,    -1.    ,  'rainfall'  ,  .true.    , .true.  ,  ' yearly'  , ''       , ''         ,  ''    ) 
    731732         sn_sblm = FLD_N( 'sublim_1m'    ,    -1.    ,  'sublim'    ,  .true.    , .true.  ,  ' yearly'  , ''       , ''         ,  ''    ) 
    732733         sn_top1 = FLD_N( 'topmeltn1_1m' ,    -1.    ,  'topmeltn1' ,  .true.    , .true.  ,  ' yearly'  , ''       , ''         ,  ''    ) 
     
    754755         slf_i(jp_bot2) = sn_bot2   ;   slf_i(jp_bot3) = sn_bot3   ;   slf_i(jp_bot4) = sn_bot4 
    755756         slf_i(jp_bot5) = sn_bot5 
    756           
     757 
    757758         ! set sf structure 
    758759         ALLOCATE( sf(jpfld), STAT=ierror ) 
     
    792793      ! control print (if less than 100 time-step asked) 
    793794      IF( nitend-nit000 <= 100 .AND. lwp ) THEN 
    794          WRITE(numout,*)  
     795         WRITE(numout,*) 
    795796         WRITE(numout,*) '        read forcing fluxes for CICE OK' 
    796797         CALL FLUSH(numout) 
     
    802803      !!--------------------------------------------------------------------- 
    803804      !!                    ***  ROUTINE nemo2cice  *** 
    804       !! ** Purpose :   Transfer field in NEMO array to field in CICE array.   
     805      !! ** Purpose :   Transfer field in NEMO array to field in CICE array. 
    805806#if defined key_nemocice_decomp 
    806       !!              
     807      !! 
    807808      !!                NEMO and CICE PE sub domains are identical, hence 
    808       !!                there is no need to gather or scatter data from  
     809      !!                there is no need to gather or scatter data from 
    809810      !!                one PE configuration to another. 
    810811#else 
    811       !!                Automatically gather/scatter between  
     812      !!                Automatically gather/scatter between 
    812813      !!                different processors and blocks 
    813814      !! ** Method :    A. Ensure all haloes are filled in NEMO field (pn) 
    814815      !!                B. Gather pn into global array (png) 
    815816      !!                C. Map png into CICE global array (pcg) 
    816       !!                D. Scatter pcg to CICE blocks (pc) + update haloes   
     817      !!                D. Scatter pcg to CICE blocks (pc) + update haloes 
    817818#endif 
    818819      !!--------------------------------------------------------------------- 
     
    858859      IF( jpnij > 1) THEN 
    859860         CALL mppsync 
    860          CALL mppgather (pn,0,png)  
     861         CALL mppgather (pn,0,png) 
    861862         CALL mppsync 
    862863      ELSE 
     
    869870! (may be OK but not 100% sure) 
    870871 
    871       IF(nproc==0) THEN      
     872      IF(nproc==0) THEN 
    872873!        pcg(:,:)=0.0 
    873874         DO jn=1,jpnij 
     
    890891         CASE ( 'T' ) 
    891892            grid_loc=field_loc_center 
    892          CASE ( 'F' )                               
     893         CASE ( 'F' ) 
    893894            grid_loc=field_loc_NEcorner 
    894895      END SELECT 
     
    897898         CASE ( -1 ) 
    898899            field_type=field_type_vector 
    899          CASE ( 1 )                               
     900         CASE ( 1 ) 
    900901            field_type=field_type_scalar 
    901902      END SELECT 
     
    916917      !! ** Purpose :   Transfer field in CICE array to field in NEMO array. 
    917918#if defined key_nemocice_decomp 
    918       !!              
     919      !! 
    919920      !!                NEMO and CICE PE sub domains are identical, hence 
    920       !!                there is no need to gather or scatter data from  
     921      !!                there is no need to gather or scatter data from 
    921922      !!                one PE configuration to another. 
    922 #else  
     923#else 
    923924      !!                Automatically deal with scatter/gather between 
    924925      !!                different processors and blocks 
     
    926927      !!                B. Map pcg into NEMO global array (png) 
    927928      !!                C. Scatter png into NEMO field (pn) for each processor 
    928       !!                D. Ensure all haloes are filled in pn  
     929      !!                D. Ensure all haloes are filled in pn 
    929930#endif 
    930931      !!--------------------------------------------------------------------- 
     
    958959         CASE ( 'T' ) 
    959960            grid_loc=field_loc_center 
    960          CASE ( 'F' )                               
     961         CASE ( 'F' ) 
    961962            grid_loc=field_loc_NEcorner 
    962963      END SELECT 
     
    965966         CASE ( -1 ) 
    966967            field_type=field_type_vector 
    967          CASE ( 1 )                               
     968         CASE ( 1 ) 
    968969            field_type=field_type_scalar 
    969970      END SELECT 
     
    979980#else 
    980981 
    981 !      A. Gather CICE blocks (pc) into global array (pcg)  
     982!      A. Gather CICE blocks (pc) into global array (pcg) 
    982983 
    983984      CALL gather_global(pcg, pc, 0, distrb_info) 
     
    10051006      IF( jpnij > 1) THEN 
    10061007         CALL mppsync 
    1007          CALL mppscatter (png,0,pn)  
     1008         CALL mppscatter (png,0,pn) 
    10081009         CALL mppsync 
    10091010      ELSE 
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/TRA/traatfLF.F90

    r12581 r12583  
    103103      ENDIF 
    104104 
    105       ! Update after tracer on domain lateral boundaries 
    106       ! 
    107 #if defined key_agrif 
    108       CALL Agrif_tra                     ! AGRIF zoom boundaries 
    109 #endif 
    110       !                                              ! local domain boundaries  (T-point, unchanged sign) 
    111       CALL lbc_lnk_multi( 'traatfLF', pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1. ) 
    112       ! 
    113       IF( ln_bdy )   CALL bdy_tra( kt, Kbb, pts, Kaa )  ! BDY open boundaries 
     105!       ! Update after tracer on domain lateral boundaries 
     106!       ! 
     107! #if defined key_agrif 
     108!       CALL Agrif_tra                     ! AGRIF zoom boundaries 
     109! #endif 
     110!       !                                              ! local domain boundaries  (T-point, unchanged sign) 
     111!       CALL lbc_lnk_multi( 'traatfLF', pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1. ) 
     112!       ! 
     113!       IF( ln_bdy )   CALL bdy_tra( kt, Kbb, pts, Kaa )  ! BDY open boundaries 
    114114 
    115115      ! set time step size (Euler/Leapfrog) 
  • NEMO/branches/2020/dev_r12377_KERNEL-06_techene_e3/src/OCE/steplf.F90

    r12581 r12583  
    4343   USE traatfLF          ! time filtering                   (tra_atf_lf routine) 
    4444   USE dynatfLF          ! time filtering                   (dyn_atf_lf routine) 
     45   USE dynspg_ts         ! surface pressure gradient: split-explicit scheme (define un_adv) 
     46   USE bdydyn            ! ocean open boundary conditions (define bdy_dyn) 
    4547 
    4648   IMPLICIT NONE 
     
    288290!! 
    289291!!jc2: dynnxt must be the latest call. e3t(:,:,:,Nbb) are indeed updated in that routine 
     292                         CALL zdyn_ts       ( Nnn, Naa, e3u, e3v, uu, vv )                   ! barotrope ajustment 
     293                         CALL finalize_sbc  ( kstp, Nbb, Naa, uu, vv, ts )                   ! boundary condifions 
    290294                         CALL ssh_atf       ( kstp, Nbb, Nnn, Naa, ssh )                     ! time filtering of "now" sea surface height 
    291295                         CALL dom_qe_r3c    ( ssh(:,:,Nnn), r3t_f, r3u_f, r3v_f )            ! "now" ssh/h_0 ratio from filtrered ssh 
     
    352356      ! 
    353357   END SUBROUTINE stplf 
     358 
     359 
     360   SUBROUTINE zdyn_ts (Kmm, Kaa, pe3u, pe3v, puu, pvv) 
     361      !!---------------------------------------------------------------------- 
     362      !!                  ***  ROUTINE zdyn_ts  *** 
     363      !! 
     364      !! ** Purpose :   Finalize after horizontal velocity. 
     365      !! 
     366      !! ** Method  : * Ensure after velocities transport matches time splitting 
     367      !!             estimate (ln_dynspg_ts=T) 
     368      !! 
     369      !! ** Action :   puu(Kmm),pvv(Kmm),puu(Kaa),pvv(Kaa)   now and after horizontal velocity 
     370      !!---------------------------------------------------------------------- 
     371      INTEGER                             , INTENT(in   ) :: Kmm, Kaa    ! before and after time level indices 
     372      REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv         ! velocities 
     373      REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(in   ) :: pe3u, pe3v   ! scale factors 
     374      ! 
     375      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   zue, zve 
     376      ! 
     377      INTEGER  ::   jk   ! dummy loop indices 
     378      !!---------------------------------------------------------------------- 
     379 
     380      IF ( ln_dynspg_ts ) THEN 
     381         ALLOCATE( zue(jpi,jpj)     , zve(jpi,jpj)     ) 
     382         ! Ensure below that barotropic velocities match time splitting estimate 
     383         ! Compute actual transport and replace it with ts estimate at "after" time step 
     384         zue(:,:) = pe3u(:,:,1,Kaa) * puu(:,:,1,Kaa) * umask(:,:,1) 
     385         zve(:,:) = pe3v(:,:,1,Kaa) * pvv(:,:,1,Kaa) * vmask(:,:,1) 
     386         DO jk = 2, jpkm1 
     387            zue(:,:) = zue(:,:) + pe3u(:,:,jk,Kaa) * puu(:,:,jk,Kaa) * umask(:,:,jk) 
     388            zve(:,:) = zve(:,:) + pe3v(:,:,jk,Kaa) * pvv(:,:,jk,Kaa) * vmask(:,:,jk) 
     389         END DO 
     390         DO jk = 1, jpkm1 
     391            puu(:,:,jk,Kaa) = ( puu(:,:,jk,Kaa) - zue(:,:) * r1_hu(:,:,Kaa) + uu_b(:,:,Kaa) ) * umask(:,:,jk) 
     392            pvv(:,:,jk,Kaa) = ( pvv(:,:,jk,Kaa) - zve(:,:) * r1_hv(:,:,Kaa) + vv_b(:,:,Kaa) ) * vmask(:,:,jk) 
     393         END DO 
     394         ! 
     395         IF( .NOT.ln_bt_fw ) THEN 
     396            ! Remove advective velocity from "now velocities" 
     397            ! prior to asselin filtering 
     398            ! In the forward case, this is done below after asselin filtering 
     399            ! so that asselin contribution is removed at the same time 
     400            DO jk = 1, jpkm1 
     401               puu(:,:,jk,Kmm) = ( puu(:,:,jk,Kmm) - un_adv(:,:)*r1_hu(:,:,Kmm) + uu_b(:,:,Kmm) )*umask(:,:,jk) 
     402               pvv(:,:,jk,Kmm) = ( pvv(:,:,jk,Kmm) - vn_adv(:,:)*r1_hv(:,:,Kmm) + vv_b(:,:,Kmm) )*vmask(:,:,jk) 
     403            END DO 
     404         ENDIF 
     405         ! 
     406         DEALLOCATE( zue, zve ) 
     407         ! 
     408      ENDIF 
     409      ! 
     410   END SUBROUTINE zdyn_ts 
     411 
     412 
     413   SUBROUTINE finalize_sbc (kt, Kbb, Kaa, puu, pvv, pts) 
     414      !!---------------------------------------------------------------------- 
     415      !!                  ***  ROUTINE finalize_sbc  *** 
     416      !! 
     417      !! ** Purpose :   Apply the boundary condition on the after velocity 
     418      !! 
     419      !! ** Method  : * Apply lateral boundary conditions on after velocity 
     420      !!             at the local domain boundaries through lbc_lnk call, 
     421      !!             at the one-way open boundaries (ln_bdy=T), 
     422      !!             at the AGRIF zoom   boundaries (lk_agrif=T) 
     423      !! 
     424      !! ** Action :   puu(Kaa),pvv(Kaa)   after horizontal velocity and tracers 
     425      !!---------------------------------------------------------------------- 
     426      INTEGER                             , INTENT(in   ) :: kt               ! ocean time-step index 
     427      INTEGER                             , INTENT(in   ) :: Kbb, Kaa    ! before and after time level indices 
     428      REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv         ! velocities to be time filtered 
     429      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts            ! active tracers 
     430      ! 
     431      ! Update after tracer and velocity on domain lateral boundaries 
     432      ! 
     433#if defined key_agrif 
     434            CALL Agrif_tra                     ! AGRIF zoom boundaries 
     435            CALL Agrif_dyn( kt )               !* AGRIF zoom boundaries 
     436#endif 
     437      !                                        ! local domain boundaries  (T-point, unchanged sign) 
     438      CALL lbc_lnk_multi( 'finalize_sbc', puu(:,:,:,       Kaa), 'U', -1., pvv(:,:,:       ,Kaa), 'V', -1.   & 
     439                       &                , pts(:,:,:,jp_tem,Kaa), 'T',  1., pts(:,:,:,jp_sal,Kaa), 'T',  1. )     !* local domain boundaries 
     440      ! 
     441      !                                        !* BDY open boundaries 
     442      IF( ln_bdy )   THEN 
     443                               CALL bdy_tra( kt, Kbb, pts,      Kaa ) 
     444         IF( ln_dynspg_exp )   CALL bdy_dyn( kt, Kbb, puu, pvv, Kaa ) 
     445         IF( ln_dynspg_ts  )   CALL bdy_dyn( kt, Kbb, puu, pvv, Kaa, dyn3d_only=.true. ) 
     446      ENDIF 
     447      ! 
     448   END SUBROUTINE finalize_sbc 
     449 
    354450   ! 
    355451   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.