source: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90 @ 11403

MODULE isfstp
   !!======================================================================
   !!                       ***  MODULE  isfstp  ***
   !! Surface module :  compute iceshelf load, melt and heat flux
   !!======================================================================
   !! History :  3.2  !  2011-02  (C.Harris  ) Original code isf cav
   !!            X.X  !  2006-02  (C. Wang   ) Original code bg03
   !!            3.4  !  2013-03  (P. Mathiot) Merging + parametrization
   !!            4.1  !  2019-09  (P. Mathiot) Split param/explicit ice shelf and re-organisation
   !!----------------------------------------------------------------------

   !!----------------------------------------------------------------------
   !!   isfstp       : compute iceshelf melt and heat flux
   !!----------------------------------------------------------------------
   USE oce            ! ocean dynamics and tracers
   USE dom_oce        ! ocean space and time domain
   USE phycst         ! physical constants
   USE eosbn2         ! equation of state
   USE sbc_oce        ! surface boundary condition: ocean fields
   USE zdfdrg         ! vertical physics: top/bottom drag coef.
   !
   USE in_out_manager ! I/O manager
   USE iom            ! I/O library
   USE fldread        ! read input field at current time step
   USE lbclnk         !
   USE lib_fortran    ! glob_sum
   !
   USE isfrst         ! iceshelf restart
   USE isftbl         ! ice shelf boundary layer
   USE isfpar         ! ice shelf parametrisation
   USE isfcav         ! ice shelf cavity
   USE isfload        ! ice shelf load
   USE isf            ! isf variables

   IMPLICIT NONE

   PRIVATE

   PUBLIC   isf_stp, isf_stp_init  ! routine called in sbcmod and divhor

   !!----------------------------------------------------------------------
   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
   !! $Id$
   !! Software governed by the CeCILL license (see ./LICENSE)
   !!----------------------------------------------------------------------
CONTAINS
 
  SUBROUTINE isf_stp( kt )
      !!---------------------------------------------------------------------
      !!                  ***  ROUTINE isf_stp  ***
      !!
      !! ** Purpose : compute total heat flux and total fwf due to ice shelf melt
      !!
      !! ** Method  : For each case (parametrisation or explicity cavity) :
      !!              - define the before fields
      !!              - compute top boundary layer properties
      !!                (in case of parametrisation, this is the 
      !!                 depth range model array used to compute mean far fields properties)
      !!              - compute fluxes
      !!              - write restart variables
      !!
      !!----------------------------------------------------------------------
      INTEGER, INTENT(in) ::   kt   ! ocean time step
      !!---------------------------------------------------------------------
      !
      IF ( ln_isfcav_mlt ) THEN
         !
         ! before time step 
         IF ( kt /= nit000 ) THEN 
            risf_cav_tsc_b (:,:,:) = risf_cav_tsc (:,:,:)
            fwfisf_cav_b(:,:)      = fwfisf_cav(:,:)
         END IF
         !
         ! compute tbl lvl/h
         CALL isf_tbl_lvl(ht_n, e3t_n, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav)
         !
         ! compute ice shelf melt
         CALL isf_cav( kt, risf_cav_tsc, fwfisf_cav)
         !
         ! write restart variables (qoceisf, qhcisf, fwfisf for now and before)
         IF (lrst_oce) CALL isfrst_write(kt, 'cav', risf_cav_tsc, fwfisf_cav)
         !
      END IF
      ! 
      IF ( ln_isfpar_mlt ) THEN
         !
         ! before time step 
         IF ( kt /= nit000 ) THEN 
            risf_par_tsc_b (:,:) = risf_par_tsc (:,:)
            fwfisf_par_b(:,:)    = fwfisf_par(:,:)
         END IF
         !
         ! compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl)
         CALL isf_tbl_lvl(ht_n, e3t_n, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par)
         !
         ! compute ice shelf melt
         CALL isf_par( kt, risf_par_tsc, fwfisf_par)
         !
         ! write restart variables (qoceisf, qhcisf, fwfisf for now and before)
         IF (lrst_oce) CALL isfrst_write(kt, 'par', risf_par_tsc, fwfisf_par)
         !
      END IF
      !
   END SUBROUTINE isf_stp

   SUBROUTINE isf_stp_init
      !!---------------------------------------------------------------------
      !!                  ***  ROUTINE isfstp_init  ***
      !!
      !! ** Purpose :   Initialisation of the ice shelf public variables
      !!
      !! ** Method  :   Read the namsbc namelist and set derived parameters
      !!                Call init routines for all other SBC modules that have one
      !!
      !! ** Action  : - read namsbc parameters
      !!              - allocate memory
      !!              - call cav/param init routine
      !!----------------------------------------------------------------------
      INTEGER               :: inum, ierror
      INTEGER               :: ios                  ! Local integer output status for namelist read
      INTEGER               :: ikt, ikb
      INTEGER               :: ji, jj
      !!----------------------------------------------------------------------
      NAMELIST/namisf/ ln_isfcav_mlt, cn_isfcav_mlt, cn_gammablk, rn_gammat0, rn_gammas0, rn_htbl, sn_isfcav_fwf,  &
         &             ln_isfpar_mlt, cn_isfpar_mlt, sn_isfpar_fwf, sn_isfpar_zmin, sn_isfpar_zmax, sn_isfpar_Leff
      !!----------------------------------------------------------------------
      !
      ! Allocate public array
      CALL isf_alloc()
      !
      riceload(:,:)       = 0.0_wp
      fwfisf_cpl(:,:)     = 0.0_wp
      fwfisf_par(:,:)     = 0.0_wp    ; fwfisf_par_b(:,:)     = 0.0_wp
      fwfisf_cav(:,:)     = 0.0_wp    ; fwfisf_cav_b(:,:)     = 0.0_wp
      risf_cav_tsc(:,:,:) = 0.0_wp    ; risf_cav_tsc_b(:,:,:) = 0.0_wp
      risf_par_tsc(:,:,:) = 0.0_wp    ; risf_par_tsc_b(:,:,:) = 0.0_wp
      !
      ! terminate routine now if no ice shelf melt formulation specify
      IF ( .NOT. ln_isf ) RETURN
      !
      REWIND( numnam_ref )              ! Namelist namsbc_rnf in reference namelist : Runoffs 
      READ  ( numnam_ref, namisf, IOSTAT = ios, ERR = 901)
