diff --git a/coupled_AM2_LM3_SIS2/AM2_SIS2_MOM6i_1deg/MOM_parameter_doc.all b/coupled_AM2_LM3_SIS2/AM2_SIS2_MOM6i_1deg/MOM_parameter_doc.all
index 30a767ede..f98321d1c 100644
--- a/coupled_AM2_LM3_SIS2/AM2_SIS2_MOM6i_1deg/MOM_parameter_doc.all
+++ b/coupled_AM2_LM3_SIS2/AM2_SIS2_MOM6i_1deg/MOM_parameter_doc.all
@@ -1157,6 +1157,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/coupled_AM2_LM3_SIS2/Concurrent_ice_1deg/MOM_parameter_doc.all b/coupled_AM2_LM3_SIS2/Concurrent_ice_1deg/MOM_parameter_doc.all
index a002ee261..47272b9f9 100644
--- a/coupled_AM2_LM3_SIS2/Concurrent_ice_1deg/MOM_parameter_doc.all
+++ b/coupled_AM2_LM3_SIS2/Concurrent_ice_1deg/MOM_parameter_doc.all
@@ -1157,6 +1157,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/coupled_AM2_LM3_SIS2/Intersperse_ice_1deg/MOM_parameter_doc.all b/coupled_AM2_LM3_SIS2/Intersperse_ice_1deg/MOM_parameter_doc.all
index 45ee60e08..df1821afd 100644
--- a/coupled_AM2_LM3_SIS2/Intersperse_ice_1deg/MOM_parameter_doc.all
+++ b/coupled_AM2_LM3_SIS2/Intersperse_ice_1deg/MOM_parameter_doc.all
@@ -1157,6 +1157,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/Baltic/MOM_parameter_doc.all b/ice_ocean_SIS2/Baltic/MOM_parameter_doc.all
index 0604cd82f..5f1254171 100644
--- a/ice_ocean_SIS2/Baltic/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/Baltic/MOM_parameter_doc.all
@@ -1157,6 +1157,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/Baltic_ALE_z_offline_tracers/MOM_parameter_doc.all b/ice_ocean_SIS2/Baltic_ALE_z_offline_tracers/MOM_parameter_doc.all
index 72ea067fc..4498d0bb1 100644
--- a/ice_ocean_SIS2/Baltic_ALE_z_offline_tracers/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/Baltic_ALE_z_offline_tracers/MOM_parameter_doc.all
@@ -1368,6 +1368,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/Baltic_OM4_025/MOM_parameter_doc.all b/ice_ocean_SIS2/Baltic_OM4_025/MOM_parameter_doc.all
index de2626dd8..05386b1de 100644
--- a/ice_ocean_SIS2/Baltic_OM4_025/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/Baltic_OM4_025/MOM_parameter_doc.all
@@ -1392,6 +1392,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/Baltic_OM4_05/MOM_parameter_doc.all b/ice_ocean_SIS2/Baltic_OM4_05/MOM_parameter_doc.all
index 259096e04..cf831828f 100644
--- a/ice_ocean_SIS2/Baltic_OM4_05/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/Baltic_OM4_05/MOM_parameter_doc.all
@@ -1411,6 +1411,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/OM4_025/MOM_parameter_doc.all b/ice_ocean_SIS2/OM4_025/MOM_parameter_doc.all
index d30ceb9b3..62b2530d8 100644
--- a/ice_ocean_SIS2/OM4_025/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/OM4_025/MOM_parameter_doc.all
@@ -1392,6 +1392,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/OM4_05/MOM_parameter_doc.all b/ice_ocean_SIS2/OM4_05/MOM_parameter_doc.all
index a2599794b..e30269e98 100644
--- a/ice_ocean_SIS2/OM4_05/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/OM4_05/MOM_parameter_doc.all
@@ -1411,6 +1411,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/OM_1deg/MOM_parameter_doc.all b/ice_ocean_SIS2/OM_1deg/MOM_parameter_doc.all
index 0055575fd..8aba3a63d 100644
--- a/ice_ocean_SIS2/OM_1deg/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/OM_1deg/MOM_parameter_doc.all
@@ -1385,6 +1385,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/SIS2/MOM_parameter_doc.all b/ice_ocean_SIS2/SIS2/MOM_parameter_doc.all
index 2457784fc..f81687885 100644
--- a/ice_ocean_SIS2/SIS2/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/SIS2/MOM_parameter_doc.all
@@ -1123,6 +1123,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/SIS2_bergs_cgrid/MOM_parameter_doc.all b/ice_ocean_SIS2/SIS2_bergs_cgrid/MOM_parameter_doc.all
index c94878600..ade1d1b99 100644
--- a/ice_ocean_SIS2/SIS2_bergs_cgrid/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/SIS2_bergs_cgrid/MOM_parameter_doc.all
@@ -1123,6 +1123,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ice_ocean_SIS2/SIS2_cgrid/MOM_parameter_doc.all b/ice_ocean_SIS2/SIS2_cgrid/MOM_parameter_doc.all
index 9359c45a1..d14a4ae72 100644
--- a/ice_ocean_SIS2/SIS2_cgrid/MOM_parameter_doc.all
+++ b/ice_ocean_SIS2/SIS2_cgrid/MOM_parameter_doc.all
@@ -1123,6 +1123,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ocean_only/circle_obcs/MOM_parameter_doc.all b/ocean_only/circle_obcs/MOM_parameter_doc.all
index 3b954adee..15352aeb4 100644
--- a/ocean_only/circle_obcs/MOM_parameter_doc.all
+++ b/ocean_only/circle_obcs/MOM_parameter_doc.all
@@ -991,6 +991,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ocean_only/global_ALE/hycom/MOM_parameter_doc.all b/ocean_only/global_ALE/hycom/MOM_parameter_doc.all
index 5e9f79ea0..6fc4640c7 100644
--- a/ocean_only/global_ALE/hycom/MOM_parameter_doc.all
+++ b/ocean_only/global_ALE/hycom/MOM_parameter_doc.all
@@ -1418,6 +1418,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ocean_only/global_ALE/layer/MOM_parameter_doc.all b/ocean_only/global_ALE/layer/MOM_parameter_doc.all
index 5467f21c9..06152829f 100644
--- a/ocean_only/global_ALE/layer/MOM_parameter_doc.all
+++ b/ocean_only/global_ALE/layer/MOM_parameter_doc.all
@@ -1313,6 +1313,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ocean_only/global_ALE/z/MOM_parameter_doc.all b/ocean_only/global_ALE/z/MOM_parameter_doc.all
index d7b121611..285e31a7b 100644
--- a/ocean_only/global_ALE/z/MOM_parameter_doc.all
+++ b/ocean_only/global_ALE/z/MOM_parameter_doc.all
@@ -1365,6 +1365,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/ocean_only/nonBous_global/MOM_parameter_doc.all b/ocean_only/nonBous_global/MOM_parameter_doc.all
index cd7ae40b3..cb8cab5d4 100644
--- a/ocean_only/nonBous_global/MOM_parameter_doc.all
+++ b/ocean_only/nonBous_global/MOM_parameter_doc.all
@@ -1136,6 +1136,13 @@ BETTER_BOUND_AH = True          !   [Boolean] default = True
 RE_AH = 0.0                     !   [nondim] default = 0.0
                                 ! If nonzero, the biharmonic coefficient is scaled so that the biharmonic
                                 ! Reynolds number is equal to this.
+BACKSCATTER_UNDERBOUND = True   !   [Boolean] default = True
+                                ! If true, the bounds on the biharmonic viscosity are allowed to increase where
+                                ! the Laplacian viscosity is negative (due to backscatter parameterizations)
+                                ! beyond the largest timestep-dependent stable values of biharmonic viscosity
+                                ! when no Laplacian viscosity is applied.  The default is true for historical
+                                ! reasons, but this option probably should not be used because it can contribute
+                                ! to numerical instabilities.
 SMAG_BI_CONST = 0.06            !   [nondim] default = 0.0
                                 ! The nondimensional biharmonic Smagorinsky constant, typically 0.015 - 0.06.
 BOUND_CORIOLIS_BIHARM = True    !   [Boolean] default = True
diff --git a/src/MOM6 b/src/MOM6
index 09774d410..36dda7e4e 160000
--- a/src/MOM6
+++ b/src/MOM6
@@ -1 +1 @@
-Subproject commit 09774d410007b98e4ebcf938ed39d138213c5ab7
+Subproject commit 36dda7e4ec1a1ca7b17b4d83feda336803a9b2f1