diff --git a/config_src/drivers/solo_driver/MOM_driver.F90 b/config_src/drivers/solo_driver/MOM_driver.F90 index 1b88f1ce36..c2cd0a248c 100644 --- a/config_src/drivers/solo_driver/MOM_driver.F90 +++ b/config_src/drivers/solo_driver/MOM_driver.F90 @@ -26,6 +26,7 @@ program MOM_main use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end use MOM_cpu_clock, only : CLOCK_COMPONENT + use MOM_data_override, only : data_override_init use MOM_diag_mediator, only : enable_averaging, disable_averaging, diag_mediator_end use MOM_diag_mediator, only : diag_ctrl, diag_mediator_close_registration use MOM, only : initialize_MOM, step_MOM, MOM_control_struct, MOM_end @@ -47,6 +48,7 @@ program MOM_main use MOM_ice_shelf, only : initialize_ice_shelf, ice_shelf_end, ice_shelf_CS use MOM_ice_shelf, only : shelf_calc_flux, add_shelf_forces, ice_shelf_save_restart use MOM_ice_shelf, only : initialize_ice_shelf_fluxes, initialize_ice_shelf_forces + use MOM_ice_shelf, only : ice_shelf_query use MOM_interpolate, only : time_interp_external_init use MOM_io, only : file_exists, open_ASCII_file, close_file use MOM_io, only : check_nml_error, io_infra_init, io_infra_end @@ -176,6 +178,8 @@ program MOM_main type(surface_forcing_CS), pointer :: surface_forcing_CSp => NULL() type(write_cputime_CS), pointer :: write_CPU_CSp => NULL() type(ice_shelf_CS), pointer :: ice_shelf_CSp => NULL() + logical :: override_shelf_fluxes !< If true, and shelf dynamics are active, + !! the data_override feature is enabled (only for MOSAIC grid types) type(wave_parameters_cs), pointer :: waves_CSp => NULL() type(MOM_restart_CS), pointer :: & restart_CSp => NULL() !< A pointer to the restart control structure @@ -302,6 +306,8 @@ program MOM_main ! when using an ice shelf call initialize_ice_shelf_fluxes(ice_shelf_CSp, grid, US, fluxes) call initialize_ice_shelf_forces(ice_shelf_CSp, grid, US, forces) + call ice_shelf_query(ice_shelf_CSp, grid, data_override_shelf_fluxes=override_shelf_fluxes) + if (override_shelf_fluxes) call data_override_init(Ocean_Domain_in=grid%domain%mpp_domain) endif diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90 index be40f7bcb7..23b8b9da7c 100644 --- a/src/core/MOM_forcing_type.F90 +++ b/src/core/MOM_forcing_type.F90 @@ -170,6 +170,8 @@ module MOM_forcing_type !! exactly 0 away from shelves or on land. real, pointer, dimension(:,:) :: iceshelf_melt => NULL() !< Ice shelf melt rate (positive) !! or freezing (negative) [R Z T-1 ~> kg m-2 s-1] + real, pointer, dimension(:,:) :: shelf_sfc_mass_flux => NULL() !< Ice shelf surface mass flux + !! deposition from the atmosphere. [R Z T-1 ~> kg m-2 s-1] ! Scalars set by surface forcing modules real :: vPrecGlobalAdj = 0. !< adjustment to restoring vprec to zero out global net [kg m-2 s-1] @@ -377,7 +379,6 @@ module MOM_forcing_type ! Iceberg + Ice shelf diagnostic handles integer :: id_ustar_ice_cover = -1 integer :: id_frac_ice_cover = -1 - ! wave forcing diagnostics handles. integer :: id_lamult = -1 !>@} @@ -2099,6 +2100,12 @@ subroutine fluxes_accumulate(flux_tmp, fluxes, G, wt2, forces) fluxes%iceshelf_melt(i,j) = flux_tmp%iceshelf_melt(i,j) enddo ; enddo endif + if (associated(fluxes%shelf_sfc_mass_flux) & + .and. associated(flux_tmp%shelf_sfc_mass_flux)) then + do i=isd,ied ; do j=jsd,jed + fluxes%shelf_sfc_mass_flux(i,j) = flux_tmp%shelf_sfc_mass_flux(i,j) + enddo ; enddo + endif if (associated(fluxes%frac_shelf_h) .and. associated(flux_tmp%frac_shelf_h)) then do i=isd,ied ; do j=jsd,jed fluxes%frac_shelf_h(i,j) = flux_tmp%frac_shelf_h(i,j) @@ -2928,7 +2935,8 @@ end subroutine forcing_diagnostics !> Conditionally allocate fields within the forcing type subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, & - shelf, iceberg, salt, fix_accum_bug, cfc, waves) + shelf, iceberg, salt, fix_accum_bug, cfc, waves, & + shelf_sfc_accumulation) type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(forcing), intent(inout) :: fluxes !< A structure containing thermodynamic forcing fields logical, optional, intent(in) :: water !< If present and true, allocate water fluxes @@ -2942,14 +2950,21 @@ subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, & !! accumulation of ustar_gustless logical, optional, intent(in) :: cfc !< If present and true, allocate cfc fluxes logical, optional, intent(in) :: waves !< If present and true, allocate wave fields + logical, optional, intent(in) :: shelf_sfc_accumulation !< If present and true, and shelf is true, + !! then allocate surface flux deposition from the atmosphere + !! over ice shelves and ice sheets. ! Local variables integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB logical :: heat_water + logical :: shelf_sfc_acc isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB + shelf_sfc_acc=.false. + if (present(shelf_sfc_accumulation)) shelf_sfc_acc=shelf_sfc_accumulation + call myAlloc(fluxes%ustar,isd,ied,jsd,jed, ustar) call myAlloc(fluxes%ustar_gustless,isd,ied,jsd,jed, ustar) @@ -2987,9 +3002,13 @@ subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, & call myAlloc(fluxes%p_surf,isd,ied,jsd,jed, press) - call myAlloc(fluxes%frac_shelf_h,isd,ied,jsd,jed, shelf) - call myAlloc(fluxes%ustar_shelf,isd,ied,jsd,jed, shelf) - call myAlloc(fluxes%iceshelf_melt,isd,ied,jsd,jed, shelf) + ! These fields should only be allocated if ice shelf is enabled. + if (present(shelf)) then; if (shelf) then + call myAlloc(fluxes%frac_shelf_h,isd,ied,jsd,jed, shelf) + call myAlloc(fluxes%ustar_shelf,isd,ied,jsd,jed, shelf) + call myAlloc(fluxes%iceshelf_melt,isd,ied,jsd,jed, shelf) + if (shelf_sfc_acc) call myAlloc(fluxes%shelf_sfc_mass_flux,isd,ied,jsd,jed, shelf_sfc_acc) + endif; endif !These fields should only on allocated when iceberg area is being passed through the coupler. call myAlloc(fluxes%ustar_berg,isd,ied,jsd,jed, iceberg) @@ -3257,6 +3276,8 @@ subroutine deallocate_forcing_type(fluxes) if (associated(fluxes%ustar_tidal)) deallocate(fluxes%ustar_tidal) if (associated(fluxes%ustar_shelf)) deallocate(fluxes%ustar_shelf) if (associated(fluxes%iceshelf_melt)) deallocate(fluxes%iceshelf_melt) + if (associated(fluxes%shelf_sfc_mass_flux)) & + deallocate(fluxes%shelf_sfc_mass_flux) if (associated(fluxes%frac_shelf_h)) deallocate(fluxes%frac_shelf_h) if (associated(fluxes%ustar_berg)) deallocate(fluxes%ustar_berg) if (associated(fluxes%area_berg)) deallocate(fluxes%area_berg) @@ -3355,6 +3376,7 @@ subroutine rotate_forcing(fluxes_in, fluxes, turns) call rotate_array(fluxes_in%frac_shelf_h, turns, fluxes%frac_shelf_h) call rotate_array(fluxes_in%ustar_shelf, turns, fluxes%ustar_shelf) call rotate_array(fluxes_in%iceshelf_melt, turns, fluxes%iceshelf_melt) + call rotate_array(fluxes_in%shelf_sfc_mass_flux, turns, fluxes%shelf_sfc_mass_flux) endif if (do_iceberg) then @@ -3603,6 +3625,7 @@ subroutine homogenize_forcing(fluxes, G, GV, US) call homogenize_field_t(fluxes%frac_shelf_h, G) call homogenize_field_t(fluxes%ustar_shelf, G, tmp_scale=US%Z_to_m*US%s_to_T) call homogenize_field_t(fluxes%iceshelf_melt, G, tmp_scale=US%RZ_T_to_kg_m2s) + call homogenize_field_t(fluxes%shelf_sfc_mass_flux, G, tmp_scale=US%RZ_T_to_kg_m2s) endif if (do_iceberg) then diff --git a/src/ice_shelf/MOM_ice_shelf.F90 b/src/ice_shelf/MOM_ice_shelf.F90 index 13af5a936a..c1e5392d5e 100644 --- a/src/ice_shelf/MOM_ice_shelf.F90 +++ b/src/ice_shelf/MOM_ice_shelf.F90 @@ -9,6 +9,7 @@ module MOM_ice_shelf use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end use MOM_cpu_clock, only : CLOCK_COMPONENT, CLOCK_ROUTINE use MOM_coms, only : num_PEs +use MOM_data_override, only : data_override use MOM_diag_mediator, only : MOM_diag_ctrl=>diag_ctrl use MOM_IS_diag_mediator, only : post_data=>post_IS_data use MOM_IS_diag_mediator, only : register_diag_field=>register_MOM_IS_diag_field, safe_alloc_ptr @@ -160,6 +161,8 @@ module MOM_ice_shelf !! equation of state to use. logical :: active_shelf_dynamics !< True if the ice shelf mass changes as a result !! the dynamic ice-shelf model. + logical :: data_override_shelf_fluxes !< True if the ice shelf surface mass fluxes can be + !! written using the data_override feature (only for MOSAIC grids) logical :: override_shelf_movement !< If true, user code specifies the shelf movement !! instead of using the dynamic ice-shelf mode. logical :: isthermo !< True if the ice shelf can exchange heat and @@ -188,7 +191,8 @@ module MOM_ice_shelf id_h_shelf = -1, id_h_mask = -1, & id_surf_elev = -1, id_bathym = -1, & id_area_shelf_h = -1, & - id_ustar_shelf = -1, id_shelf_mass = -1, id_mass_flux = -1 + id_ustar_shelf = -1, id_shelf_mass = -1, id_mass_flux = -1, & + id_shelf_sfc_mass_flux = -1 !>@} integer :: id_read_mass !< An integer handle used in time interpolation of @@ -321,6 +325,11 @@ subroutine shelf_calc_flux(sfc_state_in, fluxes_in, Time, time_step, CS) G => CS%grid ; US => CS%US ISS => CS%ISS + if (CS%data_override_shelf_fluxes .and. CS%active_shelf_dynamics) then + call data_override(G%Domain, 'shelf_sfc_mass_flux', fluxes_in%shelf_sfc_mass_flux, CS%Time, & + scale=US%kg_m2s_to_RZ_T) + endif + if (CS%rotate_index) then allocate(sfc_state) call rotate_surface_state(sfc_state_in, CS%Grid_in, sfc_state, CS%Grid, CS%turns) @@ -730,6 +739,15 @@ subroutine shelf_calc_flux(sfc_state_in, fluxes_in, Time, time_step, CS) call hchksum(ISS%mass_shelf, "mass_shelf after change thickness using melt", G%HI, haloshift=0, & scale=US%RZ_to_kg_m2) endif + + call change_thickness_using_precip(CS, ISS, G, US, fluxes,US%s_to_T*time_step, Time) + + if (CS%debug) then + call hchksum(ISS%h_shelf, "h_shelf after change thickness using surf acc", G%HI, haloshift=0, scale=US%Z_to_m) + call hchksum(ISS%mass_shelf, "mass_shelf after change thickness using surf acc", G%HI, haloshift=0, & + scale=US%RZ_to_kg_m2) + endif + endif if (CS%debug) call MOM_forcing_chksum("Before add shelf flux", fluxes, G, CS%US, haloshift=0) @@ -753,6 +771,8 @@ subroutine shelf_calc_flux(sfc_state_in, fluxes_in, Time, time_step, CS) if (CS%id_shelf_mass > 0) call post_data(CS%id_shelf_mass, ISS%mass_shelf, CS%diag) if (CS%id_area_shelf_h > 0) call post_data(CS%id_area_shelf_h, ISS%area_shelf_h, CS%diag) if (CS%id_ustar_shelf > 0) call post_data(CS%id_ustar_shelf, fluxes%ustar_shelf, CS%diag) + if (CS%id_shelf_sfc_mass_flux > 0) call post_data(CS%id_shelf_sfc_mass_flux, fluxes%shelf_sfc_mass_flux, CS%diag) + if (CS%id_melt > 0) call post_data(CS%id_melt, fluxes%iceshelf_melt, CS%diag) if (CS%id_thermal_driving > 0) call post_data(CS%id_thermal_driving, (sfc_state%sst-ISS%tfreeze), CS%diag) if (CS%id_Sbdry > 0) call post_data(CS%id_Sbdry, Sbdry, CS%diag) @@ -1265,7 +1285,7 @@ subroutine initialize_ice_shelf(param_file, ocn_grid, Time, CS, diag, forces_in, allocate(CS%Grid) call MOM_domains_init(CS%Grid%domain, param_file, min_halo=wd_halos, symmetric=GRID_SYM_,& domain_name='MOM_Ice_Shelf_in') -! allocate(CS%Grid_in%HI) + !allocate(CS%Grid_in%HI) !call hor_index_init(CS%Grid%Domain, CS%Grid%HI, param_file, & ! local_indexing=.not.global_indexing) call MOM_grid_init(CS%Grid, param_file, CS%US) @@ -1354,6 +1374,11 @@ subroutine initialize_ice_shelf(param_file, ocn_grid, Time, CS, diag, forces_in, "If true, user provided code specifies the ice-shelf "//& "movement instead of the dynamic ice model.", default=.false.) CS%active_shelf_dynamics = .not.CS%override_shelf_movement + call get_param(param_file, mdl, "DATA_OVERRIDE_SHELF_FLUXES", & + CS%data_override_shelf_fluxes, & + "If true, the data override feature is used to write "//& + "the surface mass flux deposition. This option is only "//& + "available for MOSAIC grid types.", default=.false.) call get_param(param_file, mdl, "GROUNDING_LINE_INTERPOLATE", CS%GL_regularize, & "If true, regularize the floatation condition at the "//& "grounding line as in Goldberg Holland Schoof 2009.", default=.false.) @@ -1815,6 +1840,8 @@ subroutine initialize_ice_shelf(param_file, ocn_grid, Time, CS, diag, forces_in, if (CS%active_shelf_dynamics) then CS%id_h_mask = register_diag_field('ice_shelf_model', 'h_mask', CS%diag%axesT1, CS%Time, & 'ice shelf thickness mask', 'none') + CS%id_shelf_sfc_mass_flux = register_diag_field('ice_shelf_model', 'sfc_mass_flux', CS%diag%axesT1, CS%Time, & + 'ice shelf surface mass flux deposition from atmosphere', 'none', conversion=US%RZ_T_to_kg_m2s) endif call MOM_IS_diag_mediator_close_registration(CS%diag) @@ -1845,7 +1872,7 @@ subroutine initialize_ice_shelf_fluxes(CS, ocn_grid, US, fluxes_in) ! when SHELF_THERMO = True. These fluxes are necessary if one wants to ! use either ENERGETICS_SFC_PBL (ALE mode) or BULKMIXEDLAYER (layer mode). call allocate_forcing_type(CS%Grid_in, fluxes_in, ustar=.true., shelf=.true., & - press=.true., water=CS%isthermo, heat=CS%isthermo) + press=.true., water=CS%isthermo, heat=CS%isthermo, shelf_sfc_accumulation = CS%active_shelf_dynamics) else call MOM_mesg("MOM_ice_shelf.F90, initialize_ice_shelf: allocating fluxes in solo mode.") call allocate_forcing_type(CS%Grid_in, fluxes_in, ustar=.true., shelf=.true., press=.true.) @@ -1975,6 +2002,57 @@ subroutine initialize_shelf_mass(G, param_file, CS, ISS, new_sim) end select end subroutine initialize_shelf_mass +!> This subroutine applies net accumulation/ablation at the top surface to the dynamic ice shelf. +!>>acc_rate[m-s]=surf_mass_flux/density_ice is ablation/accumulation rate +!>>positive for accumulation negative for ablation +subroutine change_thickness_using_precip(CS, ISS, G, US, fluxes,time_step, Time) + type(ice_shelf_CS), intent(in) :: CS !< A pointer to the ice shelf control structure + type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure. + type(ice_shelf_state), intent(inout) :: ISS !< A structure with elements that describe + !! the ice-shelf state + type(forcing), intent(in) :: fluxes !< A structure of surface fluxes that + !! includes surface mass flux + type(time_type), intent(in) :: Time !< The current model time + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + real, intent(in) :: time_step !< The time step for this update [T ~> s]. + + ! locals + integer :: i, j + real ::I_rho_ice + + I_rho_ice = 1.0 / CS%density_ice + + !update time +! CS%Time = Time + + +! CS%time_step = time_step + ! update surface mass flux rate +! if (CS%surf_mass_flux_from_file) call update_surf_mass_flux(G, US, CS, ISS, Time) + + do j=G%jsc,G%jec ; do i=G%isc,G%iec + if ((ISS%hmask(i,j) == 1) .or. (ISS%hmask(i,j) == 2)) then + + if (-fluxes%shelf_sfc_mass_flux(i,j) * time_step < ISS%h_shelf(i,j)) then + ISS%h_shelf(i,j) = ISS%h_shelf(i,j) + fluxes%shelf_sfc_mass_flux(i,j) * time_step * I_rho_ice + else + ! the ice is about to ablate, so set thickness, area, and mask to zero + ! NOTE: this is not mass conservative should maybe scale salt & heat flux for this cell + ISS%h_shelf(i,j) = 0.0 + ISS%hmask(i,j) = 0.0 + ISS%area_shelf_h(i,j) = 0.0 + endif + ISS%mass_shelf(i,j) = ISS%h_shelf(i,j) * CS%density_ice + endif + enddo ; enddo + + call pass_var(ISS%area_shelf_h, G%domain, complete=.false.) + call pass_var(ISS%h_shelf, G%domain, complete=.false.) + call pass_var(ISS%hmask, G%domain, complete=.false.) + call pass_var(ISS%mass_shelf, G%domain, complete=.true.) + +end subroutine change_thickness_using_precip + !> Updates the ice shelf mass using data from a file. subroutine update_shelf_mass(G, US, CS, ISS, Time) @@ -2032,12 +2110,13 @@ subroutine update_shelf_mass(G, US, CS, ISS, Time) end subroutine update_shelf_mass !> Save the ice shelf restart file -subroutine ice_shelf_query(CS, G, frac_shelf_h, mass_shelf) +subroutine ice_shelf_query(CS, G, frac_shelf_h, mass_shelf, data_override_shelf_fluxes) type(ice_shelf_CS), pointer :: CS !< ice shelf control structure type(ocean_grid_type), intent(in) :: G !< A pointer to an ocean grid control structure. real, optional, dimension(SZI_(G),SZJ_(G)), intent(out) :: frac_shelf_h !< Ice shelf area fraction [nodim]. real, optional, dimension(SZI_(G),SZJ_(G)), intent(out) :: mass_shelf ! kg m-2] - + logical, optional :: data_override_shelf_fluxes !< If true, shelf fluxes can be written using + !! the data_override capability (only for MOSAIC grids) integer :: i, j @@ -2055,6 +2134,11 @@ subroutine ice_shelf_query(CS, G, frac_shelf_h, mass_shelf) enddo ; enddo endif + if (present(data_override_shelf_fluxes)) then + data_override_shelf_fluxes=.false. + if (CS%active_shelf_dynamics) data_override_shelf_fluxes = CS%data_override_shelf_fluxes + endif + end subroutine ice_shelf_query !> Save the ice shelf restart file diff --git a/src/ice_shelf/MOM_ice_shelf_dynamics.F90 b/src/ice_shelf/MOM_ice_shelf_dynamics.F90 index 8fb674e36c..5eccf92976 100644 --- a/src/ice_shelf/MOM_ice_shelf_dynamics.F90 +++ b/src/ice_shelf/MOM_ice_shelf_dynamics.F90 @@ -895,7 +895,7 @@ subroutine ice_shelf_solve_outer(CS, ISS, G, US, u_shlf, v_shlf, taudx, taudy, i enddo call calc_shelf_driving_stress(CS, ISS, G, US, taudx, taudy, CS%OD_av) -! call pass_vector(taudx, taudy, G%domain, TO_ALL, BGRID_NE) + call pass_vector(taudx, taudy, G%domain, TO_ALL, BGRID_NE) ! This is to determine which cells contain the grounding line, the criterion being that the cell ! is ice-covered, with some nodes floating and some grounded flotation condition is estimated by ! assuming topography is cellwise constant and H is bilinear in a cell; floating where @@ -1816,13 +1816,14 @@ subroutine calc_shelf_driving_stress(CS, ISS, G, US, taudx, taudy, OD) real :: dxh, dyh,Dx,Dy ! Local grid spacing [L ~> m] real :: grav ! The gravitational acceleration [L2 Z-1 T-2 ~> m s-2] - integer :: i, j, iscq, iecq, jscq, jecq, isd, jsd, is, js, iegq, jegq + integer :: i, j, iscq, iecq, jscq, jecq, isd, jsd, ied, jed, is, js, iegq, jegq integer :: giec, gjec, gisc, gjsc, cnt, isc, jsc, iec, jec integer :: i_off, j_off isc = G%isc ; jsc = G%jsc ; iec = G%iec ; jec = G%jec iscq = G%iscB ; iecq = G%iecB ; jscq = G%jscB ; jecq = G%jecB isd = G%isd ; jsd = G%jsd + isd = G%isd ; jsd = G%jsd ; ied = G%ied ; jed = G%jed iegq = G%iegB ; jegq = G%jegB ! gisc = G%domain%nihalo+1 ; gjsc = G%domain%njhalo+1 gisc = 1 ; gjsc = 1 @@ -1842,6 +1843,7 @@ subroutine calc_shelf_driving_stress(CS, ISS, G, US, taudx, taudy, OD) BASE(:,:) = -CS%bed_elev(:,:) + OD(:,:) S(:,:) = -CS%bed_elev(:,:) + ISS%h_shelf(:,:) ! check whether the ice is floating or grounded + do j=jsc-G%domain%njhalo,jec+G%domain%njhalo do i=isc-G%domain%nihalo,iec+G%domain%nihalo if (rhoi_rhow * ISS%h_shelf(i,j) - CS%bed_elev(i,j) <= 0) then @@ -1851,6 +1853,13 @@ subroutine calc_shelf_driving_stress(CS, ISS, G, US, taudx, taudy, OD) endif enddo enddo + + call pass_var(S, G%domain) + allocate(Phi(1:8,1:4,isd:ied,jsd:jed), source=0.0) + do j=jscq,jecq ; do i=iscq,iecq + call bilinear_shape_fn_grid(G, i, j, Phi(:,:,i,j)) + enddo ; enddo + do j=jsc-1,jec+1 do i=isc-1,iec+1 cnt = 0 @@ -1996,6 +2005,8 @@ subroutine calc_shelf_driving_stress(CS, ISS, G, US, taudx, taudy, OD) endif enddo enddo + + deallocate(Phi) end subroutine calc_shelf_driving_stress subroutine init_boundary_values(CS, G, time, hmask, input_flux, input_thick, new_sim) @@ -2592,7 +2603,8 @@ subroutine calc_shelf_visc(CS, ISS, G, US, u_shlf, v_shlf) allocate(Phi(1:8,1:4,isd:ied,jsd:jed), source=0.0) - do j=jsc,jec ; do i=isc,iec +! do j=jsc,jec ; do i=isc,iec + do j=jscq,jecq ; do i=iscq,iecq call bilinear_shape_fn_grid(G, i, j, Phi(:,:,i,j)) enddo ; enddo diff --git a/src/ice_shelf/MOM_ice_shelf_initialize.F90 b/src/ice_shelf/MOM_ice_shelf_initialize.F90 index 7cc3c020a3..751a6953e6 100644 --- a/src/ice_shelf/MOM_ice_shelf_initialize.F90 +++ b/src/ice_shelf/MOM_ice_shelf_initialize.F90 @@ -395,8 +395,6 @@ end subroutine initialize_ice_shelf_boundary_channel !> Initialize ice shelf flow from file -!subroutine initialize_ice_flow_from_file(bed_elev,u_shelf, v_shelf,float_cond,& -! hmask,h_shelf, G, US, PF) subroutine initialize_ice_flow_from_file(bed_elev,u_shelf, v_shelf,float_cond,& G, US, PF) type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure @@ -410,12 +408,6 @@ subroutine initialize_ice_flow_from_file(bed_elev,u_shelf, v_shelf,float_cond,& real, dimension(SZDI_(G),SZDJ_(G)), & intent(inout) :: float_cond !< An array indicating where the ice !! shelf is floating: 0 if floating, 1 if not. -! real, dimension(SZDI_(G),SZDJ_(G)), & -! intent(in) :: hmask !< A mask indicating which tracer points are - !! partly or fully covered by an ice-shelf -! real, dimension(SZDI_(G),SZDJ_(G)), & -! intent(in) :: h_shelf !< A mask indicating which tracer points are -! !! partly or fully covered by an ice-shelf type(unit_scale_type), intent(in) :: US !< A structure containing unit conversion factors type(param_file_type), intent(in) :: PF !< A structure to parse for run-time parameters @@ -508,7 +500,7 @@ subroutine initialize_ice_shelf_boundary_from_file(u_face_mask_bdry, v_face_mask character(len=200) :: filename, bc_file, inputdir, icethick_file ! Strings for file/path character(len=200) :: ufcmskbdry_varname, vfcmskbdry_varname, & ubdryv_varname, vbdryv_varname, umask_varname, vmask_varname, & - h_varname,hmsk_varname ! Variable name in file + hmsk_varname ! Variable name in file character(len=40) :: mdl = "initialize_ice_shelf_boundary_from_file" ! This subroutine's name. integer :: i, j, isc, jsc, iec, jec @@ -526,9 +518,6 @@ subroutine initialize_ice_shelf_boundary_from_file(u_face_mask_bdry, v_face_mask call get_param(PF, mdl, "ICE_THICKNESS_FILE", icethick_file, & "The file from which the ice-shelf thickness is read.", & default="ice_shelf_thick.nc") -! call get_param(PF, mdl, "ICE_THICKNESS_VARNAME", h_varname, & -! "The name of the thickness variable in ICE_THICKNESS_FILE.", & -! default="h_shelf") call get_param(PF, mdl, "ICE_THICKNESS_MASK_VARNAME", hmsk_varname, & "The name of the icethickness mask variable in ICE_THICKNESS_FILE.", & default="h_mask") @@ -565,7 +554,6 @@ subroutine initialize_ice_shelf_boundary_from_file(u_face_mask_bdry, v_face_mask call MOM_read_data(filename,trim(vmask_varname), vmask, G%Domain, position=CORNER,scale=1.) filename = trim(inputdir)//trim(icethick_file) -! call MOM_read_data(filename, trim(h_varname), h_shelf, G%Domain, scale=US%m_to_Z) call MOM_read_data(filename,trim(hmsk_varname), hmask, G%Domain, scale=1.) isc = G%isc ; jsc = G%jsc ; iec = G%iec ; jec = G%jec