Calve icebergs from ice-shelf flux

full/coupler_main.F90

@@ -278,6 +278,14 @@
 !!       the slow sea-ice processes are on the same PEs as the fast sea-ice.</td>
 !!   </tr>
 !!   <tr>
+!!     <td>calve_ice_shelf_bergs</td>
+!!     <td>logical</td>
+!!     <td>.FALSE.</td>
+!!      <td> If true, the ice sheet flux through a fixed ice-shelf front is
+!!           converted to icebergs, rather than initializing icebergs from frozen
+!!           freshwater discharge.</td>
+!!   </tr>
+!!   <tr>
 !!     <td>restart_interval</td>
 !!     <td>integer, dimension(6)</td>
 !!     <td>(/0,0,0,0,0,0/)</td>
@@ -447,9 +455,14 @@ program coupler_main
       ! slow ice are on different PEs.  call fms_mpp_set_current_pelist(Ice%pelist)
       ! is called if(.not.Ice%shared_slow_fast_PEs)
       call coupler_exchange_slow_to_fast_ice(Ice, coupler_clocks)
-      !> This call occurs all ice PEs.
-      if (concurrent_ice) call coupler_exchange_fast_to_slow_ice(Ice, coupler_clocks)
-      !> call fms_mpp_set_current_pelist(Ice%pelist) is called if(.not.Ice%shared_slow_fast_PEs)
+
+      if (concurrent_ice) then
+        !> This call occurs all ice PEs.
+        call coupler_exchange_fast_to_slow_ice(Ice, coupler_clocks)
+        if (Ice%slow_ice_pe .and. calve_ice_shelf_bergs) &
+          call coupler_unpack_ocean_ice_boundary_calved_ice_shelf_bergs(Ice, Ocean_ice_boundary, coupler_clocks)
+      endif
+
       if (Ice%fast_ice_pe) call coupler_set_ice_surface_fields(Ice, coupler_clocks)
     endif
 
@@ -605,8 +618,11 @@ program coupler_main
     if (do_ice .and. Ice%pe) then
       if (Ice%fast_ice_PE) call coupler_unpack_land_ice_boundary(Ice, Land_ice_boundary, coupler_clocks)
       !> This could be a point where the model is serialized; This calls on all ice PEs
-      if (.not.concurrent_ice) call coupler_exchange_fast_to_slow_ice(Ice, coupler_clocks, &
-                               set_ice_current_pelist=.True.)
+      if (.not.concurrent_ice) then
+        call coupler_exchange_fast_to_slow_ice(Ice, coupler_clocks, set_ice_current_pelist=.True.)
+        if (Ice%slow_ice_pe .and. calve_ice_shelf_bergs) &
+          call coupler_unpack_ocean_ice_boundary_calved_ice_shelf_bergs(Ice, Ocean_ice_boundary, coupler_clocks)
+      endif
       !> slow-ice model
       !! This call occurs on whichever PEs handle the slow ice processess.
       if (Ice%slow_ice_PE .and. .not.combined_ice_and_ocean) &

full/full_coupler_mod.F90

@@ -60,6 +60,7 @@ module full_coupler_mod
   use ice_model_mod,           only: ice_data_type_chksum, ocn_ice_bnd_type_chksum
   use ice_model_mod,           only: atm_ice_bnd_type_chksum, lnd_ice_bnd_type_chksum
   use ice_model_mod,           only: unpack_ocean_ice_boundary, exchange_slow_to_fast_ice
+  use ice_model_mod,           only: unpack_ocean_ice_boundary_calved_shelf_bergs
   use ice_model_mod,           only: ice_model_fast_cleanup, unpack_land_ice_boundary
   use ice_model_mod,           only: exchange_fast_to_slow_ice, update_ice_model_slow
 
@@ -128,6 +129,7 @@ module full_coupler_mod
 
   public :: coupler_update_land_model_slow, coupler_flux_land_to_ice
   public :: coupler_unpack_land_ice_boundary, coupler_flux_ice_to_ocean
+  public :: coupler_unpack_ocean_ice_boundary_calved_ice_shelf_bergs
   public :: coupler_update_ice_model_slow_and_stocks, coupler_update_ocean_model
 
   public :: coupler_clock_type, coupler_components_type, coupler_chksum_type
@@ -219,6 +221,9 @@ module full_coupler_mod
   logical, public :: do_debug=.FALSE.!< If .TRUE. print additional debugging messages.
   integer, public :: check_stocks = 0 !< -1: never 0: at end of run only n>0: every n coupled steps
   logical, public :: use_hyper_thread = .false.
+  logical, public :: calve_ice_shelf_bergs = .false. !< If true, the ice sheet flux through a fixed ice-shelf front is
+                                                     !! converted to icebergs, rather than initializing icebergs from
+                                                     !! frozen freshwater discharge
 
   namelist /coupler_nml/ current_date, calendar, force_date_from_namelist,         &
                          months, days, hours, minutes, seconds, dt_cpld, dt_atmos, &
