diff --git a/src/core_seaice/Registry.xml b/src/core_seaice/Registry.xml
index 2810cb33c9..a60caafb99 100644
--- a/src/core_seaice/Registry.xml
+++ b/src/core_seaice/Registry.xml
@@ -1685,6 +1685,7 @@
+
diff --git a/src/core_seaice/column/Makefile b/src/core_seaice/column/Makefile
index 822ecc1305..ca9b303e2f 100644
--- a/src/core_seaice/column/Makefile
+++ b/src/core_seaice/column/Makefile
@@ -28,7 +28,8 @@ OBJS = ice_colpkg.o \
ice_brine.o \
ice_algae.o \
ice_zbgc.o \
- ice_zbgc_shared.o
+ ice_zbgc_shared.o \
+ ice_bergs_mechred.o
all: $(OBJS)
@@ -40,7 +41,7 @@ ice_constants_colpkg.F90:
cp constants/cice/ice_constants_colpkg.F90 .
endif
-ice_colpkg.o: ice_kinds_mod.o ice_constants_colpkg.o ice_warnings.o ice_colpkg_shared.o ice_therm_shared.o ice_orbital.o ice_atmo.o ice_age.o ice_firstyear.o ice_flux_colpkg.o ice_meltpond_cesm.o ice_meltpond_lvl.o ice_meltpond_topo.o ice_therm_vertical.o ice_itd.o ice_therm_itd.o ice_shortwave.o ice_mechred.o ice_colpkg_tracers.o ice_atmo.o ice_mushy_physics.o ice_zbgc.o ice_zbgc_shared.o ice_aerosol.o ice_algae.o ice_brine.o ice_zsalinity.o
+ice_colpkg.o: ice_kinds_mod.o ice_constants_colpkg.o ice_warnings.o ice_colpkg_shared.o ice_therm_shared.o ice_orbital.o ice_atmo.o ice_age.o ice_firstyear.o ice_flux_colpkg.o ice_meltpond_cesm.o ice_meltpond_lvl.o ice_meltpond_topo.o ice_therm_vertical.o ice_itd.o ice_therm_itd.o ice_shortwave.o ice_mechred.o ice_colpkg_tracers.o ice_atmo.o ice_mushy_physics.o ice_zbgc.o ice_zbgc_shared.o ice_aerosol.o ice_algae.o ice_brine.o ice_zsalinity.o ice_bergs_mechred.o
ice_kinds_mod.o:
@@ -96,7 +97,9 @@ ice_brine.o: ice_kinds_mod.o ice_constants_colpkg.o ice_colpkg_shared.o ice_colp
ice_zbgc.o: ice_kinds_mod.o ice_constants_colpkg.o ice_colpkg_shared.o ice_colpkg_tracers.o ice_zbgc_shared.o ice_therm_shared.o ice_itd.o
-ice_zbgc_shared.o: ice_kinds_mod.o ice_constants_colpkg.o ice_colpkg_shared.o
+ice_zbgc_shared.o: ice_kinds_mod.o ice_constants_colpkg.o ice_colpkg_shared.o
+
+ice_bergs_mechred.o: ice_kinds_mod.o ice_constants_colpkg.o ice_colpkg_shared.o ice_mechred.o ice_itd.o ice_colpkg_tracers.o
.F90.o:
diff --git a/src/core_seaice/column/ice_bergs_mechred.F90 b/src/core_seaice/column/ice_bergs_mechred.F90
new file mode 100644
index 0000000000..ece890e7f3
--- /dev/null
+++ b/src/core_seaice/column/ice_bergs_mechred.F90
@@ -0,0 +1,558 @@
+! SVN:$Id: $
+!=======================================================================
+
+! Ridging of sea ice by advancing ice shelves or icebergs
+! based on ice_mechred.F90
+!
+! authors: Elizabeth C. Hunke, LANL
+
+ module ice_bergs_mechred
+
+ use ice_kinds_mod
+ use ice_constants_colpkg, only: c0, c1, c2, c10, p5, &
+ puny, Lfresh, rhoi, rhos
+ use ice_itd, only: column_sum, &
+ column_conservation_check
+ use ice_mechred, only: asum_ridging, ridge_itd, ridge_shift, ridge_prep
+ use ice_colpkg_tracers, only: nt_qice, nt_qsno, nt_fbri, nt_sice
+ use ice_warnings, only: add_warning
+
+ implicit none
+ save
+
+ private
+ public :: ridge_ice_by_bergs
+
+ logical (kind=log_kind), parameter :: &
+ l_conservation_check = .false. ! if true, check conservation
+! l_conservation_check = .true. ! if true, check conservation
+
+!=======================================================================
+
+ contains
+
+!=======================================================================
+
+! Compute changes in the sea ice thickness distribution due to icebergs.
+! Based on ice_mechred.F90
+!
+! NOTE: This subroutine operates over a single block.
+!
+! !REVISION HISTORY:
+!
+! author: William H. Lipscomb, LANL
+! 2008: Elizabeth Hunke modified for icebergs.
+
+ subroutine ridge_ice_by_bergs (dt, dice, &
+ ncat, n_aero, &
+ nilyr, nslyr, &
+ ntrcr, hin_max, &
+ aicen, btrcrn, &
+ vicen, vsnon, &
+ aice0, bice, &
+ trcr_depend, trcr_base , &
+ n_trcr_strata, &
+ nt_strata, l_stop, &
+ stop_label, &
+ krdg_partic, krdg_redist, &
+ mu_rdg, tr_brine, &
+ dardg1dt, dardg2dt, &
+ dvirdgdt, opening, &
+ fpond, &
+ fresh, fhocn, &
+ faero_ocn, &
+ aparticn, krdgn, &
+ aredistn, vredistn, &
+ dardg1ndt, dardg2ndt, &
+ dvirdgndt, &
+ araftn, vraftn)
+
+!
+! !INPUT/OUTPUT PARAMETERS:
+!
+ real (kind=dbl_kind), intent(in) :: &
+ dt , & ! time step
+ dice , & ! fractional area to be ridged (displaced)
+ bice , & ! fractional berg (and shelf) area
+ mu_rdg ! gives e-folding scale of ridged ice (m^.5)
+
+ integer (kind=int_kind), intent(in) :: &
+ ncat , & ! number of thickness categories
+ n_aero , & ! number of aerosol tracers
+ nilyr , & ! number of ice layers
+ nslyr , & ! number of snow layers
+ ! nblyr , & ! number of bio layers
+ ntrcr ! number of tracers in use
+
+ real (kind=dbl_kind), dimension(0:ncat), intent(inout) :: &
+ hin_max ! category limits (m)
+
+ real (kind=dbl_kind), dimension (ncat), intent(inout) :: &
+ aicen , & ! concentration of ice
+ vicen , & ! volume per unit area of ice (m)
+ vsnon ! volume per unit area of snow (m)
+
+ real (kind=dbl_kind), dimension (:,:), intent(inout) :: &
+ btrcrn ! temporary array for tracers
+
+ real (kind=dbl_kind), intent(inout) :: &
+ aice0 ! concentration of open water
+
+ integer (kind=int_kind), dimension (:), intent(in) :: &
+ trcr_depend, & ! = 0 for aicen tracers, 1 for vicen, 2 for vsnon
+ n_trcr_strata ! number of underlying tracer layers
+
+ real (kind=dbl_kind), dimension (:,:), intent(in) :: &
+ trcr_base ! = 0 or 1 depending on tracer dependency
+ ! argument 2: (1) aice, (2) vice, (3) vsno
+
+ integer (kind=int_kind), dimension (:,:), intent(in) :: &
+ nt_strata ! indices of underlying tracer layers
+
+ logical (kind=log_kind), intent(out) :: &
+ l_stop ! if true, abort on return
+
+ character (char_len), intent(out) :: &
+ stop_label ! diagnostic information for abort
+
+ integer (kind=int_kind), intent(in) :: &
+ krdg_partic , & ! selects participation function
+ krdg_redist ! selects redistribution function
+
+ logical (kind=log_kind), intent(in) :: &
+ tr_brine ! if .true., brine height differs from ice thickness
+
+ ! optional history fields
+ real (kind=dbl_kind), intent(inout), optional :: &
+ dardg1dt , & ! rate of fractional area loss by ridging ice (1/s)
+ dardg2dt , & ! rate of fractional area gain by new ridges (1/s)
+ dvirdgdt , & ! rate of ice volume ridged (m/s)
+ opening , & ! rate of opening due to divergence/shear (1/s)
+ fpond , & ! fresh water flux to ponds (kg/m^2/s)
+ fresh , & ! fresh water flux to ocean (kg/m^2/s)
+ fhocn ! net heat flux to ocean (W/m^2)
+
+ real (kind=dbl_kind), dimension(:), intent(inout), optional :: &
+ dardg1ndt , & ! rate of fractional area loss by ridging ice (1/s)
+ dardg2ndt , & ! rate of fractional area gain by new ridges (1/s)
+ dvirdgndt , & ! rate of ice volume ridged (m/s)
+ aparticn , & ! participation function
+ krdgn , & ! mean ridge thickness/thickness of ridging ice
+ araftn , & ! rafting ice area
+ vraftn , & ! rafting ice volume
+ aredistn , & ! redistribution function: fraction of new ridge area
+ vredistn ! redistribution function: fraction of new ridge volume
+
+ real (kind=dbl_kind), dimension(:), intent(inout), optional :: &
+ faero_ocn ! aerosol flux to ocean (kg/m^2/s)
+
+ ! local variables
+
+ real (kind=dbl_kind), dimension (ncat) :: &
+ eicen , & ! energy of melting for each ice layer (J/m^2)
+ esnon , & ! energy of melting for each snow layer (J/m^2)
+ vbrin , & ! ice volume with defined by brine height (m)
+ sicen ! Bulk salt in h ice (ppt*m)
+
+ ! variables for ridging routines
+ real (kind=dbl_kind) :: &
+ asum , & ! sum of ice and open water area
+ aksum , & ! ratio of area removed to area ridged
+ msnow_mlt , & ! mass of snow added to ocean (kg m-2)
+ esnow_mlt , & ! energy needed to melt snow in ocean (J m-2)
+ mpond , & ! mass of pond added to ocean (kg m-2)
+ closing_net, & ! net rate at which area is removed (1/s)
+ ! (ridging ice area - area of new ridges) / dt
+ divu_adv , & ! divu as implied by transport scheme (1/s)
+ opning , & ! rate of opening due to divergence/shear
+ ! opning is a local variable;
+ ! opening is the history diagnostic variable
+ ardg1 , & ! fractional area loss by ridging ice
+ ardg2 , & ! fractional area gain by new ridges
+ virdg , & ! ice volume ridged
+ aopen ! area opening due to divergence/shear
+
+ real (kind=dbl_kind), dimension (n_aero) :: &
+ maero ! aerosol mass added to ocean (kg m-2)
+
+ real (kind=dbl_kind), dimension (0:ncat) :: &
+ apartic ! participation function; fraction of ridging
+ ! and closing associated w/ category n
+
+ real (kind=dbl_kind), dimension (ncat) :: &
+ hrmin , & ! minimum ridge thickness
+ hrmax , & ! maximum ridge thickness (krdg_redist = 0)
+ hrexp , & ! ridge e-folding thickness (krdg_redist = 1)
+ krdg , & ! mean ridge thickness/thickness of ridging ice
+ ardg1n , & ! area of ice ridged
+ ardg2n , & ! area of new ridges
+ virdgn , & ! ridging ice volume
+ mraftn ! rafting ice mask
+
+ real (kind=dbl_kind) :: &
+ vice_init, vice_final, & ! ice volume summed over categories
+ vsno_init, vsno_final, & ! snow volume summed over categories
+ eice_init, eice_final, & ! ice energy summed over layers
+ vbri_init, vbri_final, & ! ice volume in fbri*vicen summed over categories
+ sice_init ,sice_final, & ! ice bulk salinity summed over categories
+ esno_init, esno_final ! snow energy summed over layers
+
+ integer (kind=int_kind), parameter :: &
+ nitermax = 20 ! max number of ridging iterations
+
+ integer (kind=int_kind) :: &
+ k , & ! layer index
+ n , & ! thickness category index
+ niter , & ! iteration counter
+ it ! tracer index
+
+ real (kind=dbl_kind) :: &
+ dti ! 1 / dt
+
+ logical (kind=log_kind) :: &
+ iterate_ridging, & ! if true, repeat the ridging
+ asum_error ! flag for asum .ne. 1
+
+ character (len=char_len) :: &
+ fieldid ! field identifier
+
+ character(len=char_len_long) :: &
+ warning ! warning message
+
+ !-----------------------------------------------------------------
+ ! Initialize
+ !-----------------------------------------------------------------
+
+ l_stop = .false.
+
+ msnow_mlt = c0
+ esnow_mlt = c0
+ maero(:) = c0
+ ardg1 = c0
+ ardg2 = c0
+ virdg = c0
+ ardg1n(:) = c0
+ ardg2n(:) = c0
+ virdgn(:) = c0
+ mpond = c0
+ aopen = c0
+
+ divu_adv = c0
+ closing_net = c0
+ opning = c0
+ asum = c0
+
+ !-----------------------------------------------------------------
+ ! Compute initial values of conserved quantities.
+ !-----------------------------------------------------------------
+
+ if (l_conservation_check) then
+
+ do n = 1, ncat
+ eicen(n) = c0
+ esnon(n) = c0
+ sicen(n) = c0
+ vbrin(n) = c0
+
+ do k = 1, nilyr
+ eicen(n) = eicen(n) + btrcrn(nt_qice+k-1,n) &
+ * vicen(n)/real(nilyr,kind=dbl_kind)
+ enddo
+
+ do k = 1, nslyr
+ esnon(n) = esnon(n) + btrcrn(nt_qsno+k-1,n) &
+ * vsnon(n)/real(nslyr,kind=dbl_kind)
+ enddo ! k
+
+ vbrin(n) = vicen(n)
+ if (tr_brine) vbrin(n) = btrcrn(nt_fbri,n) * vicen(n)
+
+ enddo ! n
+
+ call column_sum (ncat, vicen(:), vice_init)
+ call column_sum (ncat, vsnon(:), vsno_init)
+ call column_sum (ncat, eicen(:), eice_init)
+ call column_sum (ncat, esnon(:), esno_init)
+ call column_sum (ncat, sicen(:), sice_init)
+ call column_sum (ncat, vbrin(:), vbri_init)
+
+ endif ! conservation check
+
+ do niter = 1, nitermax
+
+ !-----------------------------------------------------------------
+ ! Compute the thickness distribution of ridging ice
+ ! and various quantities associated with the new ridged ice.
+ !-----------------------------------------------------------------
+
+ call ridge_itd (ncat, aice0, &
+ aicen(:), vicen(:), &
+ krdg_partic, krdg_redist, &
+ mu_rdg, &
+ aksum, apartic(:), &
+ hrmin(:), hrmax(:), &
+ hrexp(:), krdg(:), &
+ aparticn(:), krdgn(:), &
+ mraftn(:))
+
+ !-----------------------------------------------------------------
+ ! Compute the area opening and closing.
+ !-----------------------------------------------------------------
+
+ closing_net = dice/dt ! net closing rate
+ divu_adv = -closing_net
+ opning = c0
+
+ !-----------------------------------------------------------------
+ ! Redistribute area, volume, and energy.
+ !-----------------------------------------------------------------
+
+ call ridge_shift (ntrcr, dt, &
+ ncat, hin_max(:), &
+ aicen(:), btrcrn(:,:), &
+ vicen(:), vsnon(:), &
+ aice0, trcr_depend, &
+ trcr_base, n_trcr_strata, &
+ nt_strata, krdg_redist, &
+ aksum, apartic(:), &
+ hrmin(:), hrmax(:), &
+ hrexp(:), krdg(:), &
+ closing_net, opning, &
+ ardg1, ardg2, &
+ virdg, aopen, &
+ ardg1n(:), ardg2n(:), &
+ virdgn(:), &
+ nslyr, n_aero, &
+ msnow_mlt, esnow_mlt, &
+ maero(:), mpond, &
+ l_stop, stop_label, &
+ aredistn(:), vredistn(:))
+ if (l_stop) return
+
+ !-----------------------------------------------------------------
+ ! Make sure the new area = 1. If not (because the closing
+ ! and opening rates were reduced above), prepare to ridge again
+ ! with new rates.
+ !-----------------------------------------------------------------
+
+ ! dcmod - adjust with berg area
+ call asum_ridging (ncat, aicen(:), aice0, asum)
+ asum = asum + bice
+
+ ! dcmod - if less than 1 due to forced ridging, add missing area back to open water
+ if (asum < c1) aice0 = aice0 + (c1 - asum)
+
+ ! check sum again
+ call asum_ridging (ncat, aicen(:), aice0, asum)
+ asum = asum + bice
+
+ if (abs(asum - c1) < puny) then
+ iterate_ridging = .false.
+ closing_net = c0
+ opning = c0
+ else
+ iterate_ridging = .true.
+ divu_adv = (c1 - asum) / dt
+ closing_net = max(c0, -divu_adv)
+ opning = max(c0, divu_adv)
+ endif
+
+ !-----------------------------------------------------------------
+ ! If done, exit. If not, prepare to ridge again.
+ !-----------------------------------------------------------------
+
+ if (iterate_ridging) then
+ write(warning,*) 'Repeat ridging (bergs), niter =', niter
+ call add_warning(warning)
+ ! ! dcmod debugging messages
+ ! write(warning,*) 'area (inc. bergs):', asum
+ ! call add_warning(warning)
+ ! write(warning,*) 'ice:', asum - bice - aice0
+ ! call add_warning(warning)
+ ! write(warning,*) 'open water:', aice0
+ ! call add_warning(warning)
+ else
+ ! exit rdg_iteration
+ exit
+ endif
+
+ if (niter == nitermax) then
+ write(warning,*) ' '
+ call add_warning(warning)
+ write(warning,*) 'Exceeded max number of ridging iterations'
+ call add_warning(warning)
+ write(warning,*) 'max =',nitermax
+ call add_warning(warning)
+ l_stop = .true.
+ stop_label = "ridge_ice: Exceeded max number of ridging iterations"
+ return
+ endif
+
+ enddo ! niter
+
+ !-----------------------------------------------------------------
+ ! Compute final values of conserved quantities.
+ ! Check for conservation (allowing for snow thrown into ocean).
+ !-----------------------------------------------------------------
+
+ if (l_conservation_check) then
+
+ do n = 1, ncat
+ eicen(n) = c0
+ esnon(n) = c0
+ sicen(n) = c0
+ vbrin(n) = c0
+
+ do k = 1, nilyr
+ eicen(n) = eicen(n) + btrcrn(nt_qice+k-1,n) &
+ * vicen(n)/real(nilyr,kind=dbl_kind)
+ enddo
+ do k = 1, nslyr
+ esnon(n) = esnon(n) + btrcrn(nt_qsno+k-1,n) &
+ * vsnon(n)/real(nslyr,kind=dbl_kind)
+ enddo
+
+ vbrin(n) = vicen(n)
+ if (tr_brine) vbrin(n) = btrcrn(nt_fbri,n) * vbrin(n)
+
+ enddo ! ncat
+
+ call column_sum (ncat, vicen(:), vice_final)
+ call column_sum (ncat, vsnon(:), vsno_final)
+ call column_sum (ncat, eicen(:), eice_final)
+ call column_sum (ncat, esnon(:), esno_final)
+ call column_sum (ncat, sicen(:), sice_final)
+ call column_sum (ncat, vbrin(:), vbri_final)
+
+ vsno_final = vsno_final + msnow_mlt/rhos
+ esno_final = esno_final + esnow_mlt
+
+ fieldid = 'vice, ridging'
+ call column_conservation_check (fieldid, &
+ vice_init, vice_final, &
+ puny, &
+ l_stop)
+ if (l_stop) return
+
+ fieldid = 'vsno, ridging'
+ call column_conservation_check (fieldid, &
+ vsno_init, vsno_final, &
+ puny, &
+ l_stop)
+ if (l_stop) return
+
+ fieldid = 'eice, ridging'
+ call column_conservation_check (fieldid, &
+ eice_init, eice_final, &
+ puny*Lfresh*rhoi, &
+ l_stop)
+ if (l_stop) return
+
+ fieldid = 'esno, ridging'
+ call column_conservation_check (fieldid, &
+ esno_init, esno_final, &
+ puny*Lfresh*rhos, &
+ l_stop)
+ if (l_stop) return
+
+ endif ! l_conservation_check
+
+ !-----------------------------------------------------------------
+ ! Compute ridging diagnostics.
+ !-----------------------------------------------------------------
+
+ dti = c1/dt
+
+ if (present(dardg1dt)) then
+ dardg1dt = ardg1*dti
+ endif
+ if (present(dardg2dt)) then
+ dardg2dt = ardg2*dti
+ endif
+ if (present(dvirdgdt)) then
+ dvirdgdt = virdg*dti
+ endif
+ if (present(opening)) then
+ opening = aopen*dti
+ endif
+
+ if (present(dardg1ndt)) then
+ do n = 1, ncat
+ dardg1ndt(n) = ardg1n(n)*dti
+ enddo
+ endif
+ if (present(dardg2ndt)) then
+ do n = 1, ncat
+ dardg2ndt(n) = ardg2n(n)*dti
+ enddo
+ endif
+ if (present(dvirdgndt)) then
+ do n = 1, ncat
+ dvirdgndt(n) = virdgn(n)*dti
+ enddo
+ endif
+ if (present(araftn)) then
+ do n = 1, ncat
+ araftn(n) = mraftn(n)*ardg2n(n)
+! araftn(n) = mraftn(n)*ardg1n(n)*p5
+ enddo
+ endif
+ if (present(vraftn)) then
+ do n = 1, ncat
+ vraftn(n) = mraftn(n)*virdgn(n)
+ enddo
+ endif
+
+ !-----------------------------------------------------------------
+ ! Update fresh water and heat fluxes due to snow melt.
+ !-----------------------------------------------------------------
+
+ if (present(fresh)) then
+ fresh = fresh + msnow_mlt*dti
+ endif
+ if (present(fhocn)) then
+ fhocn = fhocn + esnow_mlt*dti
+ endif
+ if (present(faero_ocn)) then
+ do it = 1, n_aero
+ faero_ocn(it) = faero_ocn(it) + maero(it)*dti
+ enddo
+ endif
+ if (present(fpond)) then
+ fpond = fpond - mpond ! units change later
+ endif
+
+ !-----------------------------------------------------------------
+ ! Check for fractional ice area > 1.
+ !-----------------------------------------------------------------
+
+ call asum_ridging (ncat, aicen(:), aice0, asum)
+ ! dcmod - adjust with berg area
+ asum = asum + bice
+ if (abs(asum - c1) > puny) then
+ l_stop = .true.
+
+ write(warning,*) ' '
+ call add_warning(warning)
+ write(warning,*) 'Ridging error (bergs): total area > 1'
+ call add_warning(warning)
+ write(warning,*) 'area:', asum
+ call add_warning(warning)
+ write(warning,*) 'n, aicen:'
+ call add_warning(warning)
+ write(warning,*) 0, aice0
+ call add_warning(warning)
+ do n = 1, ncat
+ write(warning,*) n, aicen(n)
+ call add_warning(warning)
+ enddo
+ return
+ endif
+
+ end subroutine ridge_ice_by_bergs
+
+!=======================================================================
+
+ end module ice_bergs_mechred
+
+!=======================================================================
diff --git a/src/core_seaice/column/ice_mechred.F90 b/src/core_seaice/column/ice_mechred.F90
index 668621fee7..daa59f1b02 100644
--- a/src/core_seaice/column/ice_mechred.F90
+++ b/src/core_seaice/column/ice_mechred.F90
@@ -46,7 +46,7 @@ module ice_mechred
save
private
- public :: ridge_ice, asum_ridging, ridge_itd
+ public :: ridge_ice, asum_ridging, ridge_itd, ridge_shift, ridge_prep
real (kind=dbl_kind), parameter :: &
Cs = p25 , & ! fraction of shear energy contrbtng to ridging
diff --git a/src/core_seaice/shared/mpas_seaice_bergs.F b/src/core_seaice/shared/mpas_seaice_bergs.F
index bfaa170948..1876b9186a 100644
--- a/src/core_seaice/shared/mpas_seaice_bergs.F
+++ b/src/core_seaice/shared/mpas_seaice_bergs.F
@@ -219,7 +219,7 @@ end subroutine seaice_run_berg_dynamics
!
!-----------------------------------------------------------------------
- subroutine seaice_run_berg_postdynamics(domain, clock, itimestep)
+ subroutine seaice_run_berg_postdynamics(domain, clock, itimestep, ierr)
use seaice_berg_decay, only: &
seaice_run_berg_decay
@@ -233,6 +233,9 @@ subroutine seaice_run_berg_postdynamics(domain, clock, itimestep)
use seaice_berg_velocity_solver, only: &
seaice_berg_forcing_for_ice
+ use seaice_column, only: &
+ seaice_column_ridging_by_bergs
+
type(domain_type), intent(inout) :: &
domain !< Input/Output:
@@ -242,14 +245,19 @@ subroutine seaice_run_berg_postdynamics(domain, clock, itimestep)
integer, intent(in) :: &
itimestep !< Input:
+ integer, intent(inout) :: &
+ ierr !< Input/Output:
+
logical, pointer :: &
config_use_berg_decay, & ! flag to turn on / off berg decay parameterization
config_berg_update_open_water, & ! flag to turn on / off bergs affecting open water
- config_berg_seaice_interaction ! flag to turn on / off bergs affecting open water
+ config_berg_seaice_interaction, &! flag to turn on / off bergs term in sea ice momentum equation
+ config_berg_ridge_seaice ! flag to turn on / off bergs causing sea ice to ridge
call MPAS_pool_get_config(domain % configs, "config_use_berg_decay", config_use_berg_decay)
call MPAS_pool_get_config(domain % configs, "config_berg_update_open_water", config_berg_update_open_water)
call MPAS_pool_get_config(domain % configs, "config_berg_seaice_interaction", config_berg_seaice_interaction)
+ call MPAS_pool_get_config(domain % configs, "config_berg_ridge_seaice", config_berg_ridge_seaice)
! berg mass and area state variables update
call mpas_timer_start("Berg update mass")
@@ -273,7 +281,14 @@ subroutine seaice_run_berg_postdynamics(domain, clock, itimestep)
call mpas_timer_stop("Berg sea ice interaction")
! update open water
+ call mpas_timer_start("Berg update open water")
if (config_berg_update_open_water) call seaice_berg_update_open_water_area(domain)
+ call mpas_timer_stop("Berg update open water")
+
+ ! berg sea ice interaction
+ call mpas_timer_start("Berg induced ridging")
+ if (config_berg_ridge_seaice) call seaice_column_ridging_by_bergs(domain, ierr)
+ call mpas_timer_stop("Berg induced ridging")
end subroutine seaice_run_berg_postdynamics
diff --git a/src/core_seaice/shared/mpas_seaice_column.F b/src/core_seaice/shared/mpas_seaice_column.F
index 064ae259c2..1c31ff9794 100644
--- a/src/core_seaice/shared/mpas_seaice_column.F
+++ b/src/core_seaice/shared/mpas_seaice_column.F
@@ -42,7 +42,8 @@ module seaice_column
seaice_column_reinitialize_diagnostics_dynamics, &
seaice_column_reinitialize_diagnostics_bgc, &
seaice_column_coupling_prep, &
- seaice_column_finalize
+ seaice_column_finalize, &
+ seaice_column_ridging_by_bergs
! tracer object
type, private :: ciceTracerObjectType
@@ -13265,6 +13266,319 @@ subroutine column_combine_snow_ice_tracers(domain)
end subroutine column_combine_snow_ice_tracers
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! seaice_column_ridging_by_bergs
+!
+!> \brief
+!> \author Darin Comeau, LANL
+!> \date 9th April 2018
+!> \details
+!>
+!
+!-----------------------------------------------------------------------
+
+ subroutine seaice_column_ridging_by_bergs(domain, ierr)
+
+ use ice_colpkg, only: &
+ ! colpkg_step_ridge, &
+ colpkg_get_warnings, &
+ colpkg_clear_warnings
+
+ use ice_bergs_mechred, only: &
+ ridge_ice_by_bergs
+
+ type(domain_type), intent(inout) :: domain
+
+ integer, intent(inout) :: ierr
+
+ type(block_type), pointer :: block
+
+ type(MPAS_pool_type), pointer :: &
+ mesh, &
+ icestate, &
+ tracers, &
+ tracers_aggregate, &
+ ponds, &
+ ocean_fluxes, &
+ ridging, &
+ aerosols, &
+ biogeochemistry, &
+ initial, &
+ velocity_solver, &
+ berg_state, &
+ berg_tracers_aggregate
+
+ ! configs
+ logical, pointer :: &
+ config_use_column_biogeochemistry, &
+ config_use_brine
+
+ real(kind=RKIND), pointer :: &
+ config_dt
+
+ ! integer, pointer :: &
+ ! config_dynamics_subcycle_number
+
+ character(len=strKIND), pointer :: &
+ config_ridging_participation_function, &
+ config_ridging_redistribution_function
+
+ real(kind=RKIND), pointer :: &
+ config_ridiging_efolding_scale
+
+ ! dimensions
+ integer, pointer :: &
+ nCellsSolve, &
+ nCategories, &
+ nIceLayers, &
+ nSnowLayers, &
+ nAerosols, &
+ nBioLayers
+
+ ! variables
+ real(kind=RKIND), dimension(:), pointer :: &
+ pondFreshWaterFlux, &
+ oceanFreshWaterFlux, &
+ oceanSaltFlux, &
+ oceanHeatFlux, &
+ iceAreaCell, &
+ ridgeConvergence, &
+ ridgeShear, &
+ openWaterArea, &
+ areaLossRidge, &
+ areaGainRidge, &
+ iceVolumeRidged, &
+ openingRateRidge, &
+ categoryThicknessLimits, &
+ zSalinityFlux, &
+ bergDisplacedArea, &
+ bergAreaCell
+
+ real(kind=RKIND), dimension(:,:), pointer :: &
+ oceanAerosolFlux, &
+ ridgeParticipationFunction, &
+ ratioRidgeThicknessToIce, &
+ fractionNewRidgeArea, &
+ fractionNewRidgeVolume, &
+ areaLossRidgeCategory, &
+ areaGainRidgeCategory, &
+ iceVolumeRidgedCategory, &
+ raftingIceArea, &
+ raftingIceVolume, &
+ oceanBioFluxes
+
+ real(kind=RKIND), dimension(:,:,:), pointer :: &
+ iceAreaCategory, &
+ iceVolumeCategory, &
+ snowVolumeCategory
+
+ ! integer, dimension(:,:), pointer :: &
+ ! newlyFormedIce
+
+ integer, dimension(:), pointer :: &
+ indexToCellID
+
+ real(kind=RKIND), pointer :: &
+ dynamicsTimeStep
+
+ ! local
+ integer :: &
+ iCell, &
+ iCategory
+
+ ! logical, dimension(:), allocatable :: &
+ ! newlyFormedIceLogical
+
+ ! real(kind=RKIND), dimension(:,:) :: &
+ ! bergTracerArray
+
+ logical :: &
+ abortFlag, &
+ setGetPhysicsTracers, &
+ setGetBGCTracers
+
+ character(len=strKIND) :: &
+ abortMessage, &
+ abortLocation
+
+ character(len=strKINDWarnings), dimension(:), allocatable :: &
+ warnings
+
+ block => domain % blocklist
+ do while (associated(block))
+
+ call MPAS_pool_get_subpool(block % structs, "mesh", mesh)
+ call MPAS_pool_get_subpool(block % structs, "tracers", tracers)
+ call MPAS_pool_get_subpool(block % structs, "tracers_aggregate", tracers_aggregate)
+ call MPAS_pool_get_subpool(block % structs, "icestate", icestate)
+ call MPAS_pool_get_subpool(block % structs, "ponds", ponds)
+ call MPAS_pool_get_subpool(block % structs, "ocean_fluxes", ocean_fluxes)
+ call MPAS_pool_get_subpool(block % structs, "ridging", ridging)
+ call MPAS_pool_get_subpool(block % structs, "aerosols", aerosols)
+ call MPAS_pool_get_subpool(block % structs, "biogeochemistry", biogeochemistry)
+ call MPAS_pool_get_subpool(block % structs, "initial", initial)
+ call MPAS_pool_get_subpool(block % structs, "velocity_solver", velocity_solver)
+ call MPAS_pool_get_subpool(block % structs, "berg_state", berg_state)
+ call MPAS_pool_get_subpool(block % structs, "berg_tracers_aggregate", berg_tracers_aggregate)
+
+ ! call MPAS_pool_get_config(block % configs, "config_dynamics_subcycle_number", config_dynamics_subcycle_number)
+ call MPAS_pool_get_config(block % configs, "config_use_column_biogeochemistry", config_use_column_biogeochemistry)
+ call MPAS_pool_get_config(block % configs, "config_use_brine", config_use_brine)
+ call MPAS_pool_get_config(block % configs, "config_ridging_participation_function", config_ridging_participation_function)
+ call MPAS_pool_get_config(block % configs, "config_ridging_redistribution_function", config_ridging_redistribution_function)
+ call MPAS_pool_get_config(block % configs, "config_ridiging_efolding_scale", config_ridiging_efolding_scale)
+
+ call MPAS_pool_get_array(velocity_solver, "dynamicsTimeStep", dynamicsTimeStep)
+
+ call MPAS_pool_get_dimension(mesh, "nCellsSolve", nCellsSolve)
+ call MPAS_pool_get_dimension(mesh, "nCategories", nCategories)
+ call MPAS_pool_get_dimension(mesh, "nIceLayers", nIceLayers)
+ call MPAS_pool_get_dimension(mesh, "nSnowLayers", nSnowLayers)
+ call MPAS_pool_get_dimension(mesh, "nAerosols", nAerosols)
+ ! call MPAS_pool_get_dimension(block % dimensions, "nBioLayers", nBioLayers)
+
+ call MPAS_pool_get_array(mesh, "indexToCellID", indexToCellID)
+
+ call MPAS_pool_get_array(tracers_aggregate, "iceAreaCell", iceAreaCell)
+
+ call MPAS_pool_get_array(icestate, "openWaterArea", openWaterArea)
+
+ call MPAS_pool_get_array(tracers, "iceAreaCategory", iceAreaCategory, 1)
+ call MPAS_pool_get_array(tracers, "iceVolumeCategory", iceVolumeCategory, 1)
+ call MPAS_pool_get_array(tracers, "snowVolumeCategory", snowVolumeCategory, 1)
+
+ call MPAS_pool_get_array(ocean_fluxes, "oceanFreshWaterFlux", oceanFreshWaterFlux)
+ call MPAS_pool_get_array(ocean_fluxes, "oceanSaltFlux", oceanSaltFlux)
+ call MPAS_pool_get_array(ocean_fluxes, "oceanHeatFlux", oceanHeatFlux)
+
+ call MPAS_pool_get_array(ridging, "ridgeConvergence", ridgeConvergence)
+ call MPAS_pool_get_array(ridging, "ridgeShear", ridgeShear)
+ call MPAS_pool_get_array(ridging, "areaLossRidge", areaLossRidge)
+ call MPAS_pool_get_array(ridging, "areaGainRidge", areaGainRidge)
+ call MPAS_pool_get_array(ridging, "iceVolumeRidged", iceVolumeRidged)
+ call MPAS_pool_get_array(ridging, "openingRateRidge", openingRateRidge)
+ call MPAS_pool_get_array(ridging, "ridgeParticipationFunction", ridgeParticipationFunction)
+ call MPAS_pool_get_array(ridging, "ratioRidgeThicknessToIce", ratioRidgeThicknessToIce)
+ call MPAS_pool_get_array(ridging, "fractionNewRidgeArea", fractionNewRidgeArea)
+ call MPAS_pool_get_array(ridging, "fractionNewRidgeVolume", fractionNewRidgeVolume)
+ call MPAS_pool_get_array(ridging, "areaLossRidgeCategory", areaLossRidgeCategory)
+ call MPAS_pool_get_array(ridging, "areaGainRidgeCategory", areaGainRidgeCategory)
+ call MPAS_pool_get_array(ridging, "iceVolumeRidgedCategory", iceVolumeRidgedCategory)
+ call MPAS_pool_get_array(ridging, "raftingIceArea", raftingIceArea)
+ call MPAS_pool_get_array(ridging, "raftingIceVolume", raftingIceVolume)
+
+ call MPAS_pool_get_array(aerosols, "oceanAerosolFlux", oceanAerosolFlux)
+
+ call MPAS_pool_get_array(ponds, "pondFreshWaterFlux", pondFreshWaterFlux)
+
+ ! call MPAS_pool_get_array(biogeochemistry, "newlyFormedIce", newlyFormedIce)
+ ! call MPAS_pool_get_array(biogeochemistry, "oceanBioFluxes", oceanBioFluxes)
+ ! call MPAS_pool_get_array(biogeochemistry, "zSalinityFlux", zSalinityFlux)
+
+ call MPAS_pool_get_array(initial, "categoryThicknessLimits", categoryThicknessLimits)
+
+ call MPAS_pool_get_array(berg_state, "bergDisplacedArea", bergDisplacedArea)
+
+ call MPAS_pool_get_array(berg_tracers_aggregate, "bergAreaCell", bergAreaCell)
+
+ ! ! newly formed ice
+ ! allocate(newlyFormedIceLogical(nCategories))
+
+ setGetPhysicsTracers = .true.
+ setGetBGCTracers = config_use_column_biogeochemistry
+
+ do iCell = 1, nCellsSolve
+!dcmod - not needed, for bgc?
+ ! newly formed ice
+ ! do iCategory = 1, nCategories
+ ! newlyFormedIceLogical(iCategory) = (newlyFormedIce(iCategory,iCell) == 1)
+ ! enddo ! iCategory
+
+ ! set the category tracer array
+ call set_cice_tracer_array_category(block, ciceTracerObject, iCell, setGetPhysicsTracers, setGetBGCTracers)
+
+ abortFlag = .false.
+ abortMessage = ""
+!dcmod temporary so as to not update
+ ! bergTracerArray = ciceTracerObject % tracerArrayCategory
+
+ call colpkg_clear_warnings()
+ call ridge_ice_by_bergs(&
+ dynamicsTimeStep, &
+ bergDisplacedArea(iCell), &
+ nCategories, &
+ nAerosols, &
+ nIceLayers, &
+ nSnowLayers, &
+ ciceTracerObject % nTracers, &
+ categoryThicknessLimits, &
+ iceAreaCategory(1,:,iCell), &
+! bergTracerArray, & ! trcrn !dcmod is there reason this needs to be temporary array?
+ ciceTracerObject % tracerArrayCategory, & ! trcrn !dcmod is there reason this needs to be temporary array?
+ iceVolumeCategory(1,:,iCell), &
+ snowVolumeCategory(1,:,iCell), &
+ openWaterArea(iCell), &
+ bergAreaCell(iCell), &
+ ciceTracerObject % parentIndex, & ! trcr_depend
+ ciceTracerObject % firstAncestorMask, & ! trcr_base
+ ciceTracerObject % ancestorNumber, & !n_trcr_strata
+ ciceTracerObject % ancestorIndices, & !nt_strata
+ abortFlag, & ! move these at bottom
+ abortMessage, &
+ config_cice_int("config_ridging_participation_function", config_ridging_participation_function), &
+ config_cice_int("config_ridging_redistribution_function", config_ridging_redistribution_function), &
+ config_ridiging_efolding_scale, & ! dcmod change
+ config_use_brine, &
+ areaLossRidge(iCell), &
+ areaGainRidge(iCell), &
+ iceVolumeRidged(iCell), &
+ openingRateRidge(iCell), &
+ pondFreshWaterFlux(iCell), &
+ oceanFreshWaterFlux(iCell), &
+ oceanHeatFlux(iCell), &
+ oceanAerosolFlux(:,iCell), &
+ ridgeParticipationFunction(:,iCell), &
+ ratioRidgeThicknessToIce(:,iCell), &
+ fractionNewRidgeArea(:,iCell), &
+ fractionNewRidgeVolume(:,iCell), &
+ areaLossRidgeCategory(:,iCell), &
+ areaGainRidgeCategory(:,iCell), &
+ iceVolumeRidgedCategory(:,iCell), &
+ raftingIceArea(:,iCell), &
+ raftingIceVolume(:,iCell))
+
+ call colpkg_get_warnings(warnings)
+ call column_write_warnings(warnings)
+
+ ! ! update
+ ! do iCategory = 1, nCategories
+ ! newlyFormedIce(iCategory,iCell) = 0
+ ! if (newlyFormedIceLogical(iCategory)) newlyFormedIce(iCategory,iCell) = 1
+ ! enddo ! iCategory
+
+ ! ! get category tracer array
+ ! call get_cice_tracer_array_category(block, ciceTracerObject, iCell, setGetPhysicsTracers, setGetBGCTracers)
+
+ ! code abort
+ if (abortFlag) then
+ call mpas_log_write("column_ridging_by_bergs: "//trim(abortMessage) , messageType=MPAS_LOG_ERR)
+ call mpas_log_write("iCell: $i", messageType=MPAS_LOG_ERR, intArgs=(/indexToCellID(iCell)/))
+ ierr = SEAICE_ERROR_COL_RIDGING
+ call seaice_critical_error_write_block(domain, block)
+ return
+ endif
+
+ enddo ! iCell
+
+ ! ! newly formed ice
+ ! deallocate(newlyFormedIceLogical)
+
+ block => block % next
+ enddo
+
+ end subroutine seaice_column_ridging_by_bergs
+
!-----------------------------------------------------------------------
! Warning messages
!-----------------------------------------------------------------------
diff --git a/src/core_seaice/shared/mpas_seaice_time_integration.F b/src/core_seaice/shared/mpas_seaice_time_integration.F
index 8270b847c1..107fbbceea 100644
--- a/src/core_seaice/shared/mpas_seaice_time_integration.F
+++ b/src/core_seaice/shared/mpas_seaice_time_integration.F
@@ -187,7 +187,7 @@ subroutine seaice_timestep(&
! berg decay and calving check
if (config_use_bergs) then
call mpas_timer_start("Berg post-dynamics")
- call seaice_run_berg_postdynamics(domain, clock, itimestep)
+ call seaice_run_berg_postdynamics(domain, clock, itimestep, ierr)
call mpas_timer_stop("Berg post-dynamics")
endif
diff --git a/testing_and_setup/seaice/configurations/icebergs/namelist.seaice b/testing_and_setup/seaice/configurations/icebergs/namelist.seaice
index 89c9a19deb..a9ff16a1b7 100644
--- a/testing_and_setup/seaice/configurations/icebergs/namelist.seaice
+++ b/testing_and_setup/seaice/configurations/icebergs/namelist.seaice
@@ -347,6 +347,7 @@
config_berg_advection_outside_velocity_subcycling = true
config_berg_update_open_water = true
config_berg_seaice_interaction = true
+ config_berg_ridge_seaice = true
/
&AM_highFrequencyOutput
config_AM_highFrequencyOutput_enable = false
diff --git a/testing_and_setup/seaice/configurations/standard_physics/namelist.seaice b/testing_and_setup/seaice/configurations/standard_physics/namelist.seaice
index 47ad16f388..a7930b2c1b 100644
--- a/testing_and_setup/seaice/configurations/standard_physics/namelist.seaice
+++ b/testing_and_setup/seaice/configurations/standard_physics/namelist.seaice
@@ -347,6 +347,7 @@
config_berg_advection_outside_velocity_subcycling = false
config_berg_update_open_water = false
config_berg_seaice_interaction = false
+ config_berg_ridge_seaice = false
/
&AM_highFrequencyOutput
config_AM_highFrequencyOutput_enable = false