901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namisf in reference namelist', lwp )
      !
      REWIND( numnam_cfg )              ! Namelist namsbc_rnf in configuration namelist : Runoffs
      READ  ( numnam_cfg, namisf, IOSTAT = ios, ERR = 902 )
902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namisf in configuration namelist', lwp )
      IF(lwm) WRITE ( numond, namisf )
      !
      IF(lwp) WRITE(numout,*)
      IF(lwp) WRITE(numout,*) 'isf_init : ice shelf initialisation'
      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
      IF(lwp) WRITE(numout,*) '   Namelist namisf :'
      !
      IF(lwp) WRITE(numout,*) '      melt inside the cavity                  ln_isfcav_mlt   = ', ln_isfcav_mlt
      IF ( ln_isfcav ) THEN
         IF(lwp) WRITE(numout,*) '         melt formulation                        cn_isfcav_mlt   = ', cn_isfcav_mlt
         IF(lwp) WRITE(numout,*) '         thickness of the top boundary layer     rn_htbl     = ', rn_htbl
         IF(lwp) WRITE(numout,*) '         gamma formulation                       cn_gammablk = ', cn_gammablk 
         IF ( TRIM(cn_gammablk) .NE. 'spe' ) THEN 
            IF(lwp) WRITE(numout,*) '            gammat coefficient                      rn_gammat0  = ', rn_gammat0  
            IF(lwp) WRITE(numout,*) '            gammas coefficient                      rn_gammas0  = ', rn_gammas0  
            IF(lwp) WRITE(numout,*) '            top drag coef. used (from namdrg_top)   rn_Cd0      = ', r_Cdmin_top 
         END IF
      END IF
      !
      IF(lwp) WRITE(numout,*) ''
      !
      IF(lwp) WRITE(numout,*) '      ice shelf melt parametrisation          ln_isfpar_mlt    = ', ln_isfpar_mlt
      IF ( ln_isfpar_mlt ) THEN
         IF(lwp) WRITE(numout,*) '         isf parametrisation formulation         cn_isfpar_mlt   = ', cn_isfpar_mlt
      END IF
      IF(lwp) WRITE(numout,*) ''
      !
      ! sanity check
      ! melt in the cavity without cavity
      IF ( ln_isfcav_mlt .AND. (.NOT. ln_isfcav) ) &
         &   CALL ctl_stop('ice shelf melt in the cavity activated (ln_isfcav_mlt) but no cavity detected in domcfg (ln_isfcav), STOP' )
      !
      IF ( ln_cpl ) THEN
         !
         CALL ctl_stop( ' ln_ctl and ice shelf not tested' )
         !
         ! NEMO coupled to ATMO model with isf cavity need oasis method for melt computation 
         IF ( ln_isfcav_mlt .AND. TRIM(cn_isfcav_mlt) /= 'oasis' ) CALL ctl_stop( 'cn_isfcav_mlt = oasis is the only option availble with ln_cpl' )
         IF ( ln_isfpar_mlt .AND. TRIM(cn_isfpar_mlt) /= 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis is the only option availble with ln_cpl' )
         !
         ! oasis melt computation not tested (coded but not tested)
         IF ( ln_isfcav_mlt .OR. ln_isfpar_mlt ) THEN
            IF ( TRIM(cn_isfcav_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfcav_mlt = oasis not tested' )
            IF ( TRIM(cn_isfpar_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis not tested' )
         END IF
         !
         ! oasis melt computation with cavity open and cavity parametrised (not coded)
         IF ( ln_isfcav_mlt .AND. ln_isfpar_mlt ) THEN
            IF ( TRIM(cn_isfpar_mlt) == 'oasis' .AND. TRIM(cn_isfcav_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis and cn_isfcav_mlt = oasis not coded' )
         END IF
      END IF
      !
      ! initialisation ice shelf load
      IF ( ln_isfcav ) THEN
         !
         ! compute ice shelf mask
         mskisf_cav(:,:) = (1._wp - tmask(:,:,1)) * ssmask(:,:)
         !
         ! compute ice shelf load
         CALL isf_load( risfload )
         !
      END IF
      !
      ! initialisation useful variable
      r1_Lfusisf =  1._wp / rLfusisf
      !
      ! initialisation melt in the cavity
      IF ( ln_isfcav_mlt ) THEN
         !
         ! initialisation  of cav variable
         CALL isf_cav_init()
         !
         ! read cav variable from restart
         IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b)
         !
      END IF
      !
      ! initialisation parametrised melt
      IF ( ln_isfpar_mlt ) THEN
         !
         ! initialisation  of par variable
         CALL isf_par_init()
         !
         ! read par variable from restart
         IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b)
         !
      END IF
      !
  END SUBROUTINE isf_stp_init
   !!======================================================================
END MODULE isfstp