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 11852 for NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF – NEMO

Ignore:
Timestamp:
2019-11-04T19:00:27+01:00 (4 years ago)
Author:
mathiot
Message:

ENHANCE-02_ISF_nemo: fix WED025 restartability, finish removing useless USE, remove useless lbc_lnk

Location:
NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF
Files:
14 edited

Legend:

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

    r11844 r11852  
    1414   !!---------------------------------------------------------------------- 
    1515 
    16    USE in_out_manager ! I/O manager 
    17    USE lib_mpp        ! MPP library 
     16   USE in_out_manager, ONLY: wp, jpi,jpj, jpk, jpts ! I/O manager 
     17   USE lib_mpp       , ONLY: ctl_stop, mpp_sum      ! MPP library 
    1818   USE fldread        ! read input fields 
    1919 
     
    197197      ! 
    198198      CALL mpp_sum ( 'isf', ierr ) 
    199       IF( ierr /= 0 )   CALL ctl_warn('STOP','isfcpl: failed to allocate arrays.') 
     199      IF( ierr /= 0 )   CALL ctl_stop('STOP','isfcpl: failed to allocate arrays.') 
    200200      ! 
    201201   END SUBROUTINE isf_alloc_cpl 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcav.F90

    r11844 r11852  
    1313   !!   isf_cav       : update ice shelf melting under ice shelf 
    1414   !!---------------------------------------------------------------------- 
    15    USE oce            ! ocean dynamics and tracers 
    1615   USE isf            ! ice shelf public variables 
    17    USE isfutils 
    18    USE isftbl         ! ice shelf top boundary layer properties 
    19    USE isfcavmlt      ! ice shelf melt formulation 
    20    USE isfcavgam      ! ice shelf melt exchange coeficient 
    21    USE isfdiags       ! ice shelf diags 
    22    USE dom_oce        ! ocean space and time domain 
    23    USE phycst         ! physical constants 
    24    USE eosbn2         ! l_useCT 
     16   ! 
     17   USE isfrst   , ONLY: isfrst_write, isfrst_read ! ice shelf restart read/write subroutine 
     18   USE isfutils , ONLY: debug          ! ice shelf debug subroutine 
     19   USE isftbl   , ONLY: isf_tbl        ! ice shelf top boundary layer properties subroutine 
     20   USE isfcavmlt, ONLY: isfcav_mlt     ! ice shelf melt formulation subroutine 
     21   USE isfcavgam, ONLY: isfcav_gammats ! ice shelf melt exchange coeficient subroutine 
     22   USE isfdiags , ONLY: isf_diags_flx  ! ice shelf diags subroutine 
     23   ! 
     24   USE oce      , ONLY: tsn                    ! ocean tracers 
     25   USE dom_oce  , ONLY: jpi,jpj                ! ocean space and time domain 
     26   USE phycst   , ONLY: grav,rau0,r1_rau0_rcp  ! physical constants 
     27   USE eosbn2   , ONLY: l_useCT                ! l_useCT 
    2528   ! 
    2629   USE in_out_manager ! I/O manager 
     
    5861      !!--------------------------------------------------------------------- 
    5962      !!-------------------------- OUT -------------------------------------- 
    60       REAL(wp), DIMENSION(jpi,jpj)     , INTENT(inout) :: pqfwf 
    61       REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(inout) :: ptsc 
     63      REAL(wp), DIMENSION(jpi,jpj)     , INTENT(inout) :: pqfwf  ! ice shelf melt (>0 out) 
     64      REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(inout) :: ptsc   ! T & S ice shelf cavity contents 
    6265      !!-------------------------- IN  -------------------------------------- 
    6366      INTEGER, INTENT(in) ::   kt   ! ocean time step 
     
    6669      INTEGER :: nit 
    6770      REAL(wp) :: zerr 
    68       REAL(wp), DIMENSION(jpi,jpj) :: zqlat, zqoce, zqhc, zqh 
    69       REAL(wp), DIMENSION(jpi,jpj) :: zqoce_b 
    70       REAL(wp), DIMENSION(jpi,jpj) :: zgammat, zgammas 
    71       REAL(wp), DIMENSION(jpi,jpj) :: zttbl, zstbl 
     71      REAL(wp), DIMENSION(jpi,jpj) :: zqlat, zqoce, zqhc, zqh  ! heat fluxes 
     72      REAL(wp), DIMENSION(jpi,jpj) :: zqoce_b                  ! 
     73      REAL(wp), DIMENSION(jpi,jpj) :: zgammat, zgammas         ! exchange coeficient 
     74      REAL(wp), DIMENSION(jpi,jpj) :: zttbl, zstbl             ! temp. and sal. in top boundary layer 
    7275      !!--------------------------------------------------------------------- 
    7376      ! 
     
    140143      ptsc(:,:,jp_tem) = - zqh(:,:) * r1_rau0_rcp 
    141144      ! 
     145      ! write restart variables (qoceisf, qhcisf, fwfisf for now and before) 
     146      IF (lrst_oce) CALL isfrst_write(kt, 'cav', ptsc, pqfwf) 
     147      ! 
    142148      IF ( ln_isfdebug ) THEN 
    143149         CALL debug('isf_cav: ptsc T',ptsc(:,:,1)) 
     
    158164      !!--------------------------------------------------------------------- 
    159165 
    160       ! allocation isfcav gamtisf, gamsisf,  
     166      ! 0: allocation 
     167      ! 
    161168      CALL isf_alloc_cav() 
    162169      ! 
    163       ! cav 
     170      ! 1: initialisation 
     171      ! 
     172      ! top and bottom level of the 'top boundary layer' 
    164173      misfkt_cav(:,:)    = mikt(:,:) ; misfkb_cav(:,:)    = 1 
     174      ! 
     175      ! thickness of 'tbl' and fraction of bottom cell affected by 'tbl' 
    165176      rhisf_tbl_cav(:,:) = 0.0_wp    ; rfrac_tbl_cav(:,:) = 0.0_wp 
     177      ! 
     178      ! cavity mask 
     179      mskisf_cav(:,:)    = (1._wp - tmask(:,:,1)) * ssmask(:,:) 
     180      ! 
     181      ! 2: read restart 
     182      ! 
     183      ! read cav variable from restart 
     184      IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b) 
     185      ! 
     186      ! 3: specific allocation and initialisation (depending of scheme choice) 
    166187      ! 
    167188      SELECT CASE ( TRIM(cn_isfcav_mlt) ) 
     
    199220      END SELECT 
    200221      ! 
    201       ! compute mask 
    202       mskisf_cav(:,:) = (1._wp - tmask(:,:,1)) * ssmask(:,:) 
    203       ! 
    204222   END SUBROUTINE isf_cav_init 
    205223 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavgam.F90

    r11844 r11852  
    1010   !!   isfcav_gammats       : compute exchange coeficient gamma  
    1111   !!---------------------------------------------------------------------- 
    12    USE oce            ! ocean dynamics and tracers 
    1312   USE isf 
    14    USE isfutils 
    15    USE isftbl 
     13   USE isfutils, ONLY: debug 
     14   USE isftbl  , ONLY: isf_tbl 
     15 
     16   USE oce     , ONLY: un, vn, rn2         ! ocean dynamics and tracers 
     17   USE phycst  , ONLY: grav, vkarmn        ! physical constant 
     18   USE eosbn2  , ONLY: eos_rab             ! equation of state 
     19   USE zdfdrg  , ONLY: rCd0_top, r_ke0_top ! vertical physics: top/bottom drag coef. 
     20   USE iom     , ONLY: iom_put             ! 
     21   USE lib_mpp , ONLY: ctl_stop 
     22 
    1623   USE dom_oce        ! ocean space and time domain 
    17    USE phycst         ! physical constants 
    18    USE eosbn2         ! equation of state 
    19    USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    20    USE iom 
    2124   USE in_out_manager ! I/O manager 
    2225   ! 
     
    4649      !!--------------------------------------------------------------------- 
    4750      !!-------------------------- OUT ------------------------------------- 
    48       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pgt  , pgs      ! gamma t and gamma s  
     51      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pgt  , pgs      ! gamma t and gamma s  
    4952      !!-------------------------- IN  ------------------------------------- 
    5053      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pqoce, pqfwf    ! isf heat and fwf 
     
    114117      !!--------------------------------------------------------------------- 
    115118      !!-------------------------- OUT ------------------------------------- 
    116       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pgt, pgs     ! gammat and gammas [m/s] 
     119      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pgt, pgs     ! gammat and gammas [m/s] 
    117120      !!-------------------------- IN  ------------------------------------- 
    118121      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: putbl, pvtbl ! velocity in the losch top boundary layer 
     
    145148      !!--------------------------------------------------------------------- 
    146149      !!-------------------------- OUT ------------------------------------- 
    147       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pgt, pgs     ! gammat and gammas 
     150      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pgt, pgs     ! gammat and gammas 
    148151      !!-------------------------- IN  ------------------------------------- 
    149152      REAL(wp),                     INTENT(in   ) :: pke2           ! background velocity squared 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavmlt.F90

    r11844 r11852  
    1111   !!   isfcav_mlt    : update surface ocean boundary condition under ice shelf 
    1212   !!---------------------------------------------------------------------- 
    13    USE oce            ! ocean dynamics and tracers 
    14    USE isf            ! ice shelf public variables 
    15    USE isfutils       ! ice shelf debug subroutine 
    16    USE dom_oce        ! ocean space and time domain 
    17    USE phycst         ! physical constants 
    18    USE eosbn2         ! equation of state 
    19    ! 
    20    USE in_out_manager ! I/O manager 
    21    USE iom            ! I/O library 
    22    USE fldread        ! read input field at current time step 
    23    USE lib_fortran 
     13 
     14   USE isf                      ! ice shelf 
     15   USE isftbl , ONLY: isf_tbl   ! ice shelf depth average 
     16   USE isfutils,ONLY: debug     ! debug subroutine 
     17 
     18   USE dom_oce                            ! ocean space and time domain 
     19   USE oce    , ONLY: tsn                 ! ocean dynamics and tracers 
     20   USE phycst , ONLY: rcp, rau0, rau0_rcp ! physical constants 
     21   USE eosbn2 , ONLY: eos_fzp             ! equation of state 
     22 
     23   USE in_out_manager              ! I/O manager 
     24   USE iom        , ONLY: iom_put  ! I/O library 
     25   USE fldread    , ONLY: fld_read ! 
     26   USE lib_fortran, ONLY: glob_sum ! 
     27   USE lib_mpp    , ONLY: ctl_stop ! 
    2428 
    2529   IMPLICIT NONE 
     
    110114      ! 
    111115      ! read input file 
    112       CALL fld_read ( kt, nn_fsbc, sf_isfcav_fwf ) 
     116      CALL fld_read ( kt, 1, sf_isfcav_fwf ) 
    113117      ! 
    114118      ! define fwf and qoce 
     
    276280      ! 
    277281      ! read input file 
    278       CALL fld_read ( kt, nn_fsbc, sf_isfcav_fwf ) 
     282      CALL fld_read ( kt, 1, sf_isfcav_fwf ) 
    279283      ! 
    280284      ! ice shelf 2d map 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcpl.F90

    r11823 r11852  
    1212   !!   isfrst : read/write iceshelf variables in/from restart 
    1313   !!---------------------------------------------------------------------- 
     14   USE isf                             ! ice shelf variable 
     15   USE lib_mpp, ONLY: mpp_sum, mpp_max ! mpp routine 
     16   USE domvvl , ONLY: dom_vvl_zgr      ! vertical scale factor interpolation 
     17   ! 
    1418   USE oce            ! ocean dynamics and tracers 
    15    USE isf            ! ice shelf variable 
    16    USE isfutils       ! debuging 
    17    USE lib_mpp        ! mpp routine 
    18    USE domvvl         ! vertical scale factor interpolation 
    19    ! 
    2019   USE in_out_manager ! I/O manager 
    2120   USE iom            ! I/O library 
     
    4746CONTAINS 
    4847   SUBROUTINE isfcpl_init() 
     48      ! 
     49      ! start on an euler time step 
     50      neuler = 0 
    4951      !  
    5052      CALL isf_alloc_cpl() 
     
    5759      ! 
    5860      ! correction of the horizontal divergence and associated temp. and salt content flux 
     61      ! Need to : - include in the cpl cons the risfcpl_vol/tsc contribution 
     62      !           - decide how to manage thickness level change in conservation 
    5963      CALL isfcpl_vol() 
    6064      ! 
     
    6670      vn(:,:,:) = vn(:,:,:) * vmask(:,:,:) 
    6771      ! 
    68       ! Need to : - include in the cpl cons the risfcpl_vol/tsc contribution 
    69       !           - decide how to manage thickness level change in conservation 
     72      ! all before fields set to now values 
     73      tsb  (:,:,:,:) = tsn  (:,:,:,:) 
     74      ub   (:,:,:)   = un   (:,:,:) 
     75      vb   (:,:,:)   = vn   (:,:,:) 
     76      sshb (:,:)     = sshn (:,:) 
     77      e3t_b(:,:,:)   = e3t_n(:,:,:) 
     78  
     79      ! prepare writing restart 
     80      IF( lwxios ) THEN 
     81         CALL iom_set_rstw_var_active('ssmask') 
     82         CALL iom_set_rstw_var_active('tmask') 
     83         CALL iom_set_rstw_var_active('e3t_n') 
     84         CALL iom_set_rstw_var_active('e3u_n') 
     85         CALL iom_set_rstw_var_active('e3v_n') 
     86      END IF 
    7087      ! 
    7188   END SUBROUTINE isfcpl_init 
     
    92109   END SUBROUTINE isfcpl_rst_write 
    93110 
    94    SUBROUTINE isfcpl_ssh 
     111   SUBROUTINE isfcpl_ssh() 
    95112      !!----------------------------------------------------------------------  
    96113      !!                   ***  ROUTINE iscpl_ssh  *** 
     
    153170      ! 
    154171      ! recompute the vertical scale factor, depth and water thickness 
    155       CALL dom_vvl_zgr 
     172      CALL dom_vvl_zgr() 
    156173      ! 
    157174   END SUBROUTINE isfcpl_ssh 
    158175 
    159    SUBROUTINE isfcpl_tra 
     176   SUBROUTINE isfcpl_tra() 
    160177      !!----------------------------------------------------------------------  
    161178      !!                   ***  ROUTINE iscpl_tra  *** 
     
    316333   END SUBROUTINE isfcpl_tra 
    317334 
    318    SUBROUTINE isfcpl_vol 
     335   SUBROUTINE isfcpl_vol() 
    319336      !!----------------------------------------------------------------------  
    320337      !!                   ***  ROUTINE iscpl_vol  *** 
     
    402419   END SUBROUTINE isfcpl_vol 
    403420 
    404    SUBROUTINE isfcpl_cons 
     421   SUBROUTINE isfcpl_cons() 
    405422      !!----------------------------------------------------------------------  
    406423      !!                   ***  ROUTINE iscpl_cons  *** 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfhdiv.F90

    r11541 r11852  
    11MODULE isfhdiv 
    22 
    3    USE dom_oce 
    4    USE iom 
    5    USE isf 
    6    USE isfutils 
    7    USE phycst 
    8    USE in_out_manager 
     3   USE isf                    ! ice shelf 
     4   USE dom_oce                ! time and space domain 
     5   USE phycst , ONLY: r1_rau0 ! physical constant 
     6   USE in_out_manager         ! 
    97 
    108   IMPLICIT NONE 
     
    3937         ! 
    4038         ! ice sheet coupling contribution  
    41          IF ( ln_isfcpl ) THEN 
     39         IF ( ln_isfcpl .AND. kt /= 0 ) THEN 
    4240            ! 
    4341            ! correct divergence only for the first time step 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfload.F90

    r11823 r11852  
    1111   !!---------------------------------------------------------------------- 
    1212 
    13    USE isf             ! ice shelf variables 
    14    USE dom_oce         ! vertical scale factor 
     13   USE isf, ONLY: cn_isfload   ! ice shelf variables 
     14 
     15   USE dom_oce, ONLY: e3w_n, gdept_n, risfdep, mikt ! vertical scale factor 
     16   USE eosbn2 , ONLY: eos                           ! eos routine 
     17 
     18   USE lib_mpp, ONLY: ctl_stop ! ctl_stop routine 
    1519   USE in_out_manager  !  
    16    USE eosbn2          ! eos routine 
    17    USE lib_mpp         ! ctl_stop routine 
    1820 
    1921   IMPLICIT NONE 
     
    8688      ! 
    8789      !                                !- Surface value + ice shelf gradient 
    88       risfload(:,:) = 0._wp                       ! compute pressure due to ice shelf load  
     90      pisfload(:,:) = 0._wp                       ! compute pressure due to ice shelf load  
    8991      DO jj = 1, jpj                         ! (used to compute hpgi/j for all the level from 1 to miku/v) 
    9092         DO ji = 1, jpi                      ! divided by 2 later 
     
    9496               ! 
    9597               ! top layer of the ice shelf 
    96                risfload(ji,jj) = risfload(ji,jj) + (znad + zrhd(ji,jj,1) ) * e3w_n(ji,jj,1) 
     98               pisfload(ji,jj) = pisfload(ji,jj) + (znad + zrhd(ji,jj,1) ) * e3w_n(ji,jj,1) 
    9799               ! 
    98100               ! core layers of the ice shelf 
    99101               DO jk = 2, ikt-1 
    100                   risfload(ji,jj) = risfload(ji,jj) + (2._wp * znad + zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk)) * e3w_n(ji,jj,jk) 
     102                  pisfload(ji,jj) = pisfload(ji,jj) + (2._wp * znad + zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk)) * e3w_n(ji,jj,jk) 
    101103               END DO 
    102104               ! 
    103105               ! deepest part of the ice shelf (between deepest T point and ice/ocean interface 
    104                risfload(ji,jj) = risfload(ji,jj) + (2._wp * znad + zrhdtop_isf(ji,jj) + zrhd(ji,jj,ikt-1)) & 
     106               pisfload(ji,jj) = pisfload(ji,jj) + (2._wp * znad + zrhdtop_isf(ji,jj) + zrhd(ji,jj,ikt-1)) & 
    105107                  &                                              * ( risfdep(ji,jj) - gdept_n(ji,jj,ikt-1) ) 
    106108               ! 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfnxt.F90

    r11541 r11852  
    88   
    99   !!------------------------------------------------------------------------- 
    10    !!   isfnxt       : aplly correction need for the ice shelf 
     10   !!   isfnxt       : aplly correction needed for the ice shelf to ensure conservation 
    1111   !!------------------------------------------------------------------------- 
    1212 
    1313   USE isf 
    14    USE phycst 
    15    USE dom_oce 
     14 
     15   USE phycst , ONLY: r1_rau0                ! physical constant 
     16   USE dom_oce, ONLY: e3t_b, e3t_n, r1_e1e2t ! time and space domain 
     17 
    1618   USE in_out_manager 
    1719 
     
    2022   PRIVATE 
    2123 
    22    PUBLIC isf_dynnxt ! isf_tranxt 
     24   PUBLIC isf_dynnxt 
    2325 
    2426CONTAINS 
     
    3436      INTEGER ,                     INTENT(in   ) :: kt 
    3537      ! 
    36       REAL(wp),                     INTENT(in   ) :: pcoef           ! atfp * rdt 
     38      REAL(wp),                     INTENT(in   ) :: pcoef           ! atfp * rdt * r1_rau0 
    3739      !!-------------------------- IN  ------------------------------------- 
    3840      !!-------------------------------------------------------------------- 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfpar.F90

    r11541 r11852  
    1414   !!   isfpar       : compute ice shelf melt using a prametrisation of ice shelf cavities 
    1515   !!---------------------------------------------------------------------- 
    16    USE oce            ! ocean dynamics and tracers 
     16   !USE oce            ! ocean dynamics and tracers 
    1717   USE isf            ! ice shelf 
    18    USE isfutils       ! 
    19    USE isfparmlt      ! ice shelf parametrisation 
    20    USE isftbl         ! ice shelf depth average 
    21    USE isfdiags       ! ice shelf diagnostics 
    22    USE dom_oce        ! ocean space and time domain 
    23    USE phycst         ! physical constants 
    24    USE eosbn2         ! equation of state 
     18   ! 
     19   USE isfrst   , ONLY: isfrst_write, isfrst_read ! ice shelf restart read/write subroutine 
     20   USE isftbl   , ONLY: isf_tbl_ktop, isf_tbl_lvl ! ice shelf top boundary layer properties subroutine 
     21   USE isfutils , ONLY: debug, read_2dcstdta      ! ice shelf debug subroutine 
     22   USE isfparmlt, ONLY: isfpar_mlt     ! ice shelf melt formulation subroutine 
     23   USE isfdiags , ONLY: isf_diags_flx  ! ice shelf diags subroutine 
     24   ! 
     25   USE dom_oce  , ONLY: jpi,jpj        ! ocean space and time domain 
     26   USE phycst   , ONLY: r1_rau0_rcp    ! physical constants 
    2527   ! 
    2628   USE in_out_manager ! I/O manager 
    2729   USE iom            ! I/O library 
    2830   USE fldread        ! read input field at current time step 
    29    USE lbclnk         ! 
     31   USE lbclnk         ! lbc_lnk 
    3032 
    3133   IMPLICIT NONE 
     
    8082      ptsc(:,:,jp_tem) = zqh(:,:) * r1_rau0_rcp 
    8183      ! 
     84      ! write restart variables (qoceisf, qhcisf, fwfisf for now and before) 
     85      IF (lrst_oce) CALL isfrst_write(kt, 'par', ptsc, pqfwf) 
     86      ! 
     87      IF ( ln_isfdebug ) THEN 
     88         CALL debug('isf_par: ptsc T',ptsc(:,:,1)) 
     89         CALL debug('isf_par: ptsc S',ptsc(:,:,2)) 
     90         CALL debug('isf_par: pqfwf fwf',pqfwf(:,:)) 
     91      END IF 
     92      ! 
    8293   END SUBROUTINE isf_par 
    8394 
     
    100111      rhisf_tbl_par(:,:)  = 1e-20     ; rfrac_tbl_par(:,:)    = 0.0_wp 
    101112      ! 
    102       mskisf_par(:,:) = 0 
    103       ! 
    104113      ! define isf tbl tickness, top and bottom indice 
    105114      CALL read_2dcstdta(TRIM(sn_isfpar_zmax%clname), TRIM(sn_isfpar_zmax%clvar), ztblmax) 
     
    110119      ztblmin(:,:) = ztblmin(:,:) * ssmask(:,:) 
    111120      ! 
    112       ! if param used under an ice shelf overwrite ztblmax by the ice shelf draft 
     121      ! if param used under an ice shelf overwrite ztblmin by the ice shelf draft 
    113122      WHERE ( risfdep > 0._wp .AND. ztblmin > 0._wp ) 
    114123         ztblmin(:,:) = risfdep(:,:) 
     
    127136      END WHERE 
    128137      ! 
    129       ! compute misfkb_par, rhisf_tbl 
    130       rhisf_tbl_par(:,:) = rhisf0_tbl_par(:,:) 
    131       CALL isf_tbl_lvl( ht_n * mskisf_par, e3t_n, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par ) 
     138      ! read par variable from restart 
     139      IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b) 
    132140      ! 
    133141      SELECT CASE ( TRIM(cn_isfpar_mlt) ) 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfparmlt.F90

    r11541 r11852  
    88   !!---------------------------------------------------------------------- 
    99 
    10    USE oce            ! ocean dynamics and tracers 
    1110   USE isf            ! ice shelf 
    12    USE isftbl         ! ice shelf depth average 
    13    USE dom_oce        ! ocean space and time domain 
    14    USE phycst         ! physical constants 
    15    USE eosbn2         ! equation of state 
    16  
    17    USE in_out_manager ! I/O manager 
    18    USE iom            ! I/O library 
    19    USE fldread        ! 
    20    USE lib_fortran    ! 
     11   USE isftbl , ONLY: isf_tbl   ! ice shelf depth average 
     12 
     13   USE dom_oce                  ! ocean space and time domain 
     14   USE oce    , ONLY: tsn       ! ocean dynamics and tracers 
     15   USE phycst , ONLY: rcp, rau0 ! physical constants 
     16   USE eosbn2 , ONLY: eos_fzp   ! equation of state 
     17 
     18   USE in_out_manager              ! I/O manager 
     19   USE iom        , ONLY: iom_put  ! I/O library 
     20   USE fldread    , ONLY: fld_read ! 
     21   USE lib_fortran, ONLY: glob_sum ! 
     22   USE lib_mpp    , ONLY: ctl_stop ! 
    2123 
    2224   IMPLICIT NONE 
     
    9193      ! 
    9294      ! 0. ------------Read specified runoff 
    93       CALL fld_read ( kt, nn_fsbc, sf_isfpar_fwf   ) 
     95      CALL fld_read ( kt, 1, sf_isfpar_fwf   ) 
    9496      ! 
    9597      ! compute ptfrz 
     
    153155      CALL iom_put('isfthermald_par',( ztfrz(:,:) - ztavg(:,:) ) * mskisf_par(:,:)) 
    154156      ! 
    155       ! 
    156157   END SUBROUTINE isfpar_mlt_bg03 
    157158 
     
    180181      ! 
    181182      ! 0. ------------Read specified runoff 
    182       CALL fld_read ( kt, nn_fsbc, sf_isfpar_fwf   ) 
     183      CALL fld_read ( kt, 1, sf_isfpar_fwf   ) 
    183184      ! 
    184185      ! 1. ------------Mean freezing point (needed for heat content flux) 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfrst.F90

    r11423 r11852  
    1010   !!   isfrst : read/write iceshelf variables in/from restart 
    1111   !!---------------------------------------------------------------------- 
    12    USE oce            ! ocean dynamics and tracers 
    13    USE isf            ! ice shelf variable 
     12   ! 
     13   USE dom_oce, ONLY: jpi,jpj,jpk,jpts ! time and space domain 
    1414   ! 
    1515   USE in_out_manager ! I/O manager 
     
    6161      ! 
    6262      IF( lwxios ) THEN 
    63           CALL iom_set_rstw_var_active(TRIM(chc_b )) 
    64           CALL iom_set_rstw_var_active(TRIM(csc_b )) 
    65           CALL iom_set_rstw_var_active(TRIM(cfwf_b)) 
     63         CALL iom_set_rstw_var_active(TRIM(chc_b )) 
     64         CALL iom_set_rstw_var_active(TRIM(csc_b )) 
     65         CALL iom_set_rstw_var_active(TRIM(cfwf_b)) 
    6666      ENDIF 
     67 
     68      CALL FLUSH(numout) 
    6769 
    6870   END SUBROUTINE isfrst_read 
     
    7476      !!-------------------------- OUT -------------------------------------- 
    7577      !!-------------------------- IN  -------------------------------------- 
    76       INTEGER                     , INTENT(in   ) :: kt 
    77       CHARACTER(LEN=256)          , INTENT(in   ) :: cdisf 
     78      INTEGER                          , INTENT(in   ) :: kt 
     79      CHARACTER(LEN=256)               , INTENT(in   ) :: cdisf 
    7880      REAL(wp), DIMENSION(jpi,jpj)     , INTENT(in   ) :: pfwf 
    7981      REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(in   ) :: ptsc 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90

    r11823 r11852  
    1313   !!   isfstp       : compute iceshelf melt and heat flux 
    1414   !!---------------------------------------------------------------------- 
    15    USE oce            ! ocean dynamics and tracers 
    16    USE dom_oce        ! ocean space and time domain 
    17    USE domvvl, ONLY : ln_vvl_zstar 
    18    USE phycst         ! physical constants 
    19    USE eosbn2         ! equation of state 
    20    USE sbc_oce        ! surface boundary condition: ocean fields 
    21    USE zdfdrg         ! vertical physics: top/bottom drag coef. 
    2215   ! 
     16   USE isf            ! isf variables 
     17   USE isfload, ONLY: isf_load                      ! ice shelf load 
     18   USE isftbl , ONLY: isf_tbl_lvl                   ! ice shelf boundary layer 
     19   USE isfpar , ONLY: isf_par, isf_par_init         ! ice shelf parametrisation 
     20   USE isfcav , ONLY: isf_cav, isf_cav_init         ! ice shelf cavity 
     21   USE isfcpl , ONLY: isfcpl_rst_write, isfcpl_init ! isf variables 
     22 
     23   USE dom_oce, ONLY: ht_n, e3t_n, ln_isfcav, ln_linssh ! ocean space and time domain 
     24   USE domvvl,  ONLY: ln_vvl_zstar                      ! zstar logical 
     25   USE zdfdrg,  ONLY: r_Cdmin_top, r_ke0_top            ! vertical physics: top/bottom drag coef. 
     26   ! 
     27   USE lib_mpp, ONLY: ctl_stop, ctl_nam 
    2328   USE in_out_manager ! I/O manager 
    24    USE iom            ! I/O library 
    25    USE fldread        ! read input field at current time step 
    26    USE lbclnk         ! 
    27    USE lib_fortran    ! glob_sum 
    28    ! 
    29    USE isfrst         ! iceshelf restart 
    30    USE isftbl         ! ice shelf boundary layer 
    31    USE isfpar         ! ice shelf parametrisation 
    32    USE isfcav         ! ice shelf cavity 
    33    USE isfload        ! ice shelf load 
    34    USE isfcpl         ! isf variables 
    35    USE isf            ! isf variables 
    36    USE isfutils 
    3729 
    3830   IMPLICIT NONE 
     
    8274         CALL isf_cav( kt, risf_cav_tsc, fwfisf_cav) 
    8375         ! 
    84          ! write restart variables (risf_cav_tsc, fwfisf for now and before) 
    85          IF (lrst_oce) CALL isfrst_write(kt, 'cav', risf_cav_tsc, fwfisf_cav) 
    86          ! 
    8776      END IF 
    8877      !  
     
    10291         CALL isf_par( kt, risf_par_tsc, fwfisf_par) 
    10392         ! 
    104          ! write restart variables (qoceisf, qhcisf, fwfisf for now and before) 
    105          IF (lrst_oce) CALL isfrst_write(kt, 'par', risf_par_tsc, fwfisf_par) 
    106          ! 
    107       END IF 
    108  
    109       IF ( ln_isfcpl ) THEN 
    110          !  
    111          IF (lrst_oce) CALL isfcpl_rst_write(kt) 
    112          ! 
    113       END IF 
     93      END IF 
     94      ! 
     95      IF ( ln_isfcpl .AND. lrst_oce ) CALL isfcpl_rst_write(kt) 
    11496      ! 
    11597   END SUBROUTINE isf_stp 
     
    148130         ! 
    149131         !--------------------------------------------------------------------------------------------------------------------- 
    150          ! initialisation ice sheet coupling 
    151          IF( ln_isfcpl ) THEN  
    152             ! 
    153             ! start on an euler time step 
    154             neuler = 0 
    155             ! 
    156             ! ice sheet coupling: extrapolation of restart to fill new wet cell and compute divergence correction 
    157             CALL isfcpl_init() 
    158             ! 
    159             ! all before fields set to now values 
    160             tsb  (:,:,:,:) = tsn  (:,:,:,:) 
    161             ub   (:,:,:)   = un   (:,:,:) 
    162             vb   (:,:,:)   = vn   (:,:,:) 
    163             sshb (:,:)     = sshn (:,:) 
    164             e3t_b(:,:,:)   = e3t_n(:,:,:) 
    165   
    166             ! prepare writing restart 
    167             IF( lwxios ) CALL iom_set_rstw_var_active('ssmask') 
    168             IF( lwxios ) CALL iom_set_rstw_var_active('tmask') 
    169             IF( lwxios ) CALL iom_set_rstw_var_active('e3t_n') 
    170             IF( lwxios ) CALL iom_set_rstw_var_active('e3u_n') 
    171             IF( lwxios ) CALL iom_set_rstw_var_active('e3v_n') 
    172             ! 
    173          END IF 
    174          ! 
    175          !--------------------------------------------------------------------------------------------------------------------- 
    176132         ! initialisation melt in the cavity 
    177          IF ( ln_isfcav_mlt ) THEN 
    178             ! 
    179             ! initialisation  of cav variable 
    180             CALL isf_cav_init() 
    181             ! 
    182             ! read cav variable from restart 
    183             IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b) 
    184             ! 
    185          END IF 
     133         IF ( ln_isfcav_mlt ) CALL isf_cav_init() 
    186134         ! 
    187135         !--------------------------------------------------------------------------------------------------------------------- 
    188136         ! initialisation parametrised melt 
    189          IF ( ln_isfpar_mlt ) THEN 
    190             ! 
    191             ! initialisation  of par variable 
    192             CALL isf_par_init() 
    193             ! 
    194             ! read par variable from restart 
    195             IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b) 
    196             ! 
    197          END IF 
     137         IF ( ln_isfpar_mlt ) CALL isf_par_init() 
     138         ! 
     139         !--------------------------------------------------------------------------------------------------------------------- 
     140         ! initialisation ice sheet coupling 
     141         IF( ln_isfcpl ) CALL isfcpl_init() 
    198142         ! 
    199143      END IF 
     
    224168               WRITE(numout,*) '         gamma formulation                       cn_gammablk = ', TRIM(cn_gammablk)  
    225169               IF ( TRIM(cn_gammablk) .NE. 'spe' ) THEN  
    226                   WRITE(numout,*) '            gammat coefficient                      rn_gammat0  = ', rn_gammat0   
    227                   WRITE(numout,*) '            gammas coefficient                      rn_gammas0  = ', rn_gammas0   
    228                   WRITE(numout,*) '            top drag coef. used (from namdrg_top)   rn_Cd0      = ', r_Cdmin_top  
     170                  WRITE(numout,*) '            gammat coefficient                       rn_gammat0  = ', rn_gammat0   
     171                  WRITE(numout,*) '            gammas coefficient                       rn_gammas0  = ', rn_gammas0   
     172                  WRITE(numout,*) '            top drag coef.    used (from namdrg_top) rn_Cd0      = ', r_Cdmin_top 
     173                  WRITE(numout,*) '            top background ke used (from namdrg_top) rn_ke0      = ', r_ke0_top 
    229174               END IF 
    230175            END IF 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isftbl.F90

    r11844 r11852  
    1414   !!---------------------------------------------------------------------- 
    1515 
     16   USE isf     ! ice shelf variables 
     17 
    1618   USE dom_oce ! vertical scale factor 
    17    USE lbclnk  ! lbc_lnk subroutine 
    18    USE isfutils 
    19    USE isf 
    2019 
    2120   IMPLICIT NONE 
     
    6968         CALL isf_tbl_avg( miku, ikbot, zhtbl, zfrac, e3u_n, pvarin, zvarout ) 
    7069         ! 
    71          ! check if needed (probably yes) 
    72          CALL lbc_lnk('sbcisf', pvarout,'U',-1.) 
    73          ! 
    7470         ! compute tbl property at T point 
    7571         pvarout(1,:) = 0._wp 
     
    9288         CALL isf_tbl_avg( mikv, ikbot, zhtbl, zfrac, e3v_n, pvarin, zvarout ) 
    9389         ! 
    94          ! check if needed (probably yes) 
    95          CALL lbc_lnk('sbcisf', pvarout,'V',-1.) 
    96          ! 
    9790         ! pvarout is an averaging of wet point 
    9891         pvarout(:,1) = 0._wp 
     
    110103         ! 
    111104      END SELECT 
    112       ! 
    113       IF (ln_isfdebug) CALL debug('isf_tbl pvarout:',pvarout) 
    114105      ! 
    115106   END SUBROUTINE isf_tbl 
  • NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfutils.F90

    r11844 r11852  
    1212   !!---------------------------------------------------------------------- 
    1313 
    14    USE iom 
    15    USE lib_fortran 
    16    USE lib_mpp 
     14   USE iom           , ONLY: iom_open, iom_get, iom_close, jpdom_data ! read input file 
     15   USE lib_fortran   , ONLY: glob_sum, glob_min, glob_max             ! compute global value 
     16   USE dom_oce       , ONLY: jpi,jpj,jpk                              ! domain size 
     17   USE in_out_manager, ONLY: wp, lwp, numout                          ! miscelenious 
    1718 
    1819   IMPLICIT NONE 
Note: See TracChangeset for help on using the changeset viewer.