@@ -228,7 +233,8 @@ module full_coupler_mod
                          concurrent, do_concurrent_radiation, use_lag_fluxes,      &
                          check_stocks, restart_interval, do_debug, do_chksum,      &
                          use_hyper_thread, concurrent_ice, slow_ice_with_ocean,    &
-                         do_endpoint_chksum, combined_ice_and_ocean
+                         do_endpoint_chksum, combined_ice_and_ocean,               &
+                         calve_ice_shelf_bergs
 
   !> coupler_clock_type derived type consist of all clock ids that will be set and used
   !! in full coupler_main.
@@ -991,7 +997,8 @@ subroutine coupler_init(Atm, Ocean, Land, Ice, Ocean_state, Atmos_land_boundary,
 
       call fms_mpp_clock_begin(coupler_clocks%ocean_model_init)
       call ocean_model_init( Ocean, Ocean_state, Time_init, Time, &
-                             gas_fields_ocn=gas_fields_ocn  )
+                             gas_fields_ocn=gas_fields_ocn, &
+                             calve_ice_shelf_bergs=calve_ice_shelf_bergs)
       call fms_mpp_clock_end(coupler_clocks%ocean_model_init)
 
       if (concurrent) then
@@ -2305,6 +2312,23 @@ subroutine coupler_unpack_land_ice_boundary(Ice, Land_ice_boundary, coupler_cloc
 
   end subroutine coupler_unpack_land_ice_boundary
 
+  !> This subroutine calls unpack_ocean_ice_boundary_calved_shelf_bergs
+  subroutine coupler_unpack_ocean_ice_boundary_calved_ice_shelf_bergs(Ice, Ocean_ice_boundary, coupler_clocks)
+
+    implicit none
+    type(ice_data_type),          intent(inout) :: Ice                 !< Ice
+    type(ocean_ice_boundary_type), intent(inout) :: Ocean_ice_boundary !< Ocean_ice_boundary
+    type(coupler_clock_type), intent(inout) :: coupler_clocks          !< coupler_clocks
+
+    call fms_mpp_set_current_pelist(Ice%slow_pelist)
+    call fms_mpp_clock_begin(coupler_clocks%update_ice_model_slow_slow)
+
+    call unpack_ocean_ice_boundary_calved_shelf_bergs(Ice, Ocean_ice_boundary)
+
+    call fms_mpp_clock_end(coupler_clocks%update_ice_model_slow_slow)
+
+  end subroutine coupler_unpack_ocean_ice_boundary_calved_ice_shelf_bergs
+
   !> This subroutine calls update_ice_model_slow and flux_ice_to_ocean_stocks
   subroutine coupler_update_ice_model_slow_and_stocks(Ice, coupler_clocks)
 

full/ice_ocean_flux_exchange.F90

@@ -90,13 +90,17 @@ subroutine ice_ocean_flux_exchange_init(Time, Ice, Ocean, Ocean_state, ice_ocean
     !frazil and sea_level are optional, if not present they should be nullified
     allocate( ocean_ice_boundary%frazil(is:ie,js:je) )
     allocate( ocean_ice_boundary%sea_level(is:ie,js:je) )
+    allocate( ocean_ice_boundary%calving(is:ie,js:je) )
+    allocate( ocean_ice_boundary%calving_hflx(is:ie,js:je) )
     ! initialize boundary fields for override experiments (mjh)
     ocean_ice_boundary%u=0.0
     ocean_ice_boundary%v=0.0
     ocean_ice_boundary%t=273.0
     ocean_ice_boundary%s=0.0
     ocean_ice_boundary%frazil=0.0
     ocean_ice_boundary%sea_level=0.0
+    ocean_ice_boundary%calving=0.0
+    ocean_ice_boundary%calving_hflx=0.0
 
     ! allocate fields for extra tracers in ocean_ice_boundary
     if (.not.fms_coupler_type_initialized(ocean_ice_boundary%fields)) &
@@ -223,7 +227,7 @@ end subroutine ice_ocean_flux_exchange_init
   !!        fprec = mass of frozen precipitation since last
   !!                  time step (Kg/m2)
   !!       runoff = mass of runoff since last time step (Kg/m2)
-  !!       runoff = mass of calving since last time step (Kg/m2)
+  !!       calving = mass of calving since last time step (Kg/m2)
   !!       p_surf = surface pressure (Pa)
   !! </pre>
   subroutine flux_ice_to_ocean ( Ice, Ocean, Ice_Ocean_Boundary )
@@ -372,7 +376,9 @@ end subroutine flux_ice_to_ocean_finish
   !!        u_surf = zonal ocean current/ice motion (m/s)
   !!        v_surf = meridional ocean current/ice motion (m/s)
   !!        v_surf = meridional ocean current/ice motion (m/s)
-  !!       sea_lev = sea level used to drive ice accelerations (m)
+  !!        sea_lev = sea level used to drive ice accelerations (m)
+  !!        calving = ice-sheet calving flux to ocean (kg/m2/s)
+  !!        calving_hflx = heat flux associated with ice-sheet calving (W/m2)
   !! </pre>
   !!
   !! \throw FATAL, "Ocean_Ice_Boundary%xtype must be DIRECT or REDIST."
@@ -408,6 +414,23 @@ subroutine flux_ocean_to_ice ( Ocean, Ice, Ocean_Ice_Boundary )
           endif
        endif
 
+       if( ASSOCIATED(Ocean_Ice_Boundary%calving) ) then
+          if(do_area_weighted_flux) then
+             Ocean_Ice_Boundary%calving = Ocean%calving * Ocean%area
+             call divide_by_area(data=Ocean_Ice_Boundary%calving, area=Ice%area)
+           else
+             Ocean_Ice_Boundary%calving = Ocean%calving
+          endif
+       endif
+       if( ASSOCIATED(Ocean_Ice_Boundary%calving_hflx) ) then
+          if(do_area_weighted_flux) then
+             Ocean_Ice_Boundary%calving_hflx = Ocean%calving_hflx * Ocean%area
+             call divide_by_area(data=Ocean_Ice_Boundary%calving_hflx, area=Ice%area)
+          else
+             Ocean_Ice_Boundary%calving_hflx = Ocean%calving_hflx
+          endif
+       endif
+
        ! Extra fluxes
        call fms_coupler_type_copy_data(Ocean%fields, Ocean_Ice_Boundary%fields)
 
@@ -439,6 +462,35 @@ subroutine flux_ocean_to_ice ( Ocean, Ice, Ocean_Ice_Boundary )
           else
              call fms_mpp_domains_redistribute(Ocean%Domain,Ocean%frazil, Ice%slow_Domain_NH, Ocean_Ice_Boundary%frazil)
           endif
+        endif
+
+       if( ASSOCIATED(Ocean_Ice_Boundary%calving) ) then
+          if(do_area_weighted_flux) then
+             if (Ocean%is_ocean_pe) then
+               allocate(tmp(size(Ocean%area,1), size(Ocean%area,2)))
+               tmp(:,:) = Ocean%calving(:,:) * Ocean%area(:,:)
+             endif
+             call fms_mpp_domains_redistribute( Ocean%Domain, tmp, Ice%slow_Domain_NH, Ocean_Ice_Boundary%calving)
+             if (Ice%slow_ice_pe) &
+               call divide_by_area(data=Ocean_Ice_Boundary%calving, area=Ice%area)
+             if (Ocean%is_ocean_pe) deallocate(tmp)
+          else
+             call fms_mpp_domains_redistribute(Ocean%Domain, Ocean%calving, Ice%slow_Domain_NH, Ocean_Ice_Boundary%calving)
+          endif
+       endif
+       if( ASSOCIATED(Ocean_Ice_Boundary%calving_hflx) ) then
+          if(do_area_weighted_flux) then
+             if (Ocean%is_ocean_pe) then
+               allocate(tmp(size(Ocean%area,1), size(Ocean%area,2)))
+               tmp(:,:) = Ocean%calving_hflx(:,:) * Ocean%area(:,:)
+             endif
+             call fms_mpp_domains_redistribute( Ocean%Domain, tmp, Ice%slow_Domain_NH, Ocean_Ice_Boundary%calving_hflx)
+             if (Ice%slow_ice_pe) &
+               call divide_by_area(data=Ocean_Ice_Boundary%calving_hflx, area=Ice%area)
+             if (Ocean%is_ocean_pe) deallocate(tmp)
+          else
+             call fms_mpp_domains_redistribute(Ocean%Domain, Ocean%calving_hflx, Ice%slow_Domain_NH, Ocean_Ice_Boundary%calving_hflx)
+          endif
        endif
 
        ! Extra fluxes
@@ -470,6 +522,9 @@ subroutine flux_ocean_to_ice_finish( Time, Ice, Ocean_Ice_Boundary )
     call fms_data_override('ICE', 's',         Ocean_Ice_Boundary%s,         Time)
     call fms_data_override('ICE', 'frazil',    Ocean_Ice_Boundary%frazil,    Time)
     call fms_data_override('ICE', 'sea_level', Ocean_Ice_Boundary%sea_level, Time)
+    !Ice-shelf calving
+    call fms_data_override('ICE', 'IS_calving',      Ocean_Ice_Boundary%calving,      Time)
+    call fms_data_override('ICE', 'IS_calving_hflx', Ocean_Ice_Boundary%calving_hflx, Time)
     call fms_coupler_type_data_override('ICE', Ocean_Ice_Boundary%fields, Time)
 
     !  Perform diagnostic output for the ocean_ice_boundary fields