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 6488 for branches/UKMO/dev_r5518_GO6_package/NEMOGCM – NEMO

Ignore:
Timestamp:
2016-04-20T11:42:09+02:00 (8 years ago)
Author:
davestorkey
Message:

Commit changes from dev_r5518_coupling_GSI7_GSI8_landice and its ancestor branch dev_r5518_CICE_coupling_GSI7_GSI8.
Custom merge into /branches/UKMO/dev_r5518_GO6_package/NEMOGCM: r6023 cf. r5668 of /branches/UKMO/dev_r5518_coupling_GSI7_GSI8_landice/NEMOGCM@6487

Custom merge into /branches/UKMO/dev_r5518_GO6_package/NEMOGCM: r5668 cf. r5662 of /branches/UKMO/dev_r5518_CICE_coupling_GSI7_GSI8/NEMOGCM@6487

Location:
branches/UKMO/dev_r5518_GO6_package/NEMOGCM
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/CONFIG/SHARED/field_def.xml

    r5517 r6488  
    193193 
    194194         <!-- * variable related to ice shelf forcing * --> 
     195         <field id="berg_calve"   long_name="Iceberg calving"                               unit="kg/m2/s"  /> 
    195196         <field id="fwfisf"       long_name="Ice shelf melting"                             unit="kg/m2/s"  /> 
    196197         <field id="qisf"         long_name="Ice Shelf Heat Flux"                           unit="W/m2"     /> 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/CONFIG/SHARED/namelist_ref

    r6487 r6488  
    378378   ln_usecplmask = .false.   !  use a coupling mask file to merge data received from several models 
    379379                             !   -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel) 
     380   ln_coupled_iceshelf_fluxes = .false. ! If true use rate of change of mass of Greenland and Antarctic icesheets to set the  
     381                                        ! combined magnitude of the iceberg calving and iceshelf melting freshwater fluxes. 
     382   rn_greenland_calving_fraction = 0.5  ! Set fraction of total freshwater flux for iceberg calving - remainder goes to iceshelf melting. 
     383   rn_antarctica_calving_fraction = 0.5 ! Set fraction of total freshwater flux for iceberg calving - remainder goes to iceshelf melting. 
     384   rn_iceshelf_fluxes_tolerance = 1e-6  ! Fractional threshold for detecting differences in icesheet masses (must be positive definite). 
    380385/ 
    381386!----------------------------------------------------------------------- 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90

    r6486 r6488  
    6868   REAL(wp), PUBLIC ::   rcdsn    =    0.31_wp       !: thermal conductivity of snow 
    6969   REAL(wp), PUBLIC ::   cpic     = 2067.0_wp        !: specific heat for ice  
     70#if defined key_cice 
     71   REAL(wp), PUBLIC ::   lsub     =    2.835e+6_wp   !: pure ice latent heat of sublimation                   [J/kg] 
     72#else 
    7073   REAL(wp), PUBLIC ::   lsub     =    2.834e+6_wp   !: pure ice latent heat of sublimation                   [J/kg] 
     74#endif 
    7175   REAL(wp), PUBLIC ::   lfus     =    0.334e+6_wp   !: latent heat of fusion of fresh ice                    [J/kg] 
    7276   REAL(wp), PUBLIC ::   tmut     =    0.054_wp      !: decrease of seawater meltpoint with salinity 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/ICB/icbclv.F90

    r6486 r6488  
    2525   USE icbutl         ! iceberg utility routines 
    2626 
     27   USE sbc_oce        ! for icesheet freshwater input variables 
     28   USE in_out_manager 
     29   USE iom 
     30 
    2731   IMPLICIT NONE 
    2832   PRIVATE 
     
    4852      ! 
    4953      REAL(wp)                        :: zcalving_used, zdist, zfact 
     54      REAL(wp)                        :: zgreenland_calving_sum, zantarctica_calving_sum 
    5055      INTEGER                         :: jn, ji, jj                    ! loop counters 
    5156      INTEGER                         :: imx                           ! temporary integer for max berg class 
     
    5964      zfact = ( (1000._wp)**3 / ( NINT(rday) * nyear_len(1) ) ) * 850._wp 
    6065      berg_grid%calving(:,:) = src_calving(:,:) * tmask_i(:,:) * zfact 
     66 
     67      IF( lk_oasis) THEN 
     68      ! ln_coupled_iceshelf_fluxes uninitialised unless lk_oasis=true 
     69      IF( ln_coupled_iceshelf_fluxes ) THEN 
     70 
     71        ! Adjust total calving rates so that sum of iceberg calving and iceshelf melting in the northern 
     72        ! and southern hemispheres equals rate of increase of mass of greenland and antarctic ice sheets 
     73        ! to preserve total freshwater conservation in coupled models without an active ice sheet model. 
     74 
     75         zgreenland_calving_sum = SUM( berg_grid%calving(:,:) * greenland_icesheet_mask(:,:) ) 
     76         IF( lk_mpp ) CALL mpp_sum( zgreenland_calving_sum ) 
     77         WHERE( greenland_icesheet_mask(:,:) == 1.0 )                                                                                 & 
     78        &    berg_grid%calving(:,:) = berg_grid%calving(:,:) * greenland_icesheet_mass_rate_of_change * rn_greenland_calving_fraction & 
     79        &                                     / ( zgreenland_calving_sum + 1.0e-10_wp ) 
     80 
     81         ! check 
     82         IF(lwp) WRITE(numout, *) 'Greenland iceberg calving climatology (kg/s) : ',zgreenland_calving_sum 
     83         zgreenland_calving_sum = SUM( berg_grid%calving(:,:) * greenland_icesheet_mask(:,:) ) 
     84         IF( lk_mpp ) CALL mpp_sum( zgreenland_calving_sum ) 
     85         IF(lwp) WRITE(numout, *) 'Greenland iceberg calving adjusted value (kg/s) : ',zgreenland_calving_sum 
     86 
     87         zantarctica_calving_sum = SUM( berg_grid%calving(:,:) * antarctica_icesheet_mask(:,:) ) 
     88         IF( lk_mpp ) CALL mpp_sum( zantarctica_calving_sum ) 
     89         WHERE( antarctica_icesheet_mask(:,:) == 1.0 )                                                                              & 
     90         berg_grid%calving(:,:) = berg_grid%calving(:,:) * antarctica_icesheet_mass_rate_of_change * rn_antarctica_calving_fraction & 
     91        &                           / ( zantarctica_calving_sum + 1.0e-10_wp ) 
     92  
     93         ! check 
     94         IF(lwp) WRITE(numout, *) 'Antarctica iceberg calving climatology (kg/s) : ',zantarctica_calving_sum 
     95         zantarctica_calving_sum = SUM( berg_grid%calving(:,:) * antarctica_icesheet_mask(:,:) ) 
     96         IF( lk_mpp ) CALL mpp_sum( zantarctica_calving_sum ) 
     97         IF(lwp) WRITE(numout, *) 'Antarctica iceberg calving adjusted value (kg/s) : ',zantarctica_calving_sum 
     98 
     99      ENDIF 
     100      ENDIF 
     101    
     102      CALL iom_put( 'berg_calve', berg_grid%calving(:,:) ) 
    61103 
    62104      ! Heat in units of W/m2, and mask (just in case) 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/IOM/restart.F90

    r6486 r6488  
    2424   USE trdmxl_oce      ! ocean active mixed layer tracers trends variables 
    2525   USE divcur          ! hor. divergence and curl      (div & cur routines) 
     26   USE sbc_oce         ! for icesheet freshwater input variables 
    2627 
    2728   IMPLICIT NONE 
     
    145146                     CALL iom_rstput( kt, nitrst, numrow, 'rhd'    , rhd       ) 
    146147#endif 
     148                     IF( lk_oasis) THEN 
     149                     ! ln_coupled_iceshelf_fluxes uninitialised unless lk_oasis=true 
     150                     IF( ln_coupled_iceshelf_fluxes ) THEN 
     151                        CALL iom_rstput( kt, nitrst, numrow, 'greenland_icesheet_mass', greenland_icesheet_mass ) 
     152                        CALL iom_rstput( kt, nitrst, numrow, 'greenland_icesheet_timelapsed', greenland_icesheet_timelapsed ) 
     153                        CALL iom_rstput( kt, nitrst, numrow, 'greenland_icesheet_mass_roc', greenland_icesheet_mass_rate_of_change ) 
     154                        CALL iom_rstput( kt, nitrst, numrow, 'antarctica_icesheet_mass', antarctica_icesheet_mass ) 
     155                        CALL iom_rstput( kt, nitrst, numrow, 'antarctica_icesheet_timelapsed', antarctica_icesheet_timelapsed ) 
     156                        CALL iom_rstput( kt, nitrst, numrow, 'antarctica_icesheet_mass_roc', antarctica_icesheet_mass_rate_of_change ) 
     157                     ENDIF 
     158                     ENDIF 
     159 
    147160      IF( kt == nitrst ) THEN 
    148161         CALL iom_close( numrow )     ! close the restart file (only at last time step) 
     
    258271#endif 
    259272      ! 
     273      IF( iom_varid( numror, 'greenland_icesheet_mass', ldstop = .FALSE. ) > 0 )   THEN 
     274         CALL iom_get( numror, 'greenland_icesheet_mass', greenland_icesheet_mass ) 
     275         CALL iom_get( numror, 'greenland_icesheet_timelapsed', greenland_icesheet_timelapsed ) 
     276         CALL iom_get( numror, 'greenland_icesheet_mass_roc', greenland_icesheet_mass_rate_of_change ) 
     277      ELSE 
     278         greenland_icesheet_mass = 0.0  
     279         greenland_icesheet_mass_rate_of_change = 0.0  
     280         greenland_icesheet_timelapsed = 0.0 
     281      ENDIF 
     282      IF( iom_varid( numror, 'antarctica_icesheet_mass', ldstop = .FALSE. ) > 0 )   THEN 
     283         CALL iom_get( numror, 'antarctica_icesheet_mass', antarctica_icesheet_mass ) 
     284         CALL iom_get( numror, 'antarctica_icesheet_timelapsed', antarctica_icesheet_timelapsed ) 
     285         CALL iom_get( numror, 'antarctica_icesheet_mass_roc', antarctica_icesheet_mass_rate_of_change ) 
     286      ELSE 
     287         antarctica_icesheet_mass = 0.0  
     288         antarctica_icesheet_mass_rate_of_change = 0.0  
     289         antarctica_icesheet_timelapsed = 0.0 
     290      ENDIF 
     291 
    260292      IF( neuler == 0 ) THEN                                  ! Euler restart (neuler=0) 
    261293         tsb  (:,:,:,:) = tsn  (:,:,:,:)                             ! all before fields set to now values 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90

    r6486 r6488  
    101101   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr_iu              !: ice fraction at NEMO U point 
    102102   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   fr_iv              !: ice fraction at NEMO V point 
    103     
     103   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   sstfrz             !: sea surface freezing temperature 
     104   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tsfc_ice           !: sea-ice surface skin temperature (on categories) 
     105   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   kn_ice             !: sea-ice surface layer thermal conductivity (on cats) 
     106 
    104107   ! variables used in the coupled interface 
    105108   INTEGER , PUBLIC, PARAMETER ::   jpl = ncat 
    106109   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   u_ice, v_ice          ! jpi, jpj 
     110   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_p, ht_p ! Meltpond fraction and depth 
    107111#endif 
    108112    
     
    152156 
    153157#if defined key_cice 
    154       ALLOCATE( qla_ice(jpi,jpj,1)    , qlw_ice(jpi,jpj,1)    , qsr_ice(jpi,jpj,1)    , & 
     158      ALLOCATE( qla_ice(jpi,jpj,ncat) , qlw_ice(jpi,jpj,1)    , qsr_ice(jpi,jpj,1)    , & 
    155159                wndi_ice(jpi,jpj)     , tatm_ice(jpi,jpj)     , qatm_ice(jpi,jpj)     , & 
    156160                wndj_ice(jpi,jpj)     , nfrzmlt(jpi,jpj)      , ss_iou(jpi,jpj)       , & 
    157161                ss_iov(jpi,jpj)       , fr_iu(jpi,jpj)        , fr_iv(jpi,jpj)        , & 
    158162                a_i(jpi,jpj,ncat)     , topmelt(jpi,jpj,ncat) , botmelt(jpi,jpj,ncat) , & 
    159                 STAT= ierr(1) ) 
    160       IF( ln_cpl )   ALLOCATE( u_ice(jpi,jpj)        , fr1_i0(jpi,jpj)       , tn_ice (jpi,jpj,1)    , & 
     163                sstfrz(jpi,jpj)       , STAT= ierr(1) ) 
     164   ! Alex West: Allocating tn_ice with 5 categories.  When NEMO is used with CICE, this variable 
     165   ! represents top layer ice temperature, which is multi-category. 
     166      IF( ln_cpl )   ALLOCATE( u_ice(jpi,jpj)        , fr1_i0(jpi,jpj)       , tn_ice (jpi,jpj,jpl)  , & 
    161167         &                     v_ice(jpi,jpj)        , fr2_i0(jpi,jpj)       , alb_ice(jpi,jpj,1)    , & 
    162168         &                     emp_ice(jpi,jpj)      , qns_ice(jpi,jpj,1)    , dqns_ice(jpi,jpj,1)   , & 
    163          &                     STAT= ierr(2) ) 
     169         &                     a_p(jpi,jpj,jpl)      , ht_p(jpi,jpj,jpl)     , tsfc_ice(jpi,jpj,jpl) , & 
     170         &                     kn_ice(jpi,jpj,jpl) ,    STAT=ierr(2) ) 
    164171       
    165172#endif 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90

    r6486 r6488  
    125125#endif 
    126126   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xcplmask          !: coupling mask for ln_mixcpl (warning: allocated in sbccpl) 
     127   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   greenland_icesheet_mass_array, greenland_icesheet_mask 
     128   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   antarctica_icesheet_mass_array, antarctica_icesheet_mask 
    127129 
    128130   !!---------------------------------------------------------------------- 
     
    137139   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   e3t_m     !: mean (nn_fsbc time-step) sea surface layer thickness       [m] 
    138140   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   frq_m     !: mean (nn_fsbc time-step) fraction of solar net radiation absorbed in the 1st T level [-] 
     141    
     142   !!---------------------------------------------------------------------- 
     143   !!  Surface scalars of total ice sheet mass for Greenland and Antarctica,  
     144   !! passed from atmosphere to be converted to dvol and hence a freshwater  
     145   !! flux  by using old values. New values are saved in the dump, to become 
     146   !! old values next coupling timestep. Freshwater fluxes split between  
     147   !! sub iceshelf melting and iceberg calving, scalled to flux per second 
     148   !!---------------------------------------------------------------------- 
     149    
     150   REAL(wp), PUBLIC  :: greenland_icesheet_mass, greenland_icesheet_mass_rate_of_change, greenland_icesheet_timelapsed  
     151   REAL(wp), PUBLIC  :: antarctica_icesheet_mass, antarctica_icesheet_mass_rate_of_change, antarctica_icesheet_timelapsed 
     152 
     153   ! sbccpl namelist parameters associated with icesheet freshwater input code. Included here rather than in sbccpl.F90 to  
     154   ! avoid circular dependencies. 
     155   LOGICAL, PUBLIC     ::   ln_coupled_iceshelf_fluxes     ! If true use rate of change of mass of Greenland and Antarctic icesheets to set the  
     156                                                           ! combined magnitude of the iceberg calving and iceshelf melting freshwater fluxes. 
     157   REAL(wp), PUBLIC    ::   rn_greenland_calving_fraction  ! Set fraction of total freshwater flux for iceberg calving - remainder goes to iceshelf melting. 
     158   REAL(wp), PUBLIC    ::   rn_antarctica_calving_fraction ! Set fraction of total freshwater flux for iceberg calving - remainder goes to iceshelf melting. 
     159   REAL(wp), PUBLIC    ::   rn_iceshelf_fluxes_tolerance   ! Absolute tolerance for detecting differences in icesheet masses.  
    139160 
    140161   !! * Substitutions 
     
    172193         &      ssu_m  (jpi,jpj) , sst_m(jpi,jpj) , frq_m(jpi,jpj) ,      & 
    173194         &      ssv_m  (jpi,jpj) , sss_m(jpi,jpj) , ssh_m(jpi,jpj) , STAT=ierr(4) ) 
     195      ALLOCATE( greenland_icesheet_mass_array(jpi,jpj) , antarctica_icesheet_mass_array(jpi,jpj) ) 
     196      ALLOCATE( greenland_icesheet_mask(jpi,jpj) , antarctica_icesheet_mask(jpi,jpj) ) 
    174197         ! 
    175198#if defined key_vvl 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r6487 r6488  
    4646   USE p4zflx, ONLY : oce_co2 
    4747#endif 
    48 #if defined key_cice 
    49    USE ice_domain_size, only: ncat 
    50 #endif 
    5148#if defined key_lim3 
    5249   USE limthd_dh       ! for CALL lim_thd_snwblow 
     
    105102   INTEGER, PARAMETER ::   jpr_e3t1st = 41            ! first T level thickness  
    106103   INTEGER, PARAMETER ::   jpr_fraqsr = 42            ! fraction of solar net radiation absorbed in the first ocean level 
    107    INTEGER, PARAMETER ::   jprcv      = 42            ! total number of fields received 
     104   INTEGER, PARAMETER ::   jpr_ts_ice = 43            ! skin temperature of sea-ice (used for melt-ponds) 
     105   INTEGER, PARAMETER ::   jpr_grnm   = 44            ! Greenland ice mass 
     106   INTEGER, PARAMETER ::   jpr_antm   = 45            ! Antarctic ice mass 
     107   INTEGER, PARAMETER ::   jprcv      = 45            ! total number of fields received 
    108108 
    109109   INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere 
     
    135135   INTEGER, PARAMETER ::   jps_e3t1st = 27            ! first level depth (vvl) 
    136136   INTEGER, PARAMETER ::   jps_fraqsr = 28            ! fraction of solar net radiation absorbed in the first ocean level 
    137    INTEGER, PARAMETER ::   jpsnd      = 28            ! total number of fields sended 
     137   INTEGER, PARAMETER ::   jps_a_p    = 29            ! meltpond fraction   
     138   INTEGER, PARAMETER ::   jps_ht_p   = 30            ! meltpond depth (m)  
     139   INTEGER, PARAMETER ::   jps_kice   = 31            ! ice surface layer thermal conductivity 
     140   INTEGER, PARAMETER ::   jps_sstfrz = 32            ! sea-surface freezing temperature 
     141   INTEGER, PARAMETER ::   jps_fice1  = 33            ! first-order ice concentration (for time-travelling ice coupling) 
     142   INTEGER, PARAMETER ::   jpsnd      = 33            ! total number of fields sended 
    138143 
    139144   !                                                         !!** namelist namsbc_cpl ** 
     
    146151   END TYPE FLD_C 
    147152   ! Send to the atmosphere                           ! 
    148    TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2                         
     153   TYPE(FLD_C) ::   sn_snd_temp, sn_snd_alb, sn_snd_thick, sn_snd_crt, sn_snd_co2, sn_snd_cond, sn_snd_mpnd, sn_snd_sstfrz, sn_snd_thick1 
     154 
    149155   ! Received from the atmosphere                     ! 
    150156   TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf 
    151    TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2                         
     157   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_ts_ice, sn_rcv_grnm, sn_rcv_antm 
    152158   ! Other namelist parameters                        ! 
    153159   INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
     
    216222      REAL(wp), POINTER, DIMENSION(:,:) ::   zacs, zaos 
    217223      !! 
    218       NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2,      & 
    219          &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr,      & 
    220          &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal   , sn_rcv_iceflx,   & 
    221          &                  sn_rcv_co2 , nn_cplmodel  , ln_usecplmask 
     224      NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick , sn_snd_crt   , sn_snd_co2,     & 
     225         &                  sn_snd_cond, sn_snd_mpnd  , sn_snd_sstfrz, sn_snd_thick1,                 & 
     226         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau   , sn_rcv_dqnsdt, sn_rcv_qsr,     & 
     227         &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf   , sn_rcv_cal   , sn_rcv_iceflx,  & 
     228         &                  sn_rcv_co2 , sn_rcv_grnm  , sn_rcv_antm  , sn_rcv_ts_ice, nn_cplmodel  ,  & 
     229         &                  ln_usecplmask, ln_coupled_iceshelf_fluxes, rn_greenland_calving_fraction, & 
     230         &                  rn_antarctica_calving_fraction, rn_iceshelf_fluxes_tolerance 
    222231      !!--------------------------------------------------------------------- 
    223232      ! 
     
    258267         WRITE(numout,*)'      runoffs                         = ', TRIM(sn_rcv_rnf%cldes   ), ' (', TRIM(sn_rcv_rnf%clcat   ), ')' 
    259268         WRITE(numout,*)'      calving                         = ', TRIM(sn_rcv_cal%cldes   ), ' (', TRIM(sn_rcv_cal%clcat   ), ')' 
     269         WRITE(numout,*)'      Greenland ice mass              = ', TRIM(sn_rcv_grnm%cldes  ), ' (', TRIM(sn_rcv_grnm%clcat  ), ')' 
     270         WRITE(numout,*)'      Antarctica ice mass             = ', TRIM(sn_rcv_antm%cldes  ), ' (', TRIM(sn_rcv_antm%clcat  ), ')' 
    260271         WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 
    261272         WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')' 
     
    269280         WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd 
    270281         WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')' 
     282         WRITE(numout,*)'      ice effective conductivity      = ', TRIM(sn_snd_cond%cldes   ), ' (', TRIM(sn_snd_cond%clcat   ), ')' 
     283         WRITE(numout,*)'      meltponds fraction & depth      = ', TRIM(sn_snd_mpnd%cldes  ), ' (', TRIM(sn_snd_mpnd%clcat   ), ')' 
     284         WRITE(numout,*)'      sea surface freezing temp       = ', TRIM(sn_snd_sstfrz%cldes   ), ' (', TRIM(sn_snd_sstfrz%clcat   ), ')' 
     285 
    271286         WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel 
    272287         WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask 
     288         WRITE(numout,*)'  ln_coupled_iceshelf_fluxes          = ', ln_coupled_iceshelf_fluxes 
     289         WRITE(numout,*)'  rn_greenland_calving_fraction       = ', rn_greenland_calving_fraction 
     290         WRITE(numout,*)'  rn_antarctica_calving_fraction      = ', rn_antarctica_calving_fraction 
     291         WRITE(numout,*)'  rn_iceshelf_fluxes_tolerance        = ', rn_iceshelf_fluxes_tolerance 
    273292      ENDIF 
    274293 
     
    383402      srcv(jpr_snow)%clname = 'OTotSnow'      ! Snow = solid precipitation 
    384403      srcv(jpr_tevp)%clname = 'OTotEvap'      ! total evaporation (over oce + ice sublimation) 
    385       srcv(jpr_ievp)%clname = 'OIceEvap'      ! evaporation over ice = sublimation 
     404      srcv(jpr_ievp)%clname = 'OIceEvp'      ! evaporation over ice = sublimation 
    386405      srcv(jpr_sbpr)%clname = 'OSubMPre'      ! sublimation - liquid precipitation - solid precipitation  
    387406      srcv(jpr_semp)%clname = 'OISubMSn'      ! ice solid water budget = sublimation - solid precipitation 
     
    396415      CASE default              ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_emp%cldes' ) 
    397416      END SELECT 
    398  
     417      !Set the number of categories for coupling of sublimation 
     418      IF ( TRIM( sn_rcv_emp%clcat ) == 'yes' ) srcv(jpr_ievp)%nct = jpl 
     419      ! 
    399420      !                                                      ! ------------------------- ! 
    400421      !                                                      !     Runoffs & Calving     !    
     
    410431      ! 
    411432      srcv(jpr_cal   )%clname = 'OCalving'   ;   IF( TRIM( sn_rcv_cal%cldes ) == 'coupled' )   srcv(jpr_cal)%laction = .TRUE. 
     433      srcv(jpr_grnm  )%clname = 'OGrnmass'   ;   IF( TRIM( sn_rcv_grnm%cldes ) == 'coupled' )   srcv(jpr_grnm)%laction = .TRUE. 
     434      srcv(jpr_antm  )%clname = 'OAntmass'   ;   IF( TRIM( sn_rcv_antm%cldes ) == 'coupled' )   srcv(jpr_antm)%laction = .TRUE. 
     435 
    412436 
    413437      !                                                      ! ------------------------- ! 
     
    483507         srcv(jpr_topm:jpr_botm)%laction = .TRUE. 
    484508      ENDIF 
     509       
     510#if defined key_cice && ! defined key_cice4 
     511      !                                                      ! ----------------------------- ! 
     512      !                                                      !  sea-ice skin temperature     !    
     513      !                                                      !  used in meltpond scheme      ! 
     514      !                                                      !  May be calculated in Atm     ! 
     515      !                                                      ! ----------------------------- ! 
     516      srcv(jpr_ts_ice)%clname = 'OTsfIce' 
     517      IF ( TRIM( sn_rcv_ts_ice%cldes ) == 'ice' ) srcv(jpr_ts_ice)%laction = .TRUE. 
     518      IF ( TRIM( sn_rcv_ts_ice%clcat ) == 'yes' ) srcv(jpr_ts_ice)%nct = jpl 
     519      !TODO: Should there be a consistency check here? 
     520#endif 
     521 
    485522      !                                                      ! ------------------------------- ! 
    486523      !                                                      !   OPA-SAS coupling - rcv by opa !    
     
    600637      !                                                      ! ------------------------- ! 
    601638      ssnd(jps_toce)%clname = 'O_SSTSST' 
    602       ssnd(jps_tice)%clname = 'O_TepIce' 
     639      ssnd(jps_tice)%clname = 'OTepIce' 
    603640      ssnd(jps_tmix)%clname = 'O_TepMix' 
    604641      SELECT CASE( TRIM( sn_snd_temp%cldes ) ) 
    605642      CASE( 'none'                                 )       ! nothing to do 
    606643      CASE( 'oce only'                             )   ;   ssnd( jps_toce )%laction = .TRUE. 
    607       CASE( 'oce and ice' , 'weighted oce and ice' ) 
     644      CASE( 'oce and ice' , 'weighted oce and ice' , 'oce and weighted ice') 
    608645         ssnd( (/jps_toce, jps_tice/) )%laction = .TRUE. 
    609646         IF ( TRIM( sn_snd_temp%clcat ) == 'yes' )  ssnd(jps_tice)%nct = jpl 
     
    634671 
    635672      !                                                      ! ------------------------- ! 
    636       !                                                      !  Ice fraction & Thickness !  
     673      !                                                      !  Ice fraction & Thickness  
    637674      !                                                      ! ------------------------- ! 
    638675      ssnd(jps_fice)%clname = 'OIceFrc' 
    639676      ssnd(jps_hice)%clname = 'OIceTck' 
    640677      ssnd(jps_hsnw)%clname = 'OSnwTck' 
     678      ssnd(jps_a_p)%clname  = 'OPndFrc' 
     679      ssnd(jps_ht_p)%clname = 'OPndTck' 
     680      ssnd(jps_fice1)%clname = 'OIceFrd' 
    641681      IF( k_ice /= 0 ) THEN 
    642682         ssnd(jps_fice)%laction = .TRUE.                  ! if ice treated in the ocean (even in climato case) 
     683         ssnd(jps_fice1)%laction = .TRUE.                 ! First-order regridded ice concentration, to be used 
     684                                                     ! in producing atmos-to-ice fluxes 
    643685! Currently no namelist entry to determine sending of multi-category ice fraction so use the thickness entry for now 
    644686         IF ( TRIM( sn_snd_thick%clcat ) == 'yes' ) ssnd(jps_fice)%nct = jpl 
     687         IF ( TRIM( sn_snd_thick1%clcat ) == 'yes' ) ssnd(jps_fice1)%nct = jpl 
    645688      ENDIF 
    646689       
     
    657700      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_thick%cldes' ) 
    658701      END SELECT 
     702 
     703      !                                                      ! ------------------------- ! 
     704      !                                                      ! Ice Meltponds             ! 
     705      !                                                      ! ------------------------- ! 
     706#if defined key_cice && ! defined key_cice4 
     707      ! Meltponds only CICE5  
     708      ssnd(jps_a_p)%clname = 'OPndFrc'    
     709      ssnd(jps_ht_p)%clname = 'OPndTck'    
     710      SELECT CASE ( TRIM( sn_snd_mpnd%cldes ) ) 
     711      CASE ( 'none' ) 
     712         ssnd(jps_a_p)%laction = .FALSE. 
     713         ssnd(jps_ht_p)%laction = .FALSE. 
     714      CASE ( 'ice only' )  
     715         ssnd(jps_a_p)%laction = .TRUE. 
     716         ssnd(jps_ht_p)%laction = .TRUE. 
     717         IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN 
     718            ssnd(jps_a_p)%nct = jpl 
     719            ssnd(jps_ht_p)%nct = jpl 
     720         ELSE 
     721            IF ( jpl > 1 ) THEN 
     722               CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_mpnd%cldes if not exchanging category fields' ) 
     723            ENDIF 
     724         ENDIF 
     725      CASE ( 'weighted ice' )  
     726         ssnd(jps_a_p)%laction = .TRUE. 
     727         ssnd(jps_ht_p)%laction = .TRUE. 
     728         IF ( TRIM( sn_snd_mpnd%clcat ) == 'yes' ) THEN 
     729            ssnd(jps_a_p)%nct = jpl  
     730            ssnd(jps_ht_p)%nct = jpl  
     731         ENDIF 
     732      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_mpnd%cldes' ) 
     733      END SELECT 
     734#else 
     735      IF( TRIM( sn_snd_mpnd%cldes /= 'none' ) THEN 
     736         CALL ctl_stop('Meltponds can only be used with CICEv5') 
     737      ENDIF 
     738#endif 
    659739 
    660740      !                                                      ! ------------------------- ! 
     
    689769      !                                                      ! ------------------------- ! 
    690770      ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE. 
     771      ! 
     772       
     773      !                                                      ! ------------------------- ! 
     774      !                                                      ! Sea surface freezing temp ! 
     775      !                                                      ! ------------------------- ! 
     776      ssnd(jps_sstfrz)%clname = 'O_SSTFrz' ; IF( TRIM(sn_snd_sstfrz%cldes) == 'coupled' )  ssnd(jps_sstfrz)%laction = .TRUE. 
     777      ! 
     778      !                                                      ! ------------------------- ! 
     779      !                                                      !    Ice conductivity       ! 
     780      !                                                      ! ------------------------- ! 
     781      ! Note that ultimately we will move to passing an ocean effective conductivity as well so there 
     782      ! will be some changes to the parts of the code which currently relate only to ice conductivity 
     783      ssnd(jps_kice )%clname = 'OIceKn' 
     784      SELECT CASE ( TRIM( sn_snd_cond%cldes ) ) 
     785      CASE ( 'none' ) 
     786         ssnd(jps_kice)%laction = .FALSE. 
     787      CASE ( 'ice only' ) 
     788         ssnd(jps_kice)%laction = .TRUE. 
     789         IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) THEN 
     790            ssnd(jps_kice)%nct = jpl 
     791         ELSE 
     792            IF ( jpl > 1 ) THEN 
     793               CALL ctl_stop( 'sbc_cpl_init: use weighted ice option for sn_snd_cond%cldes if not exchanging category fields' ) 
     794            ENDIF 
     795         ENDIF 
     796      CASE ( 'weighted ice' ) 
     797         ssnd(jps_kice)%laction = .TRUE. 
     798         IF ( TRIM( sn_snd_cond%clcat ) == 'yes' ) ssnd(jps_kice)%nct = jpl 
     799      CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_cond%cldes' ) 
     800      END SELECT 
     801      ! 
     802       
    691803 
    692804      !                                                      ! ------------------------------- ! 
     
    785897      ncpl_qsr_freq = 86400 / ncpl_qsr_freq 
    786898 
     899      IF( ln_coupled_iceshelf_fluxes ) THEN 
     900          ! Crude masks to separate the Antarctic and Greenland icesheets. Obviously something 
     901          ! more complicated could be done if required. 
     902          greenland_icesheet_mask = 0.0 
     903          WHERE( gphit >= 0.0 ) greenland_icesheet_mask = 1.0 
     904          antarctica_icesheet_mask = 0.0 
     905          WHERE( gphit < 0.0 ) antarctica_icesheet_mask = 1.0 
     906 
     907          ! initialise other variables 
     908          greenland_icesheet_mass_array(:,:) = 0.0 
     909          antarctica_icesheet_mass_array(:,:) = 0.0 
     910 
     911          IF( .not. ln_rstart ) THEN 
     912             greenland_icesheet_mass = 0.0  
     913             greenland_icesheet_mass_rate_of_change = 0.0  
     914             greenland_icesheet_timelapsed = 0.0 
     915             antarctica_icesheet_mass = 0.0  
     916             antarctica_icesheet_mass_rate_of_change = 0.0  
     917             antarctica_icesheet_timelapsed = 0.0 
     918          ENDIF 
     919 
     920      ENDIF 
     921 
    787922      CALL wrk_dealloc( jpi,jpj, zacs, zaos ) 
    788923      ! 
     
    843978      !! 
    844979      LOGICAL  ::   llnewtx, llnewtau      ! update wind stress components and module?? 
    845       INTEGER  ::   ji, jj, jn             ! dummy loop indices 
     980      INTEGER  ::   ji, jj, jl, jn         ! dummy loop indices 
    846981      INTEGER  ::   isec                   ! number of seconds since nit000 (assuming rdttra did not change since nit000) 
    847982      REAL(wp) ::   zcumulneg, zcumulpos   ! temporary scalars      
     983      REAL(wp) ::   zgreenland_icesheet_mass_in, zantarctica_icesheet_mass_in 
     984      REAL(wp) ::   zgreenland_icesheet_mass_b, zantarctica_icesheet_mass_b 
     985      REAL(wp) ::   zmask_sum, zepsilon       
    848986      REAL(wp) ::   zcoef                  ! temporary scalar 
    849987      REAL(wp) ::   zrhoa  = 1.22          ! Air density kg/m3 
     
    9971135#endif 
    9981136 
     1137#if defined key_cice && ! defined key_cice4 
     1138      !  ! Sea ice surface skin temp: 
     1139      IF( srcv(jpr_ts_ice)%laction ) THEN 
     1140        DO jl = 1, jpl 
     1141          DO jj = 1, jpj 
     1142            DO ji = 1, jpi 
     1143              IF (frcv(jpr_ts_ice)%z3(ji,jj,jl) > 0.0) THEN 
     1144                tsfc_ice(ji,jj,jl) = 0.0 
     1145              ELSE IF (frcv(jpr_ts_ice)%z3(ji,jj,jl) < -60.0) THEN 
     1146                tsfc_ice(ji,jj,jl) = -60.0 
     1147              ELSE 
     1148                tsfc_ice(ji,jj,jl) = frcv(jpr_ts_ice)%z3(ji,jj,jl) 
     1149              ENDIF 
     1150            END DO 
     1151          END DO 
     1152        END DO 
     1153      ENDIF 
     1154#endif 
     1155 
    9991156      !  Fields received by SAS when OASIS coupling 
    10001157      !  (arrays no more filled at sbcssm stage) 
     
    11121269 
    11131270      ENDIF 
     1271       
     1272      !                                                        ! land ice masses : Greenland 
     1273      zepsilon = rn_iceshelf_fluxes_tolerance 
     1274 
     1275      IF( srcv(jpr_grnm)%laction ) THEN 
     1276         greenland_icesheet_mass_array(:,:) = frcv(jpr_grnm)%z3(:,:,1) 
     1277         ! take average over ocean points of input array to avoid cumulative error over time 
     1278         zgreenland_icesheet_mass_in = SUM( greenland_icesheet_mass_array(:,:) * tmask(:,:,1) ) 
     1279         IF(lk_mpp) CALL mpp_sum( zgreenland_icesheet_mass_in ) 
     1280         zmask_sum = SUM( tmask(:,:,1) ) 
     1281         IF(lk_mpp) CALL mpp_sum( zmask_sum )  
     1282         zgreenland_icesheet_mass_in = zgreenland_icesheet_mass_in / zmask_sum 
     1283         greenland_icesheet_timelapsed = greenland_icesheet_timelapsed + rdt          
     1284         IF( ABS( zgreenland_icesheet_mass_in - greenland_icesheet_mass ) > zepsilon ) THEN 
     1285            zgreenland_icesheet_mass_b = greenland_icesheet_mass 
     1286             
     1287            ! Only update the mass if it has increased 
     1288            IF ( (zgreenland_icesheet_mass_in - greenland_icesheet_mass) > 0.0 ) THEN 
     1289               greenland_icesheet_mass = zgreenland_icesheet_mass_in 
     1290            ENDIF 
     1291             
     1292            IF( zgreenland_icesheet_mass_b /= 0.0 ) & 
     1293           &     greenland_icesheet_mass_rate_of_change = ( greenland_icesheet_mass - zgreenland_icesheet_mass_b ) / greenland_icesheet_timelapsed  
     1294            greenland_icesheet_timelapsed = 0.0_wp        
     1295         ENDIF 
     1296         IF(lwp) WRITE(numout,*) 'Greenland icesheet mass (kg) read in is ', zgreenland_icesheet_mass_in 
     1297         IF(lwp) WRITE(numout,*) 'Greenland icesheet mass (kg) used is    ', greenland_icesheet_mass 
     1298         IF(lwp) WRITE(numout,*) 'Greenland icesheet mass rate of change (kg/s) is ', greenland_icesheet_mass_rate_of_change 
     1299         IF(lwp) WRITE(numout,*) 'Greenland icesheet seconds lapsed since last change is ', greenland_icesheet_timelapsed 
     1300      ENDIF 
     1301 
     1302      !                                                        ! land ice masses : Antarctica 
     1303      IF( srcv(jpr_antm)%laction ) THEN 
     1304         antarctica_icesheet_mass_array(:,:) = frcv(jpr_antm)%z3(:,:,1) 
     1305         ! take average over ocean points of input array to avoid cumulative error from rounding errors over time 
     1306         zantarctica_icesheet_mass_in = SUM( antarctica_icesheet_mass_array(:,:) * tmask(:,:,1) ) 
     1307         IF(lk_mpp) CALL mpp_sum( zantarctica_icesheet_mass_in ) 
     1308         zmask_sum = SUM( tmask(:,:,1) ) 
     1309         IF(lk_mpp) CALL mpp_sum( zmask_sum )  
     1310         zantarctica_icesheet_mass_in = zantarctica_icesheet_mass_in / zmask_sum 
     1311         antarctica_icesheet_timelapsed = antarctica_icesheet_timelapsed + rdt          
     1312         IF( ABS( zantarctica_icesheet_mass_in - antarctica_icesheet_mass ) > zepsilon ) THEN 
     1313            zantarctica_icesheet_mass_b = antarctica_icesheet_mass 
     1314             
     1315            ! Only update the mass if it has increased 
     1316            IF ( (zantarctica_icesheet_mass_in - antarctica_icesheet_mass) > 0.0 ) THEN 
     1317               antarctica_icesheet_mass = zantarctica_icesheet_mass_in 
     1318            END IF 
     1319             
     1320            IF( zantarctica_icesheet_mass_b /= 0.0 ) & 
     1321          &      antarctica_icesheet_mass_rate_of_change = ( antarctica_icesheet_mass - zantarctica_icesheet_mass_b ) / antarctica_icesheet_timelapsed  
     1322            antarctica_icesheet_timelapsed = 0.0_wp        
     1323         ENDIF 
     1324         IF(lwp) WRITE(numout,*) 'Antarctica icesheet mass (kg) read in is ', zantarctica_icesheet_mass_in 
     1325         IF(lwp) WRITE(numout,*) 'Antarctica icesheet mass (kg) used is    ', antarctica_icesheet_mass 
     1326         IF(lwp) WRITE(numout,*) 'Antarctica icesheet mass rate of change (kg/s) is ', antarctica_icesheet_mass_rate_of_change 
     1327         IF(lwp) WRITE(numout,*) 'Antarctica icesheet seconds lapsed since last change is ', antarctica_icesheet_timelapsed 
     1328      ENDIF 
     1329 
    11141330      ! 
    11151331      CALL wrk_dealloc( jpi,jpj, ztx, zty, zmsk, zemp, zqns, zqsr ) 
     
    14051621         zsprecip(:,:) = frcv(jpr_snow)%z3(:,:,1)                  ! May need to ensure positive here 
    14061622         ztprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:)  ! May need to ensure positive here 
    1407          zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:) 
     1623         zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:)          
     1624#if defined key_cice 
     1625         IF ( TRIM(sn_rcv_emp%clcat) == 'yes' ) THEN 
     1626            ! zemp_ice is the sum of frcv(jpr_ievp)%z3(:,:,1) over all layers - snow 
     1627            zemp_ice(:,:) = - frcv(jpr_snow)%z3(:,:,1) 
     1628            DO jl=1,jpl 
     1629               zemp_ice(:,:   ) = zemp_ice(:,:) + frcv(jpr_ievp)%z3(:,:,jl) 
     1630            ENDDO 
     1631            ! latent heat coupled for each category in CICE 
     1632            qla_ice(:,:,1:jpl) = - frcv(jpr_ievp)%z3(:,:,1:jpl) * lsub 
     1633         ELSE 
     1634            ! If CICE has multicategories it still expects coupling fields for 
     1635            ! each even if we treat as a single field 
     1636            ! The latent heat flux is split between the ice categories according 
     1637            ! to the fraction of the ice in each category 
     1638            zemp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) 
     1639            WHERE ( zicefr(:,:) /= 0._wp )  
     1640               ztmp(:,:) = 1./zicefr(:,:) 
     1641            ELSEWHERE  
     1642               ztmp(:,:) = 0.e0 
     1643            END WHERE   
     1644            DO jl=1,jpl 
     1645               qla_ice(:,:,jl) = - a_i(:,:,jl) * ztmp(:,:) * frcv(jpr_ievp)%z3(:,:,1) * lsub  
     1646            END DO 
     1647            WHERE ( zicefr(:,:) == 0._wp )  qla_ice(:,:,1) = -frcv(jpr_ievp)%z3(:,:,1) * lsub  
     1648         ENDIF 
     1649#else          
    14081650         zemp_ice(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) 
     1651#endif                   
    14091652            CALL iom_put( 'rain'         , frcv(jpr_rain)%z3(:,:,1)              )   ! liquid precipitation  
    14101653         IF( iom_use('hflx_rain_cea') )   & 
     
    17602003               CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' ) 
    17612004               END SELECT 
     2005            CASE( 'oce and weighted ice' )   ;   ztmp1(:,:) =   tsn(:,:,1,jp_tem) + rt0  
     2006               SELECT CASE( sn_snd_temp%clcat ) 
     2007               CASE( 'yes' )    
     2008           ztmp3(:,:,1:jpl) = tn_ice(:,:,1:jpl) * a_i(:,:,1:jpl) 
     2009               CASE( 'no' ) 
     2010           ztmp3(:,:,:) = 0.0 
     2011           DO jl=1,jpl 
     2012                     ztmp3(:,:,1) = ztmp3(:,:,1) + tn_ice(:,:,jl) * a_i(:,:,jl) 
     2013           ENDDO 
     2014               CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' ) 
     2015               END SELECT 
    17622016            CASE( 'mixed oce-ice'        )    
    17632017               ztmp1(:,:) = ( ztmp1(:,:) + rt0 ) * zfr_l(:,:)  
     
    18302084         END SELECT 
    18312085         IF( ssnd(jps_fice)%laction )   CALL cpl_snd( jps_fice, isec, ztmp3, info ) 
     2086      ENDIF 
     2087       
     2088      ! Send ice fraction field (first order interpolation), for weighting UM fluxes to be passed to NEMO 
     2089      IF (ssnd(jps_fice1)%laction) THEN 
     2090         SELECT CASE (sn_snd_thick1%clcat) 
     2091         CASE( 'yes' )   ;   ztmp3(:,:,1:jpl) =  a_i(:,:,1:jpl) 
     2092         CASE( 'no' )    ;   ztmp3(:,:,1) = fr_i(:,:) 
     2093         CASE default    ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick1%clcat' ) 
     2094    END SELECT 
     2095         CALL cpl_snd (jps_fice1, isec, ztmp3, info) 
    18322096      ENDIF 
    18332097       
     
    18752139         IF( ssnd(jps_hsnw)%laction )   CALL cpl_snd( jps_hsnw, isec, ztmp4, info ) 
    18762140      ENDIF 
     2141      ! 
     2142      ! Send meltpond fields  
     2143      IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN 
     2144         SELECT CASE( sn_snd_mpnd%cldes)  
     2145         CASE( 'weighted ice' )  
     2146            SELECT CASE( sn_snd_mpnd%clcat )  
     2147            CASE( 'yes' )  
     2148               ztmp3(:,:,1:jpl) =  a_p(:,:,1:jpl) * a_i(:,:,1:jpl)  
     2149               ztmp4(:,:,1:jpl) =  ht_p(:,:,1:jpl) * a_i(:,:,1:jpl)  
     2150            CASE( 'no' )  
     2151               ztmp3(:,:,:) = 0.0  
     2152               ztmp4(:,:,:) = 0.0  
     2153               DO jl=1,jpl  
     2154                 ztmp3(:,:,1) = ztmp3(:,:,1) + a_p(:,:,jpl) * a_i(:,:,jpl)  
     2155                 ztmp4(:,:,1) = ztmp4(:,:,1) + ht_p(:,:,jpl) * a_i(:,:,jpl)  
     2156               ENDDO  
     2157            CASE default    ;   CALL ctl_stop( 'sbc_cpl_mpd: wrong definition of sn_snd_mpnd%clcat' )  
     2158            END SELECT  
     2159         CASE( 'ice only' )     
     2160            ztmp3(:,:,1:jpl) = a_p(:,:,1:jpl)  
     2161            ztmp4(:,:,1:jpl) = ht_p(:,:,1:jpl)  
     2162         END SELECT  
     2163         IF( ssnd(jps_a_p)%laction )   CALL cpl_snd( jps_a_p, isec, ztmp3, info )     
     2164         IF( ssnd(jps_ht_p)%laction )   CALL cpl_snd( jps_ht_p, isec, ztmp4, info )     
     2165         ! 
     2166         ! Send ice effective conductivity 
     2167         SELECT CASE( sn_snd_cond%cldes) 
     2168         CASE( 'weighted ice' )    
     2169            SELECT CASE( sn_snd_cond%clcat ) 
     2170            CASE( 'yes' )    
     2171               ztmp3(:,:,1:jpl) =  kn_ice(:,:,1:jpl) * a_i(:,:,1:jpl) 
     2172            CASE( 'no' ) 
     2173               ztmp3(:,:,:) = 0.0 
     2174               DO jl=1,jpl 
     2175                 ztmp3(:,:,1) = ztmp3(:,:,1) + kn_ice(:,:,jl) * a_i(:,:,jl) 
     2176               ENDDO 
     2177            CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_cond%clcat' ) 
     2178            END SELECT 
     2179         CASE( 'ice only' )    
     2180           ztmp3(:,:,1:jpl) = kn_ice(:,:,1:jpl) 
     2181         END SELECT 
     2182         IF( ssnd(jps_kice)%laction )   CALL cpl_snd( jps_kice, isec, ztmp3, info ) 
     2183      ENDIF 
     2184      ! 
    18772185      ! 
    18782186#if defined key_cpl_carbon_cycle 
     
    20542362      IF( ssnd(jps_rnf   )%laction )  CALL cpl_snd( jps_rnf   , isec, RESHAPE ( rnf , (/jpi,jpj,1/) ), info ) 
    20552363      IF( ssnd(jps_taum  )%laction )  CALL cpl_snd( jps_taum  , isec, RESHAPE ( taum, (/jpi,jpj,1/) ), info ) 
    2056  
     2364       
     2365      ztmp1(:,:) = sstfrz(:,:) + rt0 
     2366      IF( ssnd(jps_sstfrz)%laction )  CALL cpl_snd( jps_sstfrz, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info ) 
     2367      ! 
    20572368      CALL wrk_dealloc( jpi,jpj, zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 ) 
    20582369      CALL wrk_dealloc( jpi,jpj,jpl, ztmp3, ztmp4 ) 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90

    r6486 r6488  
    1515   USE dom_oce         ! ocean space and time domain 
    1616   USE domvvl 
    17    USE phycst, only : rcp, rau0, r1_rau0, rhosn, rhoic 
     17   USE eosbn2, only : eos_fzp ! Function to calculate freezing point of seawater 
     18   USE phycst, only : rcp, rau0, r1_rau0, rhosn, rhoic, rt0 
    1819   USE in_out_manager  ! I/O manager 
    1920   USE iom, ONLY : iom_put,iom_use              ! I/O manager library !!Joakim edit 
     
    3738   USE ice_gather_scatter 
    3839   USE ice_calendar, only: dt 
     40# if defined key_cice4 
    3941   USE ice_state, only: aice,aicen,uvel,vvel,vsno,vsnon,vice,vicen 
    40 # if defined key_cice4 
    4142   USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  & 
    4243                strocnxT,strocnyT,                               &  
     
    4546                flatn_f,fsurfn_f,fcondtopn_f,                    & 
    4647                uatm,vatm,wind,fsw,flw,Tair,potT,Qa,rhoa,zlvl,   & 
    47                 swvdr,swvdf,swidr,swidf 
     48                swvdr,swvdf,swidr,swidf,Tf 
    4849   USE ice_therm_vertical, only: calc_Tsfc 
    4950#else 
     51   USE ice_state, only: aice,aicen,uvel,nt_hpnd,trcrn,vvel,vsno,& 
     52                vsnon,vice,vicen,nt_Tsfc 
    5053   USE ice_flux, only: strax,stray,strocnx,strocny,frain,fsnow,  & 
    5154                strocnxT,strocnyT,                               &  
    52                 sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_ai,     & 
    53                 fresh_ai,fhocn_ai,fswthru_ai,frzmlt,          & 
     55                sst,sss,uocn,vocn,ss_tltx,ss_tlty,fsalt_ai,      & 
     56                fresh_ai,fhocn_ai,fswthru_ai,frzmlt,             & 
    5457                flatn_f,fsurfn_f,fcondtopn_f,                    & 
    5558                uatm,vatm,wind,fsw,flw,Tair,potT,Qa,rhoa,zlvl,   & 
    56                 swvdr,swvdf,swidr,swidf 
    57    USE ice_therm_shared, only: calc_Tsfc 
     59                swvdr,swvdf,swidr,swidf,Tf,                      & 
     60      !! When using NEMO with CICE, this change requires use of  
     61      !! one of the following two CICE branches: 
     62      !! - at CICE5.0,   hadax/r1015_GSI8_with_GSI7 
     63      !! - at CICE5.1.2, hadax/vn5.1.2_GSI8 
     64                keffn_top,Tn_top 
     65 
     66   USE ice_therm_shared, only: calc_Tsfc, heat_capacity 
     67   USE ice_shortwave, only: apeffn 
    5868#endif 
    5969   USE ice_forcing, only: frcvdr,frcvdf,frcidr,frcidf 
     
    161171      INTEGER, INTENT( in  ) ::   ksbc                ! surface forcing type 
    162172      REAL(wp), DIMENSION(:,:), POINTER :: ztmp1, ztmp2 
    163       REAL(wp) ::   zcoefu, zcoefv, zcoeff            ! local scalar 
     173      REAL(wp), DIMENSION(:,:,:), POINTER :: ztfrz3d 
    164174      INTEGER  ::   ji, jj, jl, jk                    ! dummy loop indices 
    165175      !!--------------------------------------------------------------------- 
     
    174184      jj_off = INT ( (jpjglo - ny_global) / 2 ) 
    175185 
    176 #if defined key_nemocice_decomp 
    177       ! Pass initial SST from NEMO to CICE so ice is initialised correctly if 
    178       ! there is no restart file. 
    179       ! Values from a CICE restart file would overwrite this 
    180       IF ( .NOT. ln_rstart ) THEN     
    181          CALL nemo2cice( tsn(:,:,1,jp_tem) , sst , 'T' , 1.)  
    182       ENDIF   
    183 #endif 
    184  
    185 ! Initialize CICE 
     186      ! Initialize CICE 
    186187      CALL CICE_Initialize 
    187188 
    188 ! Do some CICE consistency checks 
     189      ! Do some CICE consistency checks 
    189190      IF ( (ksbc == jp_flx) .OR. (ksbc == jp_purecpl) ) THEN 
    190191         IF ( calc_strair .OR. calc_Tsfc ) THEN 
     
    198199 
    199200 
    200 ! allocate sbc_ice and sbc_cice arrays 
    201       IF( sbc_ice_alloc()      /= 0 )   CALL ctl_stop( 'STOP', 'sbc_ice_cice_alloc : unable to allocate arrays' ) 
     201      ! allocate sbc_ice and sbc_cice arrays 
     202      IF( sbc_ice_alloc()      /= 0 )   CALL ctl_stop( 'STOP', 'sbc_ice_alloc : unable to allocate arrays' ) 
    202203      IF( sbc_ice_cice_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_ice_cice_alloc : unable to allocate cice arrays' ) 
    203204 
    204 ! Ensure ocean temperatures are nowhere below freezing if not a NEMO restart 
     205      ! Ensure that no temperature points are below freezing if not a NEMO restart 
    205206      IF( .NOT. ln_rstart ) THEN 
    206          tsn(:,:,:,jp_tem) = MAX (tsn(:,:,:,jp_tem),Tocnfrz) 
     207 
     208         CALL wrk_alloc( jpi,jpj,jpk, ztfrz3d )  
     209         DO jk=1,jpk 
     210            ztfrz3d(:,:,jk) = eos_fzp( tsn(:,:,jk,jp_sal), fsdept_n(:,:,jk) ) 
     211         ENDDO 
     212         tsn(:,:,:,jp_tem) = MAX( tsn(:,:,:,jp_tem), ztfrz3d ) 
    207213         tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) 
    208       ENDIF 
    209  
    210       fr_iu(:,:)=0.0 
    211       fr_iv(:,:)=0.0 
     214         CALL wrk_dealloc( jpi,jpj,jpk, ztfrz3d )  
     215 
     216#if defined key_nemocice_decomp 
     217         ! Pass initial SST from NEMO to CICE so ice is initialised correctly if 
     218         ! there is no restart file. 
     219         ! Values from a CICE restart file would overwrite this 
     220         CALL nemo2cice( tsn(:,:,1,jp_tem) , sst , 'T' , 1.)  
     221#endif 
     222 
     223      ENDIF   
     224 
     225      ! calculate surface freezing temperature and send to CICE 
     226      sstfrz(:,:) = eos_fzp(sss_m(:,:), fsdept_n(:,:,1))  
     227      CALL nemo2cice(sstfrz,Tf, 'T', 1. ) 
    212228 
    213229      CALL cice2nemo(aice,fr_i, 'T', 1. ) 
     
    220236! T point to U point 
    221237! T point to V point 
     238      fr_iu(:,:)=0.0 
     239      fr_iv(:,:)=0.0 
    222240      DO jj=1,jpjm1 
    223241         DO ji=1,jpim1 
     
    343361         CALL nemo2cice(ztmp,stray,'F', -1. ) 
    344362 
     363 
     364! Alex West: From configuration GSI8 onwards, when NEMO is used with CICE in 
     365! HadGEM3 the 'time-travelling ice' coupling approach is used, whereby  
     366! atmosphere-ice fluxes are passed as pseudo-local values, formed by dividing 
     367! gridbox mean fluxes in the UM by future ice concentration obtained through   
     368! OASIS.  This allows for a much more realistic apportionment of energy through 
     369! the ice - and conserves energy. 
     370! Therefore the fluxes are now divided by ice concentration in the coupled 
     371! formulation (jp_purecpl) as well as for jp_flx.  This NEMO branch should only 
     372! be used at UM10.2 onwards (unless an explicit GSI8 UM branch is included), at 
     373! which point the GSI8 UM changes were committed. 
     374 
    345375! Surface downward latent heat flux (CI_5) 
    346376         IF (ksbc == jp_flx) THEN 
     
    348378               ztmpn(:,:,jl)=qla_ice(:,:,1)*a_i(:,:,jl) 
    349379            ENDDO 
    350          ELSE 
    351 ! emp_ice is set in sbc_cpl_ice_flx as sublimation-snow 
    352             qla_ice(:,:,1)= - ( emp_ice(:,:)+sprecip(:,:) ) * Lsub 
    353 ! End of temporary code 
    354             DO jj=1,jpj 
    355                DO ji=1,jpi 
    356                   IF (fr_i(ji,jj).eq.0.0) THEN 
    357                      DO jl=1,ncat 
    358                         ztmpn(ji,jj,jl)=0.0 
    359                      ENDDO 
    360                      ! This will then be conserved in CICE 
    361                      ztmpn(ji,jj,1)=qla_ice(ji,jj,1) 
    362                   ELSE 
    363                      DO jl=1,ncat 
    364                         ztmpn(ji,jj,jl)=qla_ice(ji,jj,1)*a_i(ji,jj,jl)/fr_i(ji,jj) 
    365                      ENDDO 
    366                   ENDIF 
    367                ENDDO 
     380         ELSE IF (ksbc == jp_purecpl) THEN 
     381            DO jl=1,ncat 
     382               ztmpn(:,:,jl)=qla_ice(:,:,jl)*a_i(:,:,jl) 
    368383            ENDDO 
     384    ELSE 
     385           !In coupled mode - qla_ice calculated in sbc_cpl for each category 
     386           ztmpn(:,:,1:ncat)=qla_ice(:,:,1:ncat) 
    369387         ENDIF 
     388 
    370389         DO jl=1,ncat 
    371390            CALL nemo2cice(ztmpn(:,:,jl),flatn_f(:,:,jl,:),'T', 1. ) 
     
    373392! GBM conductive flux through ice (CI_6) 
    374393!  Convert to GBM 
    375             IF (ksbc == jp_flx) THEN 
     394            IF (ksbc == jp_flx .OR. ksbc == jp_purecpl) THEN 
    376395               ztmp(:,:) = botmelt(:,:,jl)*a_i(:,:,jl) 
    377396            ELSE 
     
    382401! GBM surface heat flux (CI_7) 
    383402!  Convert to GBM 
    384             IF (ksbc == jp_flx) THEN 
     403            IF (ksbc == jp_flx .OR. ksbc == jp_purecpl) THEN 
    385404               ztmp(:,:) = (topmelt(:,:,jl)+botmelt(:,:,jl))*a_i(:,:,jl)  
    386405            ELSE 
     
    442461      CALL nemo2cice(ztmp,frain,'T', 1. )  
    443462 
     463! Recalculate freezing temperature and send to CICE  
     464      sstfrz(:,:)=eos_fzp(sss_m(:,:), fsdept_n(:,:,1))  
     465      CALL nemo2cice(sstfrz,Tf,'T', 1. ) 
     466 
    444467! Freezing/melting potential 
    445468! Calculated over NEMO leapfrog timestep (hence 2*dt) 
    446       nfrzmlt(:,:)=rau0*rcp*fse3t_m(:,:)*(Tocnfrz-sst_m(:,:))/(2.0*dt) 
    447  
    448       ztmp(:,:) = nfrzmlt(:,:) 
    449       CALL nemo2cice(ztmp,frzmlt,'T', 1. ) 
     469      nfrzmlt(:,:)=rau0*rcp*fse3t_m(:,:)*(sstfrz(:,:)-sst_m(:,:))/(2.0*dt)  
     470      CALL nemo2cice(nfrzmlt,frzmlt,'T', 1. ) 
    450471 
    451472! SST  and SSS 
     
    453474      CALL nemo2cice(sst_m,sst,'T', 1. ) 
    454475      CALL nemo2cice(sss_m,sss,'T', 1. ) 
     476 
     477      IF( ksbc == jp_purecpl ) THEN 
     478! Sea ice surface skin temperature 
     479         DO jl=1,ncat 
     480           CALL nemo2cice(tsfc_ice(:,:,jl), trcrn(:,:,nt_tsfc,jl,:),'T',1.) 
     481         ENDDO  
     482      ENDIF 
    455483 
    456484! x comp and y comp of surface ocean current 
     
    730758         CALL cice2nemo(vicen(:,:,jl,:),ht_i(:,:,jl),'T', 1. ) 
    731759      ENDDO 
     760 
     761#if ! defined key_cice4 
     762! Meltpond fraction and depth 
     763      DO jl = 1,ncat 
     764         CALL cice2nemo(apeffn(:,:,jl,:),a_p(:,:,jl),'T', 1. ) 
     765         CALL cice2nemo(trcrn(:,:,nt_hpnd,jl,:),ht_p(:,:,jl),'T', 1. ) 
     766      ENDDO 
     767#endif 
     768 
     769 
     770! If using multilayers thermodynamics in CICE then get top layer temperature 
     771! and effective conductivity        
     772!! When using NEMO with CICE, this change requires use of  
     773!! one of the following two CICE branches: 
     774!! - at CICE5.0,   hadax/r1015_GSI8_with_GSI7 
     775!! - at CICE5.1.2, hadax/vn5.1.2_GSI8 
     776      IF (heat_capacity) THEN 
     777         DO jl = 1,ncat 
     778            CALL cice2nemo(Tn_top(:,:,jl,:),tn_ice(:,:,jl),'T', 1. ) 
     779            CALL cice2nemo(keffn_top(:,:,jl,:),kn_ice(:,:,jl),'T', 1. ) 
     780         ENDDO 
     781! Convert surface temperature to Kelvin 
     782         tn_ice(:,:,:)=tn_ice(:,:,:)+rt0 
     783      ELSE 
     784         tn_ice(:,:,:) = 0.0 
     785         kn_ice(:,:,:) = 0.0 
     786      ENDIF        
     787 
    732788      ! 
    733789      IF( nn_timing == 1 )  CALL timing_stop('cice_sbc_hadgam') 
  • branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90

    r6487 r6488  
    2727   USE fldread         ! read input field at current time step 
    2828 
    29  
    30  
    3129   IMPLICIT NONE 
    3230   PRIVATE 
     
    8482    INTEGER                      ::   ji, jj, jk, ijkmin, inum, ierror 
    8583    INTEGER                      ::   ikt, ikb   ! top and bottom level of the isf boundary layer 
     84    REAL(wp)                     ::   zgreenland_fwfisf_sum, zantarctica_fwfisf_sum 
    8685    REAL(wp)                     ::   rmin 
    8786    REAL(wp)                     ::   zhk 
     
    250249            CALL fld_read ( kt, nn_fsbc, sf_rnfisf   ) 
    251250            fwfisf(:,:) = - sf_rnfisf(1)%fnow(:,:,1)         ! fresh water flux from the isf (fwfisf <0 mean melting)  
     251 
     252            IF( lk_oasis) THEN 
     253            ! ln_coupled_iceshelf_fluxes uninitialised unless lk_oasis=true 
     254            IF( ln_coupled_iceshelf_fluxes ) THEN 
     255 
     256              ! Adjust total iceshelf melt rates so that sum of iceberg calving and iceshelf melting in the northern 
     257              ! and southern hemispheres equals rate of increase of mass of greenland and antarctic ice sheets 
     258              ! to preserve total freshwater conservation in coupled models without an active ice sheet model. 
     259 
     260               zgreenland_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * greenland_icesheet_mask(:,:) ) 
     261               IF( lk_mpp ) CALL mpp_sum( zgreenland_fwfisf_sum ) 
     262               ! use ABS function because we need to preserve the sign of fwfisf 
     263               WHERE( greenland_icesheet_mask(:,:) == 1.0 )                                                                  & 
     264              &    fwfisf(:,:) = fwfisf(:,:)  * ABS( greenland_icesheet_mass_rate_of_change * (1.0-rn_greenland_calving_fraction) & 
     265              &                           / ( zgreenland_fwfisf_sum + 1.0e-10_wp ) ) 
     266 
     267               ! check 
     268               IF(lwp) WRITE(numout, *) 'Greenland iceshelf melting climatology (kg/s) : ',zgreenland_fwfisf_sum 
     269               zgreenland_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * greenland_icesheet_mask(:,:) ) 
     270               IF( lk_mpp ) CALL mpp_sum( zgreenland_fwfisf_sum ) 
     271               IF(lwp) WRITE(numout, *) 'Greenland iceshelf melting adjusted value (kg/s) : ',zgreenland_fwfisf_sum 
     272 
     273               zantarctica_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * antarctica_icesheet_mask(:,:) ) 
     274               IF( lk_mpp ) CALL mpp_sum( zantarctica_fwfisf_sum ) 
     275               ! use ABS function because we need to preserve the sign of fwfisf 
     276               WHERE( antarctica_icesheet_mask(:,:) == 1.0 ) & 
     277              &    fwfisf(:,:) = fwfisf(:,:)  * ABS( antarctica_icesheet_mass_rate_of_change * (1.0-rn_antarctica_calving_fraction) & 
     278              &                           / ( zantarctica_fwfisf_sum + 1.0e-10_wp ) ) 
     279       
     280               ! check 
     281               IF(lwp) WRITE(numout, *) 'Antarctica iceshelf melting climatology (kg/s) : ',zantarctica_fwfisf_sum 
     282               zantarctica_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * antarctica_icesheet_mask(:,:) ) 
     283               IF( lk_mpp ) CALL mpp_sum( zantarctica_fwfisf_sum ) 
     284               IF(lwp) WRITE(numout, *) 'Antarctica iceshelf melting adjusted value (kg/s) : ',zantarctica_fwfisf_sum 
     285 
     286            ENDIF 
     287            ENDIF 
     288 
    252289            qisf(:,:)   = fwfisf(:,:) * lfusisf              ! heat        flux 
    253290            stbl(:,:)   = soce 
     
    258295            !CALL fld_read ( kt, nn_fsbc, sf_qisf   ) 
    259296            fwfisf(:,:) = sf_fwfisf(1)%fnow(:,:,1)            ! fwf 
     297 
     298            IF( lk_oasis) THEN 
     299            ! ln_coupled_iceshelf_fluxes uninitialised unless lk_oasis=true 
     300            IF( ln_coupled_iceshelf_fluxes ) THEN 
     301 
     302              ! Adjust total iceshelf melt rates so that sum of iceberg calving and iceshelf melting in the northern 
     303              ! and southern hemispheres equals rate of increase of mass of greenland and antarctic ice sheets 
     304              ! to preserve total freshwater conservation in coupled models without an active ice sheet model. 
     305 
     306               zgreenland_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * greenland_icesheet_mask(:,:) ) 
     307               IF( lk_mpp ) CALL mpp_sum( zgreenland_fwfisf_sum ) 
     308               ! use ABS function because we need to preserve the sign of fwfisf 
     309               WHERE( greenland_icesheet_mask(:,:) == 1.0 )                                                                  & 
     310              &    fwfisf(:,:) = fwfisf(:,:)  * ABS( greenland_icesheet_mass_rate_of_change * (1.0-rn_greenland_calving_fraction) & 
     311              &                           / ( zgreenland_fwfisf_sum + 1.0e-10_wp ) ) 
     312 
     313               ! check 
     314               IF(lwp) WRITE(numout, *) 'Greenland iceshelf melting climatology (kg/s) : ',zgreenland_fwfisf_sum 
     315               zgreenland_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * greenland_icesheet_mask(:,:) ) 
     316               IF( lk_mpp ) CALL mpp_sum( zgreenland_fwfisf_sum ) 
     317               IF(lwp) WRITE(numout, *) 'Greenland iceshelf melting adjusted value (kg/s) : ',zgreenland_fwfisf_sum 
     318 
     319               zantarctica_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * antarctica_icesheet_mask(:,:) ) 
     320               IF( lk_mpp ) CALL mpp_sum( zantarctica_fwfisf_sum ) 
     321               ! use ABS function because we need to preserve the sign of fwfisf 
     322               WHERE( antarctica_icesheet_mask(:,:) == 1.0 ) & 
     323              &    fwfisf(:,:) = fwfisf(:,:)  * ABS( antarctica_icesheet_mass_rate_of_change * (1.0-rn_antarctica_calving_fraction) & 
     324              &                           / ( zantarctica_fwfisf_sum + 1.0e-10_wp ) ) 
     325       
     326               ! check 
     327               IF(lwp) WRITE(numout, *) 'Antarctica iceshelf melting climatology (kg/s) : ',zantarctica_fwfisf_sum 
     328               zantarctica_fwfisf_sum = SUM( fwfisf(:,:) * e1t(:,:) * e2t(:,:) * antarctica_icesheet_mask(:,:) ) 
     329               IF( lk_mpp ) CALL mpp_sum( zantarctica_fwfisf_sum ) 
     330               IF(lwp) WRITE(numout, *) 'Antarctica iceshelf melting adjusted value (kg/s) : ',zantarctica_fwfisf_sum 
     331 
     332            ENDIF 
     333            ENDIF 
     334 
    260335            qisf(:,:)   = fwfisf(:,:) * lfusisf              ! heat        flux 
    261336            !qisf(:,:)   = sf_qisf(1)%fnow(:,:,1)              ! heat flux 
Note: See TracChangeset for help on using the changeset viewer.