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 11931 – NEMO

Changeset 11931


Ignore:
Timestamp:
2019-11-19T18:30:57+01:00 (5 years ago)
Author:
mathiot
Message:

ENHANCE-02_ISF_nemo: add comments, improve memory usage of ln_isfcpl_cons option, fix issue in ISOMIP+ configuration

Location:
NEMO/branches/2019/ENHANCE-02_ISF_nemo
Files:
21 edited
1 moved

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynnxt.F90

    r11541 r11931  
    2626   !!------------------------------------------------------------------------- 
    2727   USE oce            ! ocean dynamics and tracers 
    28    USE isf 
    2928   USE dom_oce        ! ocean space and time domain 
    3029   USE sbc_oce        ! Surface boundary condition: ocean fields 
    3130   USE sbcrnf         ! river runoffs 
    32    USE isfnxt 
    3331   USE phycst         ! physical constants 
    3432   USE dynadv         ! dynamics: vector invariant versus flux form 
     
    4240   USE trddyn         ! trend manager: dynamics 
    4341   USE trdken         ! trend manager: kinetic energy 
     42   USE isf       , ONLY: ln_isf     ! ice shelf 
     43   USE isfdynnxt , ONLY: isf_dynnxt ! ice shelf  
    4444   ! 
    4545   USE in_out_manager ! I/O manager 
     
    243243            END IF 
    244244            ! 
    245             ! ice shelf melting 
     245            ! ice shelf melting (deal separatly as it can be in depth) 
     246            ! PM: we could probably define a generic subroutine to do the in depth correction 
     247            !     to manage rnf, isf and possibly in the futur icb, tide water glacier (...) 
    246248            IF ( ln_isf ) CALL isf_dynnxt( kt, atfp * rdt ) 
    247249            ! 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/sshwzv.F90

    r11541 r11931  
    260260               &                             + fwfisf_cav_b(:,:) - fwfisf_cav(:,:)   & 
    261261               &                             + fwfisf_par_b(:,:) - fwfisf_par(:,:)   ) * ssmask(:,:) 
     262 
     263            ! ice sheet coupling 
     264            IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) sshb(:,:) = sshb(:,:) - atfp * rdt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 
     265 
    262266         ENDIF 
    263  
    264          ! ice sheet coupling 
    265          IF ( ln_isf .AND. ln_isfcpl .AND. kt == nit000+1) sshb(:,:) = sshb(:,:) - atfp * rdt * ( risfcpl_ssh(:,:) - 0.0 ) * ssmask(:,:) 
    266267 
    267268         sshn(:,:) = ssha(:,:)                              ! now <-- after 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isf.F90

    r11876 r11931  
    1414   !!---------------------------------------------------------------------- 
    1515 
    16    USE in_out_manager, ONLY: wp, jpi,jpj, jpk, jpts ! I/O manager 
     16   USE par_oce       , ONLY: jpi, jpj, jpk 
     17   USE in_out_manager, ONLY: wp, jpts ! I/O manager 
    1718   USE lib_mpp       , ONLY: ctl_stop, mpp_sum      ! MPP library 
    1819   USE fldread        ! read input fields 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcav.F90

    r11852 r11931  
    2323   ! 
    2424   USE oce      , ONLY: tsn                    ! ocean tracers 
    25    USE dom_oce  , ONLY: jpi,jpj                ! ocean space and time domain 
     25   USE par_oce  , ONLY: jpi,jpj                ! ocean space and time domain 
    2626   USE phycst   , ONLY: grav,rau0,r1_rau0_rcp  ! physical constants 
    2727   USE eosbn2   , ONLY: l_useCT                ! l_useCT 
     
    9090      DO WHILE ( lit )    ! maybe just a constant number of iteration as in blk_core is fine 
    9191         ! 
    92          ! compute gammat every where (2d) 
     92         ! compute gammat everywhere (2d) 
    9393         ! useless if melt specified 
    9494         IF ( TRIM(cn_isfcav_mlt) .NE. 'spe' ) THEN 
     
    101101            &                         zqhc   , zqoce, pqfwf  ) 
    102102         ! 
    103          ! define if we need to iterate (nn_gammablk 0/1 do not need iteration) 
     103         ! define if we need to iterate 
    104104         SELECT CASE ( cn_gammablk ) 
    105105         CASE ( 'spe','ad15' ) 
     
    123123      END DO 
    124124      ! 
    125       ! compute heat and water flux  (change signe directly in the melt subroutine) 
     125      ! compute heat and water flux ( > 0 out ) 
    126126      pqfwf(:,:) = pqfwf(:,:) * mskisf_cav(:,:) 
    127127      zqoce(:,:) = zqoce(:,:) * mskisf_cav(:,:) 
    128128      zqhc (:,:) = zqhc(:,:)  * mskisf_cav(:,:) 
    129129      ! 
    130       ! compute heat content flux 
    131       zqlat(:,:) = - pqfwf(:,:) * rLfusisf    ! 2d latent heat flux (W/m2) ( > 0 out ) 
     130      ! compute heat content flux ( > 0 out ) 
     131      zqlat(:,:) = - pqfwf(:,:) * rLfusisf    ! 2d latent heat flux (W/m2) 
    132132      ! 
    133133      ! total heat flux ( >0 out ) 
     
    163163      INTEGER :: ierr 
    164164      !!--------------------------------------------------------------------- 
    165  
     165      ! 
     166      !============== 
    166167      ! 0: allocation 
     168      !============== 
    167169      ! 
    168170      CALL isf_alloc_cav() 
    169171      ! 
     172      !================== 
    170173      ! 1: initialisation 
     174      !================== 
    171175      ! 
    172176      ! top and bottom level of the 'top boundary layer' 
     
    179183      mskisf_cav(:,:)    = (1._wp - tmask(:,:,1)) * ssmask(:,:) 
    180184      ! 
     185      !================ 
    181186      ! 2: read restart 
     187      !================ 
    182188      ! 
    183189      ! read cav variable from restart 
    184190      IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b) 
    185191      ! 
     192      !========================================== 
    186193      ! 3: specific allocation and initialisation (depending of scheme choice) 
     194      !========================================== 
    187195      ! 
    188196      SELECT CASE ( TRIM(cn_isfcav_mlt) ) 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavgam.F90

    r11852 r11931  
    5757      !!--------------------------------------------------------------------- 
    5858      ! 
    59       ! compute velocity in the tbl if needed 
     59      !========================================== 
     60      ! 1.: compute velocity in the tbl if needed 
     61      !========================================== 
     62      ! 
    6063      SELECT CASE ( cn_gammablk ) 
    6164      CASE ( 'spe'  )  
     
    7881      END SELECT 
    7982      !  
    80       ! compute gamma 
     83      !========================================== 
     84      ! 2.: compute gamma 
     85      !========================================== 
     86      ! 
    8187      SELECT CASE ( cn_gammablk ) 
    8288      CASE ( 'spe'  ) ! gamma is constant (specified in namelist) 
     
    9197      END SELECT 
    9298      ! 
    93       ! ouput exchange coeficient and tbl velocity 
     99      !========================================== 
     100      ! 3.: output and debug 
     101      !========================================== 
     102      ! 
    94103      CALL iom_put('isfgammat', pgt(:,:)) 
    95104      CALL iom_put('isfgammas', pgs(:,:)) 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavmlt.F90

    r11876 r11931  
    33   !!                       ***  MODULE  isfcavmlt  *** 
    44   !! ice shelf module :  update surface ocean boundary condition under ice 
    5    !!                   shelf 
     5   !!                   shelves 
    66   !!====================================================================== 
    77   !! History :  4.0  !  2019-09  (P. Mathiot) Original code 
     
    99 
    1010   !!---------------------------------------------------------------------- 
    11    !!   isfcav_mlt    : update surface ocean boundary condition under ice shelf 
     11   !!   isfcav_mlt    : compute or read ice shelf fwf/heat fluxes in the ice shelf cavity 
    1212   !!---------------------------------------------------------------------- 
    1313 
     
    2323   USE in_out_manager              ! I/O manager 
    2424   USE iom        , ONLY: iom_put  ! I/O library 
    25    USE fldread    , ONLY: fld_read ! 
     25   USE fldread    , ONLY: fld_read, FLD, FLD_N ! 
    2626   USE lib_fortran, ONLY: glob_sum ! 
    2727   USE lib_mpp    , ONLY: ctl_stop ! 
     
    5151      !! ** Purpose    : compute or read ice shelf fwf/heat fluxes in the ice shelf cavity 
    5252      !! 
    53       !!--------------------------------------------------------------------- 
    5453      !!-------------------------- OUT ------------------------------------- 
    5554      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pqhc, pqoce, pqfwf  ! heat and fwf fluxes 
     
    5857      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pgt  , pgs    ! gamma t and gamma s 
    5958      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pttbl, pstbl  ! top boundary layer tracer 
    60       !!--------------------------------------------------------------------- 
    6159      !!--------------------------------------------------------------------- 
    6260      ! 
     
    110108      !!-------------------------------------------------------------------- 
    111109      ! 
    112       ! Calculate freezing temperature 
     110      ! Compute freezing temperature 
    113111      CALL eos_fzp( pstbl(:,:), ztfrz(:,:), risfdep(:,:) ) 
    114112      ! 
     
    131129      !!---------------------------------------------------------------------- 
    132130      !! 
    133       !!                          ***  ROUTINE isfcav_mlt_spe  *** 
     131      !!                          ***  ROUTINE isfcav_mlt_2eq  *** 
    134132      !! 
    135133      !! ** Purpose    : Compute ice shelf fwf/heqt fluxes using ISOMIP formulation (Hunter et al., 2006) 
     
    144142      !!                 Tech.  Rep.  June,  Antarctic  Climate  &  Ecosystems  Cooperative  Research  Centre,  available  at:   
    145143      !!                 http://staff.acecrc.org.au/~bkgalton/ISOMIP/test_cavities.pdf (last access: 21 July 2016), 2006. 
    146       !!--------------------------------------------------------------------- 
     144      !! 
    147145      !!-------------------------- OUT ------------------------------------- 
    148146      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pqhc, pqoce, pqfwf  ! hean content, ocean-ice heat and fwf fluxes 
     
    192190      !!                   MISMIP v. 3 (MISMIP +), ISOMIP v. 2 (ISOMIP +) and MISOMIP v. 1 (MISOMIP1),  
    193191      !!                   Geosci. Model Dev., 9, 2471-2497, https://doi.org/10.5194/gmd-9-2471-2016, 2016.  
    194       !!--------------------------------------------------------------------- 
     192      !! 
    195193      !!-------------------------- OUT ------------------------------------- 
    196194      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pqhc, pqoce, pqfwf  ! latent heat and fwf fluxes 
     
    308306      ! 
    309307      CALL iom_put('isftfrz_cav', ztfrz * mskisf_cav(:,:) ) 
     308      ! 
    310309   END SUBROUTINE isfcav_mlt_oasis 
    311310 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcpl.F90

    r11908 r11931  
    1212   !!   isfrst : read/write iceshelf variables in/from restart 
    1313   !!---------------------------------------------------------------------- 
    14    USE isf                             ! ice shelf variable 
     14   USE isf                              ! ice shelf variable 
    1515   USE isfutils, ONLY : debug 
    16    USE lib_mpp, ONLY: mpp_sum, mpp_max ! mpp routine 
    17    USE domvvl , ONLY: dom_vvl_zgr      ! vertical scale factor interpolation 
     16   USE lib_mpp , ONLY: mpp_sum, mpp_max ! mpp routine 
     17   USE domvvl  , ONLY: dom_vvl_zgr      ! vertical scale factor interpolation 
    1818   ! 
    1919   USE oce            ! ocean dynamics and tracers 
     
    113113   SUBROUTINE isfcpl_rst_write(kt) 
    114114      !!--------------------------------------------------------------------- 
    115       !!   isfrst_cpl_write : write icesheet coupling variables in restart 
    116       !!--------------------------------------------------------------------- 
    117       !!-------------------------- OUT -------------------------------------- 
     115      !!                   ***  ROUTINE iscpl_rst_write  *** 
     116      !!  
     117      !! ** Purpose : write icesheet coupling variables in restart 
     118      !! 
    118119      !!-------------------------- IN  -------------------------------------- 
    119120      INTEGER, INTENT(in) :: kt 
    120       !!---------------------------------------------------------------------- 
    121121      !!---------------------------------------------------------------------- 
    122122      ! 
     
    136136      !!                   ***  ROUTINE iscpl_ssh  *** 
    137137      !!  
    138       !! ** Purpose :   basic guess of ssh in new wet cell during coupling step 
     138      !! ** Purpose :   basic guess of ssh in new wet cell 
    139139      !!  
    140140      !! ** Method  :   basic extrapolation from neigbourg cells 
    141141      !! 
    142       !!---------------------------------------------------------------------- 
    143       !!-------------------------- OUT -------------------------------------- 
    144       !!-------------------------- IN  -------------------------------------- 
    145142      !!---------------------------------------------------------------------- 
    146143      !! 
     
    344341      tsn(:,:,:,jp_sal) = zts0(:,:,:,jp_sal) * tmask(:,:,:) 
    345342      ! 
    346       tsb(:,:,:,:) = tsn(:,:,:,:) 
    347       ! 
    348343      ! sanity check 
    349344      ! ----------------------------------------------------------------------------------------- 
     
    366361      !!                   ***  ROUTINE iscpl_vol  *** 
    367362      !!  
    368       !! ** Purpose :    
    369       !!                 
    370       !!  
    371       !! ** Method  :    
     363      !! ** Purpose : compute the correction of the local divergence to apply   
     364      !!              during the first time step after the coupling. 
     365      !! 
     366      !! ** Method  : - compute horizontal vol div. before/after coupling 
     367      !!              - compute vertical input 
     368      !!              - compute correction 
    372369      !!                 
    373370      !!---------------------------------------------------------------------- 
     
    376373      INTEGER :: ikb, ikt 
    377374      !! 
    378       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zqvolb, zqvoln 
    379       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztmask_b 
    380       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3u_b, ze3v_b 
    381       !!---------------------------------------------------------------------- 
    382       ! 
    383       CALL iom_get( numror, jpdom_autoglo, 'tmask'  , ztmask_b, ldxios = lrxios   ) ! need to extrapolate T/S 
     375      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zqvolb, zqvoln  ! vol flux div.         before/after coupling 
     376      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3u_b, ze3v_b  ! vertical scale factor before/after coupling 
     377      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztmask_b        ! mask                  before       coupling 
     378      !!---------------------------------------------------------------------- 
     379      ! 
     380      CALL iom_get( numror, jpdom_autoglo, 'tmask'  , ztmask_b, ldxios = lrxios ) 
    384381      CALL iom_get( numror, jpdom_autoglo, 'e3u_n'  , ze3u_b  , ldxios = lrxios ) 
    385382      CALL iom_get( numror, jpdom_autoglo, 'e3v_n'  , ze3v_b  , ldxios = lrxios ) 
    386383      ! 
    387       ! get volume flux before coupling (>0 out) 
     384      ! 1.0: compute horizontal volume flux divergence difference before-after coupling 
     385      ! 
    388386      DO jk = 1, jpk                                 ! Horizontal slab 
     387         ! 1.1: get volume flux before coupling (>0 out) 
    389388         DO jj = 2, jpjm1 
    390389            DO ji = 2, jpim1 
     
    395394         ENDDO 
    396395         ! 
     396         ! 1.2: get volume flux after coupling (>0 out) 
    397397         ! properly mask velocity  
    398398         ! (velocity are still mask with old mask at this stage) 
    399399         un(:,:,jk) = un(:,:,jk) * umask(:,:,jk) 
    400400         vn(:,:,jk) = vn(:,:,jk) * vmask(:,:,jk) 
    401          ! 
    402          ! get volume flux after coupling (>0 out) 
     401         ! compute volume flux divergence after coupling 
    403402         DO jj = 2, jpjm1 
    404403            DO ji = 2, jpim1 
     
    409408         ENDDO 
    410409         ! 
    411          ! get 3d volume flux difference (before - after cpl) (>0 out) 
    412          ! correction to add is _b - _n 
     410         ! 1.3: get 3d volume flux difference (before - after cpl) (>0 out) 
     411         !      correction to add is _b - _n 
    413412         risfcpl_vol(:,:,jk) = zqvolb(:,:,jk) - zqvoln(:,:,jk) 
    414413      END DO 
    415414      ! 
    416       ! include the contribution of the vertical velocity in the volume flux correction 
     415      ! 2.0: include the contribution of the vertical velocity in the volume flux correction 
     416      ! 
    417417      DO jj = 2, jpjm1 
    418418         DO ji = 2, jpim1 
     
    426426      ENDDO 
    427427      ! 
    428       ! mask volume flux divergence correction 
     428      CALL lbc_lnk( 'iscpl', risfcpl_vol, 'T', 1. ) 
     429      ! 
     430      ! 3.0: set total correction (div, tra, ssh) 
     431      ! 
     432      ! 3.1: mask volume flux divergence correction 
    429433      risfcpl_vol(:,:,:) = risfcpl_vol(:,:,:) * tmask(:,:,:) 
    430434      ! 
    431       CALL lbc_lnk( 'iscpl', risfcpl_vol, 'T', 1. ) 
    432       ! 
    433       ! get 3d tra increment to apply at the first time step 
     435      ! 3.2: get 3d tra increment to apply at the first time step 
    434436      ! temperature and salt content flux computed using local tsn  
    435437      ! (very simple advection scheme) 
     
    438440      risfcpl_tsc(:,:,:,jp_sal) = -risfcpl_vol(:,:,:) * tsn(:,:,:,jp_sal) 
    439441      ! 
    440       ! ssh correction (for dynspg_ts) 
     442      ! 3.3: ssh correction (for dynspg_ts) 
    441443      risfcpl_ssh(:,:) = 0.0 
    442444      DO jk = 1,jpk 
     
    459461      !!---------------------------------------------------------------------- 
    460462      ! 
    461       TYPE(isfcons), DIMENSION(:),ALLOCATABLE :: zisfpts  
    462       ! 
    463       INTEGER  ::   ji   , jj  , jk           ! loop index 
    464       INTEGER  ::   jip1 , jim1, jjp1, jjm1 
    465       INTEGER  ::   iig  , ijg, ik 
    466       INTEGER  ::   istart, iend, jisf 
    467       INTEGER  ::   nisfg , ingb, ifind 
    468       INTEGER, DIMENSION(jpnij) :: nisfl 
     463      TYPE(isfcons), DIMENSION(:),ALLOCATABLE :: zisfpts ! list of point receiving a correction 
     464      ! 
     465      INTEGER  ::   ji   , jj  , jk  , jproc          ! loop index 
     466      INTEGER  ::   jip1 , jim1, jjp1, jjm1           ! dummy indices 
     467      INTEGER  ::   iig  , ijg, ik                    ! dummy indices 
     468      INTEGER  ::   jisf                              ! start, end and current position in the increment array 
     469      INTEGER  ::   ingb, ifind                       ! 0/1 target found or need to be found  
     470      INTEGER  ::   ingb, ifind, nisfl_area           ! global number of cell concerned by the wet->dry case  
     471      INTEGER, DIMENSION(jpnij) :: nisfl              ! local  number of cell concerned by the wet->dry case 
    469472      ! 
    470473      REAL(wp) ::   z1_sum, z1_rdtiscpl 
    471       REAL(wp) ::   zdtem, zdsal, zdvol, zratio 
    472       REAL(wp) ::   zlon , zlat 
    473       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztmask_b    !! mask before 
    474       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t_b      !! scale factor before 
    475       !!---------------------------------------------------------------------- 
     474      REAL(wp) ::   zdtem, zdsal, zdvol, zratio       ! tem, sal, vol increment 
     475      REAL(wp) ::   zlon , zlat                       ! target location   
     476      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztmask_b    ! mask before 
     477      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t_b      ! scale factor before 
     478      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zt_b      ! scale factor before 
     479      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zs_b      ! scale factor before 
     480      !!---------------------------------------------------------------------- 
     481 
     482      !============================================================================== 
     483      ! 1.0: initialisation 
     484      !============================================================================== 
    476485 
    477486      ! get restart variable 
    478       CALL iom_get( numror, jpdom_autoglo, 'tmask'  , ztmask_b, ldxios = lrxios   ) ! need to extrapolate T/S 
    479       CALL iom_get( numror, jpdom_autoglo, 'e3t_n'  , ze3t_b(:,:,:), ldxios = lrxios ) 
     487      CALL iom_get( numror, jpdom_autoglo, 'tmask'  , ztmask_b(:,:,:), ldxios = lrxios   ) ! need to extrapolate T/S 
     488      CALL iom_get( numror, jpdom_autoglo, 'e3t_n'  , ze3t_b(:,:,:)  , ldxios = lrxios ) 
     489      CALL iom_get( numror, jpdom_autoglo, 'tn'     , zt_b(:,:,:)    , ldxios = lrxios ) 
     490      CALL iom_get( numror, jpdom_autoglo, 'sn'     , zs_b(:,:,:)    , ldxios = lrxios ) 
    480491 
    481492      ! compute run length 
     
    483494      rdt_iscpl   = nstp_iscpl * rn_rdt 
    484495      z1_rdtiscpl = 1._wp / rdt_iscpl  
     496 
    485497      IF (lwp) WRITE(numout,*) '            nb of stp for cons  = ', nstp_iscpl 
    486498      IF (lwp) WRITE(numout,*) '            coupling time step  = ', rdt_iscpl 
    487499 
    488       ! mask tsn and tsb  
    489       tsb(:,:,:,jp_tem) = tsb(:,:,:,jp_tem) * ztmask_b(:,:,:) 
    490       tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) *  tmask  (:,:,:) 
    491       tsb(:,:,:,jp_sal) = tsb(:,:,:,jp_sal) * ztmask_b(:,:,:) 
    492       tsn(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) *  tmask  (:,:,:) 
    493  
    494       !============================================================================== 
    495       ! diagnose the heat, salt and volume input and compute the correction variable 
    496       !============================================================================== 
    497  
     500      ! initialisation correction 
    498501      risfcpl_cons_vol = 0.0 
    499502      risfcpl_cons_ssh = 0.0 
    500503      risfcpl_cons_tsc = 0.0 
    501504 
     505      !============================================================================== 
     506      ! 2.0: diagnose the heat, salt and volume input and compute the correction variable 
     507      !      for case where we wet a cell or cell still wet (no change in cell status) 
     508      !============================================================================== 
     509 
    502510      DO jk = 1,jpk-1 
    503511         DO jj = nldj,nlej 
     
    508516 
    509517               ! heat diff 
    510                zdtem = tsn(ji,jj,jk,jp_tem) *  e3t_n(ji,jj,jk)   & 
    511                      - tsb(ji,jj,jk,jp_tem) * ze3t_b(ji,jj,jk) 
     518               zdtem = tsn (ji,jj,jk,jp_tem) *  e3t_n(ji,jj,jk) *  tmask  (ji,jj,jk)   & 
     519                     - zt_b(ji,jj,jk)        * ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
    512520 
    513521               ! salt diff 
    514                zdsal = tsn(ji,jj,jk,jp_sal) *  e3t_n(ji,jj,jk)   & 
    515                      - tsb(ji,jj,jk,jp_sal) * ze3t_b(ji,jj,jk) 
     522               zdsal = tsn(ji,jj,jk,jp_sal) *  e3t_n(ji,jj,jk) *  tmask  (ji,jj,jk)   & 
     523                     - zs_b(ji,jj,jk)       * ze3t_b(ji,jj,jk) * ztmask_b(ji,jj,jk) 
    516524             
    517525               ! volume, heat and salt differences in each cell (>0 means correction is an outward flux) 
    518                risfcpl_cons_vol(ji,jj,jk)        =   zdvol * e1e2t(ji,jj) * z1_rdtiscpl 
    519                risfcpl_cons_tsc(ji,jj,jk,jp_sal) = - zdsal * e1e2t(ji,jj) * z1_rdtiscpl 
    520                risfcpl_cons_tsc(ji,jj,jk,jp_tem) = - zdtem * e1e2t(ji,jj) * z1_rdtiscpl 
     526               ! in addition to the geometry change unconservation, need to add the divergence correction as it is flux across the boundary 
     527               risfcpl_cons_vol(ji,jj,jk)        = (   zdvol * e1e2t(ji,jj) + risfcpl_vol(ji,jj,jk)        ) * z1_rdtiscpl 
     528               risfcpl_cons_tsc(ji,jj,jk,jp_sal) = ( - zdsal * e1e2t(ji,jj) + risfcpl_tsc(ji,jj,jk,jp_sal) ) * z1_rdtiscpl 
     529               risfcpl_cons_tsc(ji,jj,jk,jp_tem) = ( - zdtem * e1e2t(ji,jj) + risfcpl_tsc(ji,jj,jk,jp_tem) ) * z1_rdtiscpl 
    521530 
    522531            END DO 
     
    524533      END DO 
    525534      ! 
    526       ! redistribute on valid point the vol/heat/salt removed during the coupling (ie when we dry a cell) 
    527       ! where we dry a cell get the number of point 
     535      !============================================================================== 
     536      ! 3.0: diagnose the heat, salt and volume input and compute the correction variable 
     537      !      for case where we close a cell 
     538      !============================================================================== 
     539      ! 
    528540      ! compute the total number of point receiving a correction increment for each processor 
    529541      ! local 
    530       nisfl=0 
     542      nisfl(:)=0 
    531543      DO jk = 1,jpk-1 
    532544         DO jj = nldj,nlej 
     
    539551      ! 
    540552      ! global  
    541       CALL mpp_sum('isfcpl',nisfl) 
    542       nisfg  = SUM(nisfl           ) 
    543       istart = SUM(nisfl(1:narea-1)) 
    544       iend   = SUM(nisfl(1:narea  )) 
     553      CALL mpp_sum('isfcpl',nisfl  ) 
    545554      ! 
    546555      ! allocate list of point receiving correction 
    547       ALLOCATE(zisfpts(nisfg)) 
     556      ALLOCATE(zisfpts(nisfl(narea))) 
     557      ! 
    548558      zisfpts(:) = isfcons(0,0,0,-HUGE(1.0), -HUGE(1.0), -HUGE(1.0), -HUGE(1.0), -HUGE(1.0), 0) 
    549559      ! 
    550560      ! start computing the correction and fill zisfpts 
    551561      ! local 
    552       jisf = istart 
     562      jisf = 0 
    553563      DO jk = 1,jpk-1 
    554564         DO jj = nldj,nlej 
     
    593603      ! 
    594604      ! share data among all processes because for some point we need to find the closest wet point (could be on other process) 
    595       DO jisf = 1,nisfg 
    596          ! 
    597          ! indices (conversion to global indices and sharing) 
    598          iig = zisfpts(jisf)%ii       ; ijg = zisfpts(jisf)%jj       ; ik = zisfpts(jisf)%kk 
    599          CALL mpp_max('isfcpl',iig)   ; CALL mpp_max('isfcpl',ijg)   ; CALL mpp_max('isfcpl',ik) 
    600          ! 
    601          ! data 
    602          zdvol = zisfpts(jisf)%dvol   ; zdsal = zisfpts(jisf)%dsal   ; zdtem = zisfpts(jisf)%dtem 
    603          CALL mpp_max('isfcpl',zdvol) ; CALL mpp_max('isfcpl',zdsal) ; CALL mpp_max('isfcpl',zdtem) 
    604          ! 
    605          ! location 
    606          zlat = zisfpts(jisf)%lat     ; zlon = zisfpts(jisf)%lon 
    607          CALL mpp_max('isfcpl',zlat)  ; CALL mpp_max('isfcpl',zlon) 
    608          ! 
    609          ! find flag 
    610          ingb = zisfpts(jisf)%ngb 
    611          CALL mpp_max('isfcpl',ingb) 
    612          ! 
    613          ! fill the 3d correction array 
    614          CALL get_correction(iig, ijg, ik, zlon, zlat, zdvol, zdsal, zdtem, ingb) 
    615          ! 
    616       END DO 
     605      DO jproc=1,jpnij 
     606         !  
     607         ! share total number of isf point treated for proc jproc 
     608         IF (jproc==narea) THEN 
     609            nisfl_area=nisfl(jproc) 
     610         ELSE 
     611            nisfl_area=0 
     612         END IF 
     613         CALL mpp_max('isfcpl',nisfl_area) 
     614         ! 
     615         DO jisf = 1,nisfl_area 
     616            ! 
     617            IF (jproc==narea) THEN 
     618               ! indices (conversion to global indices and sharing) 
     619               iig = zisfpts(jisf)%ii       ; ijg = zisfpts(jisf)%jj       ; ik = zisfpts(jisf)%kk 
     620               ! 
     621               ! data 
     622               zdvol = zisfpts(jisf)%dvol   ; zdsal = zisfpts(jisf)%dsal   ; zdtem = zisfpts(jisf)%dtem 
     623               ! 
     624               ! location 
     625               zlat = zisfpts(jisf)%lat     ; zlon = zisfpts(jisf)%lon 
     626               ! 
     627               ! find flag 
     628               ingb = zisfpts(jisf)%ngb 
     629            ELSE 
     630               iig  =0   ; ijg  =0   ; ik   =0   
     631               zdvol=-HUGE(1.0) ; zdsal=-HUGE(1.0) ; zdtem=-HUGE(1.0) 
     632               zlat =-HUGE(1.0) ; zlon =-HUGE(1.0)    
     633               ingb = 0 
     634            END IF 
     635            ! 
     636            ! share data (need synchronisation of data as get_correction call a global com) 
     637            CALL mpp_max('isfcpl',iig)   ; CALL mpp_max('isfcpl',ijg)   ; CALL mpp_max('isfcpl',ik) 
     638            CALL mpp_max('isfcpl',zdvol) ; CALL mpp_max('isfcpl',zdsal) ; CALL mpp_max('isfcpl',zdtem) 
     639            CALL mpp_max('isfcpl',zlat)  ; CALL mpp_max('isfcpl',zlon) 
     640            CALL mpp_max('isfcpl',ingb) 
     641            ! 
     642            ! fill the 3d correction array 
     643            CALL get_correction(iig, ijg, ik, zlon, zlat, zdvol, zdsal, zdtem, ingb) 
     644         END DO 
     645      END DO 
     646      ! 
     647      !============================================================================== 
     648      ! 4.0: finalisation and compute ssh equivalent of the volume correction 
     649      !============================================================================== 
    617650      ! 
    618651      ! mask (>0 out) 
     
    622655      ! 
    623656      ! add lbclnk 
    624       CALL lbc_lnk_multi( 'iscplrst', risfcpl_cons_tsc(:,:,:,jp_tem), 'T', 1., risfcpl_cons_tsc(:,:,:,jp_sal), 'T', 1., risfcpl_cons_vol, 'T', 1.) 
     657      CALL lbc_lnk_multi( 'iscplrst', risfcpl_cons_tsc(:,:,:,jp_tem), 'T', 1., risfcpl_cons_tsc(:,:,:,jp_sal), 'T', 1., & 
     658         &                            risfcpl_cons_vol(:,:,:)       , 'T', 1.) 
    625659      ! 
    626660      ! ssh correction (for dynspg_ts) 
     
    643677      !!---------------------------------------------------------------------- 
    644678      TYPE(isfcons), DIMENSION(:), INTENT(inout) :: sisfpts 
    645       INTEGER,      INTENT(inout) :: kpts 
    646       !!---------------------------------------------------------------------- 
    647       INTEGER,      INTENT(in   ) :: ki, kj, kk 
    648       INTEGER,      INTENT(in   ), OPTIONAL :: kfind 
    649       REAL(wp),     INTENT(in   ) :: pdvol, pdsal, pdtem, pratio 
     679      INTEGER,                     INTENT(inout) :: kpts 
     680      !!---------------------------------------------------------------------- 
     681      INTEGER,      INTENT(in   )           :: ki, kj, kk                  !    target location (kfind=0)  
     682      !                                                                    ! or source location (kfind=1) 
     683      INTEGER,      INTENT(in   ), OPTIONAL :: kfind                       ! 0  target cell already found 
     684      !                                                                    ! 1  target to be determined 
     685      REAL(wp),     INTENT(in   )           :: pdvol, pdsal, pdtem, pratio ! vol/sal/tem increment  
     686      !                                                                    ! and ratio in case increment span over multiple cells. 
    650687      !!---------------------------------------------------------------------- 
    651688      INTEGER :: ifind 
     
    667704   END SUBROUTINE update_isfpts 
    668705   ! 
    669    SUBROUTINE get_correction( ki, kj, kk, zlon, zlat, pvolinc, psalinc, pteminc, kfind) 
     706   SUBROUTINE get_correction( ki, kj, kk, plon, plat, pvolinc, psalinc, pteminc, kfind) 
    670707      !!--------------------------------------------------------------------- 
    671708      !!                  ***  ROUTINE get_correction  *** 
     
    676713      !! 
    677714      !!---------------------------------------------------------------------- 
    678       INTEGER , INTENT(in)           :: ki, kj, kk, kfind        ! target point 
    679       REAL(wp), INTENT(in)           :: zlon, zlat 
    680       REAL(wp), INTENT(in)           :: pvolinc, pteminc,psalinc ! correction increment for vol/temp/salt 
     715      INTEGER , INTENT(in) :: ki, kj, kk, kfind        ! target point indices 
     716      REAL(wp), INTENT(in) :: plon, plat               ! target point lon/lat 
     717      REAL(wp), INTENT(in) :: pvolinc, pteminc,psalinc ! correction increment for vol/temp/salt 
    681718      !!---------------------------------------------------------------------- 
    682719      INTEGER :: jj, ji, iig, ijg 
     
    685722      ! define global indice of correction location 
    686723      iig = ki ; ijg = kj 
    687       IF ( kfind == 1 ) CALL dom_ngb( zlon, zlat, iig, ijg,'T', kk) 
     724      IF ( kfind == 1 ) CALL dom_ngb( plon, plat, iig, ijg,'T', kk) 
    688725      ! 
    689726      ! fill the correction array 
    690727      DO jj = mj0(ijg),mj1(ijg) 
    691728         DO ji = mi0(iig),mi1(iig) 
    692             ! correct the vol_flx in the closest cell 
     729            ! correct the vol_flx and corresponding heat/salt flx in the closest cell 
    693730            risfcpl_cons_vol(ji,jj,kk)        =  risfcpl_cons_vol(ji,jj,kk       ) + pvolinc 
    694731            risfcpl_cons_tsc(ji,jj,kk,jp_sal) =  risfcpl_cons_tsc(ji,jj,kk,jp_sal) + psalinc 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfdiags.F90

    r11521 r11931  
    4646      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: phtbl, pfrac              ! thickness of the tbl and fraction of last cell affected by the tbl 
    4747      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqfwf, pqoce, pqlat, pqhc ! 2d var to map in 3d 
    48       CHARACTER(LEN=256), INTENT(in) :: cdisf                               ! parametrisation or interactive melt 
     48      CHARACTER(LEN=3), INTENT(in) :: cdisf                                 ! parametrisation or interactive melt 
    4949      !!--------------------------------------------------------------------- 
    5050      CHARACTER(LEN=256) :: cvarqfwf  , cvarqoce  , cvarqlat  , cvarqhc 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfdynnxt.F90

    r11852 r11931  
    1 MODULE isfnxt 
     1MODULE isfdynnxt 
    22   !!========================================================================= 
    33   !!                       ***  MODULE  isfnxt  *** 
     
    3232      !! ** Purpose : compute the ice shelf volume filter correction for cavity, param, ice sheet coupling case 
    3333      !! 
    34       !!-------------------------------------------------------------------- 
    3534      !!-------------------------- OUT ------------------------------------- 
    3635      INTEGER ,                     INTENT(in   ) :: kt 
    3736      ! 
    3837      REAL(wp),                     INTENT(in   ) :: pcoef           ! atfp * rdt * r1_rau0 
    39       !!-------------------------- IN  ------------------------------------- 
    4038      !!-------------------------------------------------------------------- 
    4139      INTEGER :: jk  ! loop index 
     
    6260      !! ** Purpose : compute the ice shelf volume filter correction for cavity or param 
    6361      !! 
    64       !!-------------------------------------------------------------------- 
    65       !!-------------------------- OUT ------------------------------------- 
    6662      !!-------------------------- IN  ------------------------------------- 
    6763      INTEGER , DIMENSION(jpi,jpj), INTENT(in   ) :: ktop , kbot     ! top and bottom level of tbl 
     
    9288   END SUBROUTINE isf_dynnxt_mlt 
    9389 
    94 END MODULE isfnxt 
     90END MODULE isfdynnxt 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfhdiv.F90

    r11902 r11931  
    4646            ! conservation option 
    4747            IF ( ln_isfcpl_cons ) CALL isf_hdiv_cpl(risfcpl_cons_vol, phdiv) 
    48             ! 
    49             IF ( ln_isfdebug ) THEN 
    50                CALL debug('isfdiv: phdiv'      ,phdiv(:,:,:)) 
    51                CALL debug('isfdiv: risfcpl_vol',risfcpl_vol(:,:,:)) 
    52                CALL debug('isfdiv: fwfisf     ',fwfisf_cav(:,:)+fwfisf_cav_b(:,:)) 
    53             END IF 
    5448            ! 
    5549         END IF 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfpar.F90

    r11876 r11931  
    1414   !!   isfpar       : compute ice shelf melt using a prametrisation of ice shelf cavities 
    1515   !!---------------------------------------------------------------------- 
    16    !USE oce            ! ocean dynamics and tracers 
    1716   USE isf            ! ice shelf 
    1817   ! 
     
    2322   USE isfdiags , ONLY: isf_diags_flx  ! ice shelf diags subroutine 
    2423   ! 
    25    USE dom_oce  , ONLY: jpi,jpj, bathy ! ocean space and time domain 
     24   USE dom_oce  , ONLY: bathy          ! ocean space and time domain 
     25   USE par_oce  , ONLY: jpi,jpj        ! ocean space and time domain 
    2626   USE phycst   , ONLY: r1_rau0_rcp    ! physical constants 
    2727   ! 
     
    6262      CALL isfpar_mlt( kt, zqhc, zqoce, pqfwf  ) 
    6363      ! 
    64       ! compute heat and water flux  (change signe directly in the melt subroutine) 
     64      ! compute heat and water flux ( > 0 out ) 
    6565      pqfwf(:,:) = pqfwf(:,:) * mskisf_par(:,:) 
    6666      zqoce(:,:) = zqoce(:,:) * mskisf_par(:,:) 
    6767      zqhc (:,:) = zqhc(:,:)  * mskisf_par(:,:) 
    6868      ! 
    69       ! compute heat content flux 
     69      ! compute heat content flux ( > 0 out ) 
    7070      zqlat(:,:) = pqfwf(:,:) * rLfusisf    ! 2d latent heat flux (W/m2) 
    7171      ! 
    72       ! total heat flux 
     72      ! total heat flux ( > 0 out ) 
    7373      zqh(:,:) = ( zqhc (:,:) + zqoce(:,:) ) 
    7474      ! 
     
    107107      CALL isf_alloc_par() 
    108108      ! 
    109       ! par 
     109      ! initialisation 
    110110      misfkt_par(:,:)     = 1         ; misfkb_par(:,:)       = 1          
    111111      rhisf_tbl_par(:,:)  = 1e-20     ; rfrac_tbl_par(:,:)    = 0.0_wp 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfparmlt.F90

    r11876 r11931  
    1818   USE in_out_manager              ! I/O manager 
    1919   USE iom        , ONLY: iom_put  ! I/O library 
    20    USE fldread    , ONLY: fld_read ! 
     20   USE fldread    , ONLY: fld_read, FLD, FLD_N ! 
    2121   USE lib_fortran, ONLY: glob_sum ! 
    2222   USE lib_mpp    , ONLY: ctl_stop ! 
     
    4949      !!                        1 : Specified melt flux 
    5050      !!                        2 : Beckmann & Goose parameterization 
    51       !!---------------------------------------------------------------------- 
     51      !! 
    5252      !!-------------------------- OUT ------------------------------------- 
    5353      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pqfwf, pqoce, pqhc  ! fresh water, ice-ocean heat and heat content fluxes 
     
    8181      !!              data read into a forcing files. 
    8282      !! 
    83       !!---------------------------------------------------------------------- 
    8483      !!-------------------------- OUT ------------------------------------- 
    8584      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pqhc, pqfwf, pqoce  ! fresh water and ice-ocean heat fluxes 
     
    105104      pqoce(:,:) =   pqfwf(:,:) * rLfusisf             ! ocean/ice shelf flux assume to be equal to latent heat flux 
    106105      pqhc (:,:) =   pqfwf(:,:) * ztfrz(:,:) * rcp     ! heat content flux  
     106      ! 
    107107      CALL iom_put('isftfrz_par', ztfrz ) 
    108108      ! 
     
    155155      ! output thermal driving 
    156156      CALL iom_put('isfthermald_par',( ztfrz(:,:) - ztavg(:,:) ) * mskisf_par(:,:)) 
     157      ! 
     158      ! output freezing point used to define the thermal driving and heat content fluxes 
    157159      CALL iom_put('isftfrz_par', ztfrz ) 
    158160      ! 
     
    206208      zfwf(:,:) = zfwf(:,:) * zfwf_oasis / zfwf_fld 
    207209      !  
    208       ! i3. -----------Define fwf and qoce 
     210      ! 3. -----------Define fwf and qoce 
    209211      ! ocean heat flux is assume to be equal to the latent heat 
    210212      pqfwf(:,:) =   zfwf(:,:)                         ! fwf                ( >0 out ) 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfrst.F90

    r11852 r11931  
    1111   !!---------------------------------------------------------------------- 
    1212   ! 
    13    USE dom_oce, ONLY: jpi,jpj,jpk,jpts ! time and space domain 
     13   USE par_oce, ONLY: jpi,jpj,jpk,jpts ! time and space domain 
    1414   ! 
    1515   USE in_out_manager ! I/O manager 
     
    3131   SUBROUTINE isfrst_read(cdisf, ptsc, pfwf, ptsc_b, pfwf_b ) 
    3232      !!--------------------------------------------------------------------- 
     33      !! 
    3334      !!   isfrst_read : read iceshelf variables from restart 
    34       !!--------------------------------------------------------------------- 
     35      !! 
    3536      !!-------------------------- OUT -------------------------------------- 
    3637      REAL(wp), DIMENSION(jpi,jpj)     , INTENT(  out) :: pfwf_b 
    3738      REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(  out) :: ptsc_b 
    3839      !!-------------------------- IN  -------------------------------------- 
    39       CHARACTER(LEN=256)               , INTENT(in   ) :: cdisf 
     40      CHARACTER(LEN=3)                 , INTENT(in   ) :: cdisf 
    4041      REAL(wp), DIMENSION(jpi,jpj)     , INTENT(in   ) :: pfwf 
    4142      REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(in   ) :: ptsc 
     
    6667      ENDIF 
    6768 
    68       CALL FLUSH(numout) 
    69  
    7069   END SUBROUTINE isfrst_read 
    7170   !  
    7271   SUBROUTINE isfrst_write(kt, cdisf, ptsc, pfwf ) 
    7372      !!--------------------------------------------------------------------- 
     73      !! 
    7474      !!   isfrst_write : write iceshelf variables in restart 
    75       !!--------------------------------------------------------------------- 
    76       !!-------------------------- OUT -------------------------------------- 
     75      !! 
    7776      !!-------------------------- IN  -------------------------------------- 
    7877      INTEGER                          , INTENT(in   ) :: kt 
    79       CHARACTER(LEN=256)               , INTENT(in   ) :: cdisf 
     78      CHARACTER(LEN=3)                 , INTENT(in   ) :: cdisf 
    8079      REAL(wp), DIMENSION(jpi,jpj)     , INTENT(in   ) :: pfwf 
    8180      REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(in   ) :: ptsc 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90

    r11876 r11931  
    2626   ! 
    2727   USE lib_mpp, ONLY: ctl_stop, ctl_nam 
     28   USE fldread, ONLY: FLD, FLD_N 
    2829   USE in_out_manager ! I/O manager 
    2930   USE timing 
     
    6263      IF( ln_timing )   CALL timing_start('isf') 
    6364      ! 
     65      !======================================================================= 
     66      ! 1.: compute melt and associated heat fluxes in the ice shelf cavities 
     67      !======================================================================= 
     68      ! 
    6469      IF ( ln_isfcav_mlt ) THEN 
    6570         ! 
    66          ! before time step  
     71         ! 1.1: before time step  
    6772         IF ( kt /= nit000 ) THEN  
    6873            risf_cav_tsc_b (:,:,:) = risf_cav_tsc (:,:,:) 
     
    7075         END IF 
    7176         ! 
    72          ! compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 
     77         ! 1.2: compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 
    7378         rhisf_tbl_cav(:,:) = rn_htbl * mskisf_cav(:,:) 
    7479         CALL isf_tbl_lvl(ht_n, e3t_n, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav) 
    7580         ! 
    76          ! compute ice shelf melt 
     81         ! 1.3: compute ice shelf melt 
    7782         CALL isf_cav( kt, risf_cav_tsc, fwfisf_cav) 
    7883         ! 
    7984      END IF 
    8085      !  
     86      !================================================================================= 
     87      ! 2.: compute melt and associated heat fluxes for not resolved ice shelf cavities 
     88      !================================================================================= 
     89      ! 
    8190      IF ( ln_isfpar_mlt ) THEN 
    8291         ! 
    83          ! before time step  
     92         ! 2.1: before time step  
    8493         IF ( kt /= nit000 ) THEN  
    8594            risf_par_tsc_b(:,:,:) = risf_par_tsc(:,:,:) 
     
    8796         END IF 
    8897         ! 
    89          ! compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 
    90          ! by simplicity, we assume the top level where param applied do not change with time 
     98         ! 2.2: compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 
     99         ! by simplicity, we assume the top level where param applied do not change with time (done in init part) 
    91100         rhisf_tbl_par(:,:) = rhisf0_tbl_par(:,:) 
    92101         CALL isf_tbl_lvl(ht_n, e3t_n, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par) 
    93102         ! 
    94          ! compute ice shelf melt 
     103         ! 2.3: compute ice shelf melt 
    95104         CALL isf_par( kt, risf_par_tsc, fwfisf_par) 
    96105         ! 
    97106      END IF 
     107      ! 
     108      !================================================================================== 
     109      ! 3.: output specific restart variable in case of coupling with an ice sheet model 
     110      !================================================================================== 
    98111      ! 
    99112      IF ( ln_isfcpl .AND. lrst_oce ) CALL isfcpl_rst_write(kt) 
     
    230243      IF ( l_isfoasis .AND. ln_isf ) THEN 
    231244         ! 
    232          CALL ctl_stop( ' ln_ctl and ice shelf not tested' ) 
     245         CALL ctl_stop( ' OASIS and ice shelf not tested' ) 
    233246         ! 
    234247         ! NEMO coupled to ATMO model with isf cavity need oasis method for melt computation  
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isftbl.F90

    r11876 r11931  
    3737      !!                https://doi.org/10.1029/2007JC004368 , 2008 
    3838      !! 
    39       !!-------------------------------------------------------------------- 
    4039      !!-------------------------- OUT ------------------------------------- 
    4140      REAL(wp), DIMENSION(jpi,jpj)          , INTENT(  out) :: pvarout ! 2d average of pvarin 
     
    115114      !!              over a thickness phtbl. The bottom level is partially counted (pfrac). 
    116115      !! 
    117       !!-------------------------------------------------------------------- 
    118116      !!-------------------------- OUT ------------------------------------- 
    119117      REAL(wp), DIMENSION(jpi,jpj)    , INTENT(  out) :: pvarout      ! tbl property averaged over phtbl between level ktop and kbot 
     
    154152      !!              - fraction of the bottom level affected by the tbl 
    155153      !! 
    156       !!--------------------------------------------------------------------- 
    157154      !!-------------------------- OUT -------------------------------------- 
    158155      INTEGER,  DIMENSION(jpi,jpj)    , INTENT(  out) :: kbot   ! bottom level of the top boundary layer 
     
    209206      !! ** Purpose : compute bottom level of the isf top boundary layer 
    210207      !! 
    211       !!-------------------------------------------------------------------- 
    212208      !!-------------------------- OUT ------------------------------------- 
    213209      INTEGER,  DIMENSION(jpi,jpj)    , INTENT(  out) :: kbot   ! bottom level of the top boundary layer 
     
    246242      !! ** Purpose : compute top level of the isf top boundary layer in case of an ice shelf parametrisation 
    247243      !! 
    248       !!-------------------------------------------------------------------- 
    249244      !!-------------------------- OUT ------------------------------------- 
    250245      INTEGER,  DIMENSION(jpi,jpj), INTENT(  out) :: ktop        ! top level affected by the ice shelf parametrisation 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfutils.F90

    r11852 r11931  
    1414   USE iom           , ONLY: iom_open, iom_get, iom_close, jpdom_data ! read input file 
    1515   USE lib_fortran   , ONLY: glob_sum, glob_min, glob_max             ! compute global value 
    16    USE dom_oce       , ONLY: jpi,jpj,jpk                              ! domain size 
     16   USE par_oce       , ONLY: jpi,jpj,jpk                              ! domain size 
    1717   USE in_out_manager, ONLY: wp, lwp, numout                          ! miscelenious 
    1818 
     
    3535      !! ** Purpose : read input file 
    3636      !! 
    37       !!-------------------------------------------------------------------- 
    3837      !!-------------------------- OUT ------------------------------------- 
    3938      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pvar          ! output variable 
     
    5756      !! ** Purpose : add debug print 
    5857      !! 
    59       !!-------------------------------------------------------------------- 
    60       !!-------------------------- OUT ------------------------------------- 
    6158      !!-------------------------- IN  ------------------------------------- 
    6259      CHARACTER(LEN=256)          , INTENT(in   ) :: cdtxt 
     
    8279      !! ** Purpose : add debug print 
    8380      !! 
    84       !!-------------------------------------------------------------------- 
    85       !!-------------------------- OUT ------------------------------------- 
    8681      !!-------------------------- IN  ------------------------------------- 
    8782      CHARACTER(LEN=256)              , INTENT(in   ) :: cdtxt 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/TRA/traisf.F90

    r11908 r11931  
    6868         ! ensure 0 trend due to unconservation of the ice shelf coupling 
    6969         IF ( ln_isfcpl_cons ) CALL tra_isf_cpl(risfcpl_cons_tsc, tsa) 
    70  
    71          IF ( ln_isfdebug ) THEN 
    72             CALL debug('tra_isf: risfcpl_tsc T',risfcpl_tsc(:,:,1)) 
    73             CALL debug('tra_isf: risfcpl_tsc S',risfcpl_tsc(:,:,2)) 
    74          END IF 
    7570         ! 
    7671      END IF 
    7772      ! 
    7873      IF ( ln_isfdebug ) THEN 
    79          CALL debug('tra_isf: tsa T'        ,tsa(:,:,:,1)) 
    80          CALL debug('tra_isf: tsa S'        ,tsa(:,:,:,2)) 
     74         CALL debug('tra_isf: tsa T', tsa(:,:,:,1)) 
     75         CALL debug('tra_isf: tsa S', tsa(:,:,:,2)) 
    8176      END IF 
    8277      ! 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/TRA/tranxt.F90

    r11541 r11931  
    341341                  IF( ll_isf ) THEN 
    342342                     ! 
     343                     ! melt in the cavity 
    343344                     IF ( ln_isfcav_mlt ) THEN 
    344345                        ! level fully include in the Losch_2008 ice shelf boundary layer 
     
    357358                        END IF 
    358359                     END IF 
     360                     ! 
     361                     ! parametrised melt (cavity closed) 
    359362                     IF ( ln_isfpar_mlt ) THEN 
    360363                        ! level fully include in the Losch_2008 ice shelf boundary layer 
     
    374377                     END IF 
    375378                     ! 
    376                      IF (ln_isfcpl .AND. ln_rstart .AND. kt == nit000+1 ) THEN    ! risfcpl_vol_n = 0 and risfcpl_vol_b = risfcpl_vol 
    377                         ztc_f  = ztc_f  + zfact1 * risfcpl_tsc(ji,jj,jk,jn) * r1_e1e2t(ji,jj) 
    378                         ze3t_f = ze3t_f - zfact1 * risfcpl_vol(ji,jj,jk   ) * r1_e1e2t(ji,jj) 
     379                     ! ice sheet coupling correction 
     380                     IF ( ln_isfcpl ) THEN 
     381                        ! 
     382                        ! at kt = nit000,  risfcpl_vol_n = 0 and risfcpl_vol_b = risfcpl_vol so contribution nul 
     383                        IF ( ln_rstart .AND. kt == nit000+1 ) THEN 
     384                           ztc_f  = ztc_f  + zfact1 * risfcpl_tsc(ji,jj,jk,jn) * r1_e1e2t(ji,jj) 
     385                           ze3t_f = ze3t_f - zfact1 * risfcpl_vol(ji,jj,jk   ) * r1_e1e2t(ji,jj) 
     386                        END IF 
     387                        ! 
    379388                     END IF 
    380389                     ! 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/step.F90

    r11895 r11931  
    241241                         CALL tra_sbc       ( kstp )  ! surface boundary condition 
    242242      IF( ln_traqsr  )   CALL tra_qsr       ( kstp )  ! penetrative solar radiation qsr 
    243       IF( ln_isf     ) CALL tra_isf ( kstp )  ! ice shelf heat flux 
     243      IF( ln_isf     )   CALL tra_isf      ( kstp )  ! ice shelf heat flux 
    244244      IF( ln_trabbc  )   CALL tra_bbc       ( kstp )  ! bottom heat flux 
    245245      IF( ln_trabbl  )   CALL tra_bbl       ( kstp )  ! advective (and/or diffusive) bottom boundary layer scheme 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/tests/ISOMIP+/MY_SRC/isf.F90

    r11889 r11931  
    2222   PRIVATE 
    2323 
    24    PUBLIC   isf_alloc, isf_alloc_par, isf_alloc_cav, isf_alloc_cpl 
     24   PUBLIC   isf_alloc, isf_alloc_par, isf_alloc_cav, isf_alloc_cpl, isf_dealloc_cpl 
    2525   ! 
    2626   !------------------------------------------------------- 
     
    175175      !!                  ***  ROUTINE isf_alloc_cpl  *** 
    176176      !! 
    177       !! ** Purpose :  
    178       !! 
    179       !! ** Method  :  
     177      !! ** Purpose : allocate array use for the ice sheet coupling 
    180178      !! 
    181179      !!---------------------------------------------------------------------- 
     
    202200   END SUBROUTINE isf_alloc_cpl 
    203201 
     202   SUBROUTINE isf_dealloc_cpl() 
     203      !!--------------------------------------------------------------------- 
     204      !!                  ***  ROUTINE isf_dealloc_cpl  *** 
     205      !! 
     206      !! ** Purpose : de-allocate useless public 3d array used for ice sheet coupling 
     207      !! 
     208      !!---------------------------------------------------------------------- 
     209      INTEGER :: ierr, ialloc 
     210      !!---------------------------------------------------------------------- 
     211      ierr = 0 
     212      ! 
     213      DEALLOCATE( risfcpl_ssh, risfcpl_tsc, risfcpl_vol, STAT=ialloc ) 
     214      ierr = ierr + ialloc 
     215      ! 
     216      CALL mpp_sum ( 'isf', ierr ) 
     217      IF( ierr /= 0 )   CALL ctl_stop('STOP','isfcpl: failed to deallocate arrays.') 
     218      ! 
     219   END SUBROUTINE isf_dealloc_cpl 
     220 
    204221   SUBROUTINE isf_alloc() 
    205222      !!--------------------------------------------------------------------- 
    206223      !!                  ***  ROUTINE isf_alloc  *** 
    207224      !! 
    208       !! ** Purpose :  
    209       !! 
    210       !! ** Method  :  
     225      !! ** Purpose : allocate array used for the ice shelf cavity (cav and par) 
    211226      !! 
    212227      !!---------------------------------------------------------------------- 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/tests/ISOMIP+/MY_SRC/isfstp.F90

    r11908 r11931  
    9797      END IF 
    9898      ! 
    99       IF ( ln_isfcpl .AND. lrst_oce ) CALL isfcpl_rst_write(kt) 
     99      IF ( ln_isfcpl ) THEN 
     100         ! after step nit000 + 2 we do not need anymore the risfcpl_ arrays 
     101         IF ( kt == nit000 + 2 ) CALL isf_dealloc_cpl() 
     102 
     103         IF ( lrst_oce ) CALL isfcpl_rst_write(kt) 
     104      END IF 
    100105      ! 
    101106      IF( ln_timing )   CALL timing_stop('isf') 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/tests/ISOMIP+/MY_SRC/sbcfwb.F90

    r11908 r11931  
    118118         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN 
    119119            z_fwf = glob_sum( 'sbcfwb',  e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf_cav(:,:) + fwfisf_par(:,:) - snwice_fmass(:,:) ) ) 
     120            ! 
     121            ! correction for ice sheet coupling testing (ie remove the excess through the surface) 
     122            ! test impact on the melt as conservation correction made in depth 
     123            ! test conservation level as sbcfwb is conserving 
     124            ! avoid the model to blow up for large ssh drop (isomip OCEAN3 with melt switch off and uniform T/S) 
     125            IF (ln_isfcpl .AND. ln_isfcpl_cons) THEN 
     126               z_fwf = z_fwf + glob_sum( 'sbcfwb',  e1e2t(:,:) * risfcpl_cons_ssh(:,:) * rau0 ) 
     127            END IF 
     128            ! 
    120129            z_fwf = z_fwf / area 
    121130            zcoef = z_fwf * rcp 
    122             emp(:,:) = emp(:,:) - z_fwf            * tmask(:,:,1) ! (Eq. 34 AD2015) 
    123             qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) * tmask(:,:,1) ! (Eq. 35 AD2015) ! could be sst_m if we don't want any bouyancy fluxes 
    124             sfx(:,:) = sfx(:,:) + z_fwf * sss_m(:,:) * tmask(:,:,1) ! (Eq. 36 AD2015) ! could be sss_m if we don't want any bouyancy fluxes 
     131            emp(:,:) = emp(:,:) - z_fwf              * tmask(:,:,1) ! (Eq. 34 AD2015) 
     132            qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) * tmask(:,:,1) ! (Eq. 35 AD2015) ! use sst_m to avoid generation of any bouyancy fluxes 
     133            sfx(:,:) = sfx(:,:) + z_fwf * sss_m(:,:) * tmask(:,:,1) ! (Eq. 36 AD2015) ! use sss_m to avoid generation of any bouyancy fluxes 
    125134            !qns(:,:) = qns(:,:) + zcoef * ( -1.9 ) * tmask(:,:,1) ! (Eq. 35 AD2015) ! could be sst_m if we don't want any bouyancy fluxes 
    126135            !sfx(:,:) = sfx(:,:) + z_fwf * ( 33.8 ) * tmask(:,:,1) ! (Eq. 36 AD2015) ! could be sss_m if we don't want any bouyancy fluxes 
     136            !qns(:,:) = qns(:,:) + zcoef * ( -1.0 ) * tmask(:,:,1) ! use for ISOMIP+ coupling sanity check (keep ssh cst while playing with cpl conservation option) 
     137            !sfx(:,:) = sfx(:,:) + z_fwf * ( 34.2 ) * tmask(:,:,1) ! use for ISOMIP+ coupling sanity check (keep ssh cst while playing with cpl conservation option) 
    127138         ENDIF 
    128139         ! 
Note: See TracChangeset for help on using the changeset viewer.