From 1c8dc7e8f73aba6ba5109524fbd455c44b1ea2e7 Mon Sep 17 00:00:00 2001 From: Henri Drake Date: Wed, 14 Aug 2024 16:57:35 -0700 Subject: [PATCH 1/4] Attempt at separating sea ice melt/formation from liquid precipitation For unjustified legacy reasons, these have been bundled together in the past. They are not only bundled together in SIS2 variables, but also in the ice-ocean fluxes passed to the coupler and hence to MOM6. This commit attempts to address this Issue: https://github.com/NOAA-GFDL/SIS2/issues/213 However, I currently am getting a SEGFAULT error stating that `iobt%seaice_melt` has not been allocated memory when I try to do a checksum on it. --- src/SIS_slow_thermo.F90 | 9 ++------- src/SIS_sum_output.F90 | 1 + src/SIS_types.F90 | 5 ++++- src/combined_ice_ocean_driver.F90 | 10 ++++++---- src/ice_boundary_types.F90 | 6 +++--- src/ice_model.F90 | 4 +++- src/ice_type.F90 | 12 +++++++++--- 7 files changed, 28 insertions(+), 19 deletions(-) diff --git a/src/SIS_slow_thermo.F90 b/src/SIS_slow_thermo.F90 index d4bc6f0f..7f63788b 100644 --- a/src/SIS_slow_thermo.F90 +++ b/src/SIS_slow_thermo.F90 @@ -1322,7 +1322,8 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) ! With transmuted ice, the ice is non-conservatively changed to match the ocean properties. IOF%flux_salt(i,j) = IOF%flux_salt(i,j) + salt_to_ocn * (0.001*Idt_slow) - net_melt(i,j) = net_melt(i,j) + water_to_ocn * Idt_slow ! This goes to IOF%lprec_ocn_top. + net_melt(i,j) = net_melt(i,j) + water_to_ocn * Idt_slow + IOF%seaice_melt(i,j) = net_melt(i,j) IOF%Enth_mass_out_ocn(i,j) = IOF%Enth_mass_out_ocn(i,j) + heat_to_ocn ! With transmuted ice, the imbalances are stored to close the heat and salt budgets. @@ -1369,12 +1370,6 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) call disable_SIS_averaging(CS%diag) - ! Combine the liquid precipitation with the net melt of ice and snow for - ! passing to the ocean. These may later be kept separate. - do j=jsc,jec ; do i=isc,iec - IOF%lprec_ocn_top(i,j) = IOF%lprec_ocn_top(i,j) + net_melt(i,j) - enddo ; enddo - ! Make sure TrLay is no longer allocated if(allocated(TrLay)) deallocate(TrLay) end subroutine SIS2_thermodynamics diff --git a/src/SIS_sum_output.F90 b/src/SIS_sum_output.F90 index ecfdcfc3..575db44b 100644 --- a/src/SIS_sum_output.F90 +++ b/src/SIS_sum_output.F90 @@ -759,6 +759,7 @@ subroutine accumulate_bottom_input(IST, OSS, FIA, IOF, dt, G, US, IG, CS) do j=jsc,jec ; do i=isc,iec CS%water_in_col(i,j) = CS%water_in_col(i,j) - dt * & ( ((FIA%runoff(i,j) + FIA%calving(i,j)) + & + IOF%seaice_melt(i,j) + & (IOF%lprec_ocn_top(i,j) + IOF%fprec_ocn_top(i,j))) - IOF%evap_ocn_top(i,j) ) Flux_SW = 0.0 diff --git a/src/SIS_types.F90 b/src/SIS_types.F90 index 18d5809f..f941b6c9 100644 --- a/src/SIS_types.F90 +++ b/src/SIS_types.F90 @@ -358,6 +358,7 @@ module SIS_types flux_lh_ocn_top, & !< The upward flux of latent heat at the ocean surface [Q R Z T-1 ~> W m-2]. lprec_ocn_top, & !< The downward flux of liquid precipitation at the ocean surface [R Z T-1 ~> kg m-2 s-1]. fprec_ocn_top, & !< The downward flux of frozen precipitation at the ocean surface [R Z T-1 ~> kg m-2 s-1]. + seaice_melt, & !< The downward freshwater flux into the ocean due to sea ice melt [R Z T-1 ~> kg m-2 s-1]. flux_u_ocn, & !< The flux of x-momentum into the ocean at locations given by !! flux_uv_stagger [R Z L T-2 ~> Pa]. !! Note that regardless of the staggering, flux_u_ocn is allocated as though on an A-grid. @@ -925,6 +926,7 @@ subroutine alloc_ice_ocean_flux(IOF, HI, do_stress_mag, do_iceberg_fields, do_tr allocate(IOF%flux_sw_ocn(SZI_(HI), SZJ_(HI), NBANDS), source=0.0) allocate(IOF%lprec_ocn_top(SZI_(HI), SZJ_(HI)), source=0.0) allocate(IOF%fprec_ocn_top(SZI_(HI), SZJ_(HI)), source=0.0) + allocate(IOF%seaice_melt(SZI_(HI), SZJ_(HI)), source=0.0) allocate(IOF%flux_u_ocn(SZI_(HI), SZJ_(HI)), source=0.0) allocate(IOF%flux_v_ocn(SZI_(HI), SZJ_(HI)), source=0.0) if (alloc_stress_mag) then @@ -2171,7 +2173,7 @@ subroutine dealloc_ice_ocean_flux(IOF) deallocate(IOF%flux_sh_ocn_top, IOF%evap_ocn_top) deallocate(IOF%flux_lw_ocn_top, IOF%flux_lh_ocn_top) deallocate(IOF%flux_sw_ocn) - deallocate(IOF%lprec_ocn_top, IOF%fprec_ocn_top, IOF%flux_salt) + deallocate(IOF%lprec_ocn_top, IOF%fprec_ocn_top, IOF%seaice_melt, IOF%flux_salt) deallocate(IOF%flux_u_ocn, IOF%flux_v_ocn, IOF%pres_ocn_top, IOF%mass_ice_sn_p) if (allocated(IOF%stress_mag)) deallocate(IOF%stress_mag) if (allocated(IOF%transmutation_salt_flux)) deallocate(IOF%transmutation_salt_flux) @@ -2218,6 +2220,7 @@ subroutine IOF_chksum(mesg, IOF, G, US, mech_fluxes, thermo_fluxes) call hchksum(IOF%evap_ocn_top, trim(mesg)//" IOF%evap_ocn_top", G%HI, scale=US%RZ_T_to_kg_m2s) call hchksum(IOF%lprec_ocn_top, trim(mesg)//" IOF%lprec_ocn_top", G%HI, scale=US%RZ_T_to_kg_m2s) call hchksum(IOF%fprec_ocn_top, trim(mesg)//" IOF%fprec_ocn_top", G%HI, scale=US%RZ_T_to_kg_m2s) + call hchksum(IOF%seaice_melt, trim(mesg)//" IOF%seaice_melt", G%HI, scale=US%RZ_T_to_kg_m2s) endif if (do_mech) then call hchksum(IOF%flux_u_ocn, trim(mesg)//" IOF%flux_u_ocn", G%HI, scale=US%RZ_T_to_kg_m2s*US%L_T_to_m_s) diff --git a/src/combined_ice_ocean_driver.F90 b/src/combined_ice_ocean_driver.F90 index 3c159c84..410979b7 100644 --- a/src/combined_ice_ocean_driver.F90 +++ b/src/combined_ice_ocean_driver.F90 @@ -304,6 +304,7 @@ subroutine direct_flux_ice_to_IOB(Time, Ice, IOB, do_thermo) if (ASSOCIATED(IOB%lw_flux)) IOB%lw_flux(:,:) = Ice%flux_lw(:,:) if (ASSOCIATED(IOB%lprec)) IOB%lprec(:,:) = Ice%lprec(:,:) if (ASSOCIATED(IOB%fprec)) IOB%fprec(:,:) = Ice%fprec(:,:) + if (ASSOCIATED(IOB%seaice_melt)) IOB%seaice_melt(:,:) = Ice%seaice_melt(:,:) if (ASSOCIATED(IOB%runoff)) IOB%runoff(:,:) = Ice%runoff(:,:) if (ASSOCIATED(IOB%calving)) IOB%calving(:,:) = Ice%calving if (ASSOCIATED(IOB%runoff_hflx)) IOB%runoff_hflx(:,:) = Ice%runoff_hflx(:,:) @@ -326,10 +327,11 @@ subroutine direct_flux_ice_to_IOB(Time, Ice, IOB, do_thermo) call data_override('OCN', 'sw_flux_nir_dif', IOB%sw_flux_nir_dif, Time) call data_override('OCN', 'sw_flux_vis_dir', IOB%sw_flux_vis_dir, Time) call data_override('OCN', 'sw_flux_vis_dif', IOB%sw_flux_vis_dif, Time) - call data_override('OCN', 'lprec', IOB%lprec , Time) - call data_override('OCN', 'fprec', IOB%fprec , Time) - call data_override('OCN', 'runoff', IOB%runoff , Time) - call data_override('OCN', 'calving', IOB%calving , Time) + call data_override('OCN', 'lprec', IOB%lprec , Time) + call data_override('OCN', 'fprec', IOB%fprec , Time) + call data_override('OCN', 'seaice_melt', IOB%seaice_melt , Time) + call data_override('OCN', 'runoff', IOB%runoff , Time) + call data_override('OCN', 'calving', IOB%calving , Time) call data_override('OCN', 'runoff_hflx', IOB%runoff_hflx , Time) call data_override('OCN', 'calving_hflx', IOB%calving_hflx , Time) endif diff --git a/src/ice_boundary_types.F90 b/src/ice_boundary_types.F90 index b9ea3a82..ff8ef034 100644 --- a/src/ice_boundary_types.F90 +++ b/src/ice_boundary_types.F90 @@ -37,9 +37,9 @@ module ice_boundary_types frazil => NULL(), & !< The frazil heat rejected by the ocean [J m-2]. sea_level => NULL(), & !< The sea level after adjustment for any surface !! pressure that the ocean allows to be expressed [m]. - calving => NULL(), & !< The mass flux per unit area of the ice shelf to convert to - !!bergs [RZ_T ~> kg m-2 s-1]. - calving_hflx => NULL() !< Calving heat flux [Q R Z T-1 ~> W m-2]. + calving => NULL(), & !< The mass flux per unit area of the ice shelf to convert to + !!bergs [RZ_T ~> kg m-2 s-1]. + calving_hflx => NULL() !< Calving heat flux [Q R Z T-1 ~> W m-2]. real, dimension(:,:,:), pointer :: data =>NULL() !< S collective field for "named" fields above integer :: stagger = BGRID_NE !< A flag indicating how the velocities are staggered. integer :: xtype !< A flag indicating the exchange type, which may be set to diff --git a/src/ice_model.F90 b/src/ice_model.F90 index 1042270d..6eee20bf 100644 --- a/src/ice_model.F90 +++ b/src/ice_model.F90 @@ -623,6 +623,7 @@ subroutine set_ocean_top_fluxes(Ice, IST, IOF, FIA, OSS, G, US, IG, sCS) Ice%flux_sw_vis_dir(:,:) = 0.0 ; Ice%flux_sw_vis_dif(:,:) = 0.0 Ice%flux_lw(:,:) = 0.0 ; Ice%flux_lh(:,:) = 0.0 Ice%fprec(:,:) = 0.0 ; Ice%lprec(:,:) = 0.0 + Ice%seaice_melt(:,:) = 0.0 call coupler_type_rescale_data(Ice%ocean_fluxes, 0.0) i_off = LBOUND(Ice%flux_t,1) - G%isc ; j_off = LBOUND(Ice%flux_t,2) - G%jsc @@ -640,6 +641,7 @@ subroutine set_ocean_top_fluxes(Ice, IST, IOF, FIA, OSS, G, US, IG, sCS) Ice%flux_lh(i2,j2) = US%QRZ_T_to_W_m2*IOF%flux_lh_ocn_top(i,j) Ice%fprec(i2,j2) = US%RZ_T_to_kg_m2s*IOF%fprec_ocn_top(i,j) Ice%lprec(i2,j2) = US%RZ_T_to_kg_m2s*IOF%lprec_ocn_top(i,j) + Ice%seaice_melt(i2,j2) = US%RZ_T_to_kg_m2s*IOF%seaice_melt(i,j) Ice%runoff(i2,j2) = US%RZ_T_to_kg_m2s*FIA%runoff(i,j) Ice%calving(i2,j2) = US%RZ_T_to_kg_m2s*FIA%calving(i,j) Ice%runoff_hflx(i2,j2) = US%QRZ_T_to_W_m2*FIA%runoff_hflx(i,j) @@ -659,7 +661,7 @@ subroutine set_ocean_top_fluxes(Ice, IST, IOF, FIA, OSS, G, US, IG, sCS) if (allocated(IOF%melt_nudge)) then do j=jsc,jec ; do i=isc,iec i2 = i+i_off ; j2 = j+j_off! Use these to correct for indexing differences. - Ice%lprec(i2,j2) = Ice%lprec(i2,j2) + US%RZ_T_to_kg_m2s*IOF%melt_nudge(i,j) + Ice%seaice_melt(i2,j2) = Ice%seaice_melt(i2,j2) + US%RZ_T_to_kg_m2s*IOF%melt_nudge(i,j) enddo ; enddo endif diff --git a/src/ice_type.F90 b/src/ice_type.F90 index 1687de67..950c239e 100644 --- a/src/ice_type.F90 +++ b/src/ice_type.F90 @@ -97,9 +97,10 @@ module ice_type_mod flux_sw_vis_dif => NULL(), & !< The diffuse visible shortwave heat flux into the ocean [W m-2]. flux_sw_nir_dir => NULL(), & !< The direct near-infrared heat flux into the ocean [W m-2]. flux_sw_nir_dif => NULL(), & !< The diffuse near-infrared heat flux into the ocean [W m-2]. - flux_lh => NULL(), & !< The latent heat flux out of the ocean [W m-2]. - lprec => NULL(), & !< The liquid precipitation flux into the ocean [kg m-2]. - fprec => NULL(), & !< The frozen precipitation flux into the ocean [kg m-2]. + flux_lh => NULL(), & !< The latent heat flux out of the ocean [W m-2]. + lprec => NULL(), & !< The liquid precipitation flux into the ocean [kg m-2]. + fprec => NULL(), & !< The frozen precipitation flux into the ocean [kg m-2]. + seaice_melt => NULL(), & !< The sea ice meltwater flux into the ocean [kg m-2]. p_surf => NULL(), & !< The pressure at the ocean surface [Pa]. This may !! or may not include atmospheric pressure. runoff => NULL(), & !< Liquid runoff into the ocean [kg m-2]. @@ -198,6 +199,7 @@ subroutine ice_type_slow_reg_restarts(domain, CatIce, param_file, Ice, & call safe_alloc_ptr(Ice%flux_lh, isc, iec, jsc, jec) !NI call safe_alloc_ptr(Ice%lprec, isc, iec, jsc, jec) call safe_alloc_ptr(Ice%fprec, isc, iec, jsc, jec) + call safe_alloc_ptr(Ice%seaice_melt, isc, iec, jsc, jec) call safe_alloc_ptr(Ice%p_surf, isc, iec, jsc, jec) call safe_alloc_ptr(Ice%runoff, isc, iec, jsc, jec) call safe_alloc_ptr(Ice%calving, isc, iec, jsc, jec) @@ -233,6 +235,7 @@ subroutine ice_type_slow_reg_restarts(domain, CatIce, param_file, Ice, & call register_restart_field(Ice_restart, 'flux_lw', Ice%flux_lw) call register_restart_field(Ice_restart, 'lprec', Ice%lprec) call register_restart_field(Ice_restart, 'fprec', Ice%fprec) + call register_restart_field(Ice_restart, 'seaice_melt', Ice%seaice_melt) call register_restart_field(Ice_restart, 'runoff', Ice%runoff) call register_restart_field(Ice_restart, 'calving', Ice%calving) call register_restart_field(Ice_restart, 'runoff_hflx', Ice%runoff_hflx, mandatory=.false.) @@ -338,6 +341,7 @@ subroutine dealloc_Ice_arrays(Ice) if (associated(Ice%flux_lh)) deallocate(Ice%flux_lh) if (associated(Ice%lprec)) deallocate(Ice%lprec) if (associated(Ice%fprec)) deallocate(Ice%fprec) + if (associated(Ice%seaice_melt)) deallocate(Ice%seaice_melt) if (associated(Ice%p_surf)) deallocate(Ice%p_surf) if (associated(Ice%runoff)) deallocate(Ice%runoff) if (associated(Ice%calving)) deallocate(Ice%calving) @@ -419,6 +423,7 @@ subroutine Ice_public_type_chksum(mesg, Ice, check_fast, check_slow, check_rough call chksum(Ice%flux_lh, trim(mesg)//" Ice%flux_lh") call chksum(Ice%lprec, trim(mesg)//" Ice%lprec") call chksum(Ice%fprec, trim(mesg)//" Ice%fprec") + call chksum(Ice%seaice_melt, trim(mesg)//" Ice%seaice_melt") call chksum(Ice%p_surf, trim(mesg)//" Ice%p_surf") call chksum(Ice%calving, trim(mesg)//" Ice%calving") call chksum(Ice%runoff, trim(mesg)//" Ice%runoff") @@ -680,6 +685,7 @@ subroutine ice_data_type_chksum(mesg, timestep, Ice, init_call) chks = SIS_chksum(Ice%flux_lh ) ; if (root) write(outunit,100) 'ice_data_type%flux_lh ', chks chks = SIS_chksum(Ice%lprec ) ; if (root) write(outunit,100) 'ice_data_type%lprec ', chks chks = SIS_chksum(Ice%fprec ) ; if (root) write(outunit,100) 'ice_data_type%fprec ', chks + chks = SIS_chksum(Ice%seaice_melt ) ; if (root) write(outunit,100) 'ice_data_type%seaice_melt ', chks chks = SIS_chksum(Ice%p_surf ) ; if (root) write(outunit,100) 'ice_data_type%p_surf ', chks chks = SIS_chksum(Ice%runoff ) ; if (root) write(outunit,100) 'ice_data_type%runoff ', chks chks = SIS_chksum(Ice%calving ) ; if (root) write(outunit,100) 'ice_data_type%calving ', chks From ec559baa65931c2617c170e973aaaaf9f259c75f Mon Sep 17 00:00:00 2001 From: Henri Drake Date: Thu, 15 Aug 2024 12:39:03 -0700 Subject: [PATCH 2/4] Added CMOR diagnostic for sea ice melt --- src/SIS_slow_thermo.F90 | 25 ++++++++++--------------- src/SIS_types.F90 | 2 +- 2 files changed, 11 insertions(+), 16 deletions(-) diff --git a/src/SIS_slow_thermo.F90 b/src/SIS_slow_thermo.F90 index 7f63788b..a9511ea8 100644 --- a/src/SIS_slow_thermo.F90 +++ b/src/SIS_slow_thermo.F90 @@ -573,11 +573,8 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) tmp2d, & ! A temporary array for mass balance diagnostics [R Z yr-1 ~> kg m-2 yr-1] qflx_lim_ice, & ! Ice limiting heat flux [Q R Z T-1 ~> W m-2] qflx_res_ice, & ! Ice restoring heat flux [Q R Z T-1 ~> W m-2] - cool_nudge, & ! A heat flux out of the sea ice that + cool_nudge ! A heat flux out of the sea ice that ! acts to create sea-ice [Q R Z T-1 ~> W m-2]. - net_melt ! The net mass flux from the ice and snow into the - ! ocean due to melting and freezing integrated - ! across all categories [R Z T-1 ~> kg m-2 s-1]. real, dimension(SZI_(G),SZJ_(G),1:IG%CatIce) :: & heat_in, & ! The input heat [Q R Z ~> J m-2] enth_prev, & ! The previous column integrated enthalpy by category [Q R Z ~> J m-2] @@ -818,7 +815,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) call cpu_clock_begin(iceClock5) - snow_to_ice(:,:,:) = 0.0 ; net_melt(:,:) = 0.0 + snow_to_ice(:,:,:) = 0.0 ; bsnk(:,:) = 0.0 salt_change(:,:) = 0.0 h2o_change(:,:) = 0.0 @@ -876,7 +873,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) !$OMP parallel do default(none) shared(isc,iec,jsc,jec,ncat,G,US,IST,S_col0,NkIce,S_col,dt_slow, & !$OMP snow_to_ice,heat_in,I_NK,enth_prev,enth_mass_in_col,bsnk, & -!$OMP Idt_slow,salt_change,net_melt,LatHtFus,LatHtVap,IG,CS,OSS, & +!$OMP Idt_slow,salt_change,LatHtFus,LatHtVap,IG,CS,OSS, & !$OMP FIA,IOF,npassive,nb) & !$OMP private(mass_prev,enthalpy,enthalpy_ocean,Salin, & !$OMP heat_to_ocn,h2o_ice_to_ocn,h2o_ocn_to_ice, & @@ -1010,7 +1007,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) enddo IOF%flux_sw_ocn(i,j,VIS_DIF) = IOF%flux_sw_ocn(i,j,VIS_DIF) + IST%part_size(i,j,k) * & (sw_tot * FIA%sw_abs_ocn(i,j,k)) - net_melt(i,j) = net_melt(i,j) + IST%part_size(i,j,k) * & + IOF%seaice_melt(i,j) = IOF%seaice_melt(i,j) + IST%part_size(i,j,k) * & ((h2o_ice_to_ocn-h2o_ocn_to_ice)*Idt_slow) bsnk(i,j) = bsnk(i,j) - IST%part_size(i,j,k)*bablt*Idt_slow ! bot. melt. ablation @@ -1019,7 +1016,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) !$OMP parallel do default(none) shared(isc,iec,jsc,jec,ncat,npassive,G,US,IG,IST,S_col0,NkIce, & !$OMP S_col,dt_slow,snow_to_ice,heat_in,I_NK,enth_mass_in_col, & - !$OMP enth_prev,Idt_slow,bsnk,salt_change,net_melt,LatHtFus, & + !$OMP enth_prev,Idt_slow,bsnk,salt_change,LatHtFus, & !$OMP FIA,CS,OSS,IOF) & !$OMP private(mass_prev,enthalpy,enthalpy_ocean,Salin,heat_to_ocn, & !$OMP h2o_ice_to_ocn,h2o_ocn_to_ice,evap_from_ocn,salt_to_ice, & @@ -1156,7 +1153,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) ! IST%part_size(i,j,k) * enth_ocn_to_ice IOF%Enth_Mass_in_ocn(i,j) = IOF%Enth_Mass_in_ocn(i,j) + & IST%part_size(i,j,k) * (h2o_ocn_to_ice * enthalpy_ocean) - net_melt(i,j) = net_melt(i,j) - & + IOF%seaice_melt(i,j) = IOF%seaice_melt(i,j) - & (h2o_ocn_to_ice * IST%part_size(i,j,k)) * Idt_slow if (CS%column_check) then @@ -1222,7 +1219,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) IST%mH_ice(i,j,k) = frac_keep*IST%mH_ice(i,j,k) IST%mH_snow(i,j,k) = frac_keep*IST%mH_snow(i,j,k) enddo - net_melt(i,j) = net_melt(i,j) + h2o_ice_to_ocn * Idt_slow + IOF%seaice_melt(i,j) = IOF%seaice_melt(i,j) + h2o_ice_to_ocn * Idt_slow qflx_lim_ice(i,j) = enth_to_melt * Idt_slow IOF%Enth_Mass_out_ocn(i,j) = IOF%Enth_Mass_out_ocn(i,j) - enth_ice_to_ocn if (CS%ice_rel_salin > 0.0) then @@ -1322,8 +1319,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) ! With transmuted ice, the ice is non-conservatively changed to match the ocean properties. IOF%flux_salt(i,j) = IOF%flux_salt(i,j) + salt_to_ocn * (0.001*Idt_slow) - net_melt(i,j) = net_melt(i,j) + water_to_ocn * Idt_slow - IOF%seaice_melt(i,j) = net_melt(i,j) + IOF%seaice_melt(i,j) = IOF%seaice_melt(i,j) + water_to_ocn * Idt_slow IOF%Enth_mass_out_ocn(i,j) = IOF%Enth_mass_out_ocn(i,j) + heat_to_ocn ! With transmuted ice, the imbalances are stored to close the heat and salt budgets. @@ -1362,8 +1358,8 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) scale=US%RZ_T_to_kg_m2s*Idt_slow) if (CS%id_qflim>0) call post_data(CS%id_qflim, qflx_lim_ice, CS%diag) if (CS%id_qfres>0) call post_data(CS%id_qfres, qflx_res_ice, CS%diag) - if (CS%id_net_melt>0) call post_data(CS%id_net_melt, net_melt, CS%diag) - if (CS%id_CMOR_melt>0) call post_data(CS%id_CMOR_melt, net_melt, CS%diag) + if (CS%id_net_melt>0) call post_data(CS%id_net_melt, IOF%seaice_melt, CS%diag) + if (CS%id_CMOR_melt>0) call post_data(CS%id_CMOR_melt, IOF%seaice_melt, CS%diag) if (coupler_type_initialized(IOF%tr_flux_ocn_top)) & call coupler_type_send_data(IOF%tr_flux_ocn_top, CS%Time) @@ -1557,7 +1553,6 @@ subroutine SIS_slow_thermo_init(Time, G, US, IG, param_file, diag, CS, tracer_fl CS%id_CMOR_melt = register_diag_field('ice_model','fsitherm' ,diag%axesT1, Time, & 'water_flux_into_sea_water_due_to_sea_ice_thermodynamics', & 'kg m-2 s-1', conversion=US%RZ_T_to_kg_m2s, missing_value=missing) - if (CS%do_ice_restore) then CS%id_qfres = register_diag_field('ice_model', 'QFLX_RESTORE_ICE', diag%axesT1, Time, & 'Ice Restoring heat flux', 'W/m^2', conversion=US%QRZ_T_to_W_m2, missing_value=missing) diff --git a/src/SIS_types.F90 b/src/SIS_types.F90 index f941b6c9..12d87162 100644 --- a/src/SIS_types.F90 +++ b/src/SIS_types.F90 @@ -2220,7 +2220,7 @@ subroutine IOF_chksum(mesg, IOF, G, US, mech_fluxes, thermo_fluxes) call hchksum(IOF%evap_ocn_top, trim(mesg)//" IOF%evap_ocn_top", G%HI, scale=US%RZ_T_to_kg_m2s) call hchksum(IOF%lprec_ocn_top, trim(mesg)//" IOF%lprec_ocn_top", G%HI, scale=US%RZ_T_to_kg_m2s) call hchksum(IOF%fprec_ocn_top, trim(mesg)//" IOF%fprec_ocn_top", G%HI, scale=US%RZ_T_to_kg_m2s) - call hchksum(IOF%seaice_melt, trim(mesg)//" IOF%seaice_melt", G%HI, scale=US%RZ_T_to_kg_m2s) + call hchksum(IOF%seaice_melt, trim(mesg)//" IOF%seaice_melt", G%HI, scale=US%RZ_T_to_kg_m2s) endif if (do_mech) then call hchksum(IOF%flux_u_ocn, trim(mesg)//" IOF%flux_u_ocn", G%HI, scale=US%RZ_T_to_kg_m2s*US%L_T_to_m_s) From 5ecae717f0f211cba7f354fae3760232d33a0cfc Mon Sep 17 00:00:00 2001 From: Henri Drake Date: Wed, 11 Dec 2024 14:38:58 -0800 Subject: [PATCH 3/4] Reset IOF%seaice_melt to zero every timestep just like old net_melt This fixes a bug introduced in my previous two commits, in which I was unintentionally adding IOF%seaice_melt rate on top of the previous timesteps result, without realizing that it never gets set back to zero! --- src/SIS_slow_thermo.F90 | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/SIS_slow_thermo.F90 b/src/SIS_slow_thermo.F90 index a9511ea8..e74eb7e8 100644 --- a/src/SIS_slow_thermo.F90 +++ b/src/SIS_slow_thermo.F90 @@ -854,6 +854,11 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) IOF%lprec_ocn_top(i,j) = IOF%lprec_ocn_top(i,j) + IST%part_size(i,j,k) * FIA%lprec_top(i,j,k) enddo ; enddo ; enddo + ! Reset sea ice melt rate to zero + do j=jsc,jec ; do i=isc,iec + IOF%seaice_melt(i,j) = 0. + enddo ; enddo + ! Add fluxes of snow and other properties to the ocean due to recent ridging or drifting events. if (allocated(IST%snow_to_ocn)) then !$OMP do From 49a2829d4f9f6f1024b119eb4a35c1398e935f14 Mon Sep 17 00:00:00 2001 From: Henri Drake Date: Mon, 16 Dec 2024 16:08:13 -0800 Subject: [PATCH 4/4] Simpler way of resetting seaice_melt to zero --- src/SIS_slow_thermo.F90 | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/SIS_slow_thermo.F90 b/src/SIS_slow_thermo.F90 index e74eb7e8..6f0e8044 100644 --- a/src/SIS_slow_thermo.F90 +++ b/src/SIS_slow_thermo.F90 @@ -855,9 +855,7 @@ subroutine SIS2_thermodynamics(IST, dt_slow, CS, OSS, FIA, IOF, G, US, IG) enddo ; enddo ; enddo ! Reset sea ice melt rate to zero - do j=jsc,jec ; do i=isc,iec - IOF%seaice_melt(i,j) = 0. - enddo ; enddo + IOF%seaice_melt(:,:) = 0. ! Add fluxes of snow and other properties to the ocean due to recent ridging or drifting events. if (allocated(IST%snow_to_ocn)) then