diff --git a/src/core/MOM.F90 b/src/core/MOM.F90
index 25f4f27ee7..392c603de6 100644
--- a/src/core/MOM.F90
+++ b/src/core/MOM.F90
@@ -89,7 +89,7 @@ module MOM
 use MOM_forcing_type,          only : copy_common_forcing_fields, set_derived_forcing_fields
 use MOM_forcing_type,          only : homogenize_forcing, homogenize_mech_forcing
 use MOM_grid,                  only : ocean_grid_type, MOM_grid_init, MOM_grid_end
-use MOM_grid,                  only : set_first_direction, rescale_grid_bathymetry
+use MOM_grid,                  only : set_first_direction
 use MOM_hor_index,             only : hor_index_type, hor_index_init
 use MOM_hor_index,             only : rotate_hor_index
 use MOM_interface_heights,     only : find_eta, calc_derived_thermo, thickness_to_dz
diff --git a/src/core/MOM_grid.F90 b/src/core/MOM_grid.F90
index 2e413e505b..9cbf420b23 100644
--- a/src/core/MOM_grid.F90
+++ b/src/core/MOM_grid.F90
@@ -15,7 +15,7 @@ module MOM_grid
 #include <MOM_memory.h>
 
 public MOM_grid_init, MOM_grid_end, set_derived_metrics, set_first_direction
-public isPointInCell, hor_index_type, get_global_grid_size, rescale_grid_bathymetry
+public isPointInCell, hor_index_type, get_global_grid_size
 
 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
@@ -400,40 +400,6 @@ subroutine MOM_grid_init(G, param_file, US, HI, global_indexing, bathymetry_at_v
 
 end subroutine MOM_grid_init
 
-!> rescale_grid_bathymetry permits a change in the internal units for the bathymetry on the grid,
-!! both rescaling the depths and recording the new internal units.
-subroutine rescale_grid_bathymetry(G, m_in_new_units)
-  type(ocean_grid_type), intent(inout) :: G    !< The horizontal grid structure
-  real,                  intent(in)    :: m_in_new_units !< The new internal representation of 1 m depth.
-  !### It appears that this routine is never called.
-
-  ! Local variables
-  real :: rescale ! A unit rescaling factor [various combinations of units ~> 1]
-  integer :: i, j, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
-
-  isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed
-  IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB
-
-  if (m_in_new_units == 1.0) return
-  if (m_in_new_units < 0.0) &
-    call MOM_error(FATAL, "rescale_grid_bathymetry: Negative depth units are not permitted.")
-  if (m_in_new_units == 0.0) &
-    call MOM_error(FATAL, "rescale_grid_bathymetry: Zero depth units are not permitted.")
-
-  rescale = 1.0 / m_in_new_units
-  do j=jsd,jed ; do i=isd,ied
-    G%bathyT(i,j) = rescale*G%bathyT(i,j)
-  enddo ; enddo
-  if (G%bathymetry_at_vel) then ; do j=jsd,jed ; do I=IsdB,IedB
-    G%Dblock_u(I,j) = rescale*G%Dblock_u(I,j) ; G%Dopen_u(I,j) = rescale*G%Dopen_u(I,j)
-  enddo ; enddo ; endif
-  if (G%bathymetry_at_vel) then ; do J=JsdB,JedB ; do i=isd,ied
-    G%Dblock_v(i,J) = rescale*G%Dblock_v(i,J) ; G%Dopen_v(i,J) = rescale*G%Dopen_v(i,J)
-  enddo ; enddo ; endif
-  G%max_depth = rescale*G%max_depth
-
-end subroutine rescale_grid_bathymetry
-
 !> set_derived_metrics calculates metric terms that are derived from other metrics.
 subroutine set_derived_metrics(G, US)
   type(ocean_grid_type), intent(inout) :: G  !< The horizontal grid structure