diff --git a/pyomo/contrib/gdpopt/tests/test_LBB.py b/pyomo/contrib/gdpopt/tests/test_LBB.py
index 273327b02a4..8a553398fa6 100644
--- a/pyomo/contrib/gdpopt/tests/test_LBB.py
+++ b/pyomo/contrib/gdpopt/tests/test_LBB.py
@@ -59,6 +59,7 @@ def test_infeasible_GDP(self):
         self.assertIsNone(m.d.disjuncts[0].indicator_var.value)
         self.assertIsNone(m.d.disjuncts[1].indicator_var.value)
 
+    @unittest.skipUnless(z3_available, "Z3 SAT solver is not available")
     def test_infeasible_GDP_check_sat(self):
         """Test for infeasible GDP with check_sat option True."""
         m = ConcreteModel()
diff --git a/pyomo/contrib/gdpopt/tests/test_gdpopt.py b/pyomo/contrib/gdpopt/tests/test_gdpopt.py
index 005df56ced5..873bafabc76 100644
--- a/pyomo/contrib/gdpopt/tests/test_gdpopt.py
+++ b/pyomo/contrib/gdpopt/tests/test_gdpopt.py
@@ -22,7 +22,6 @@
 from pyomo.common.collections import Bunch
 from pyomo.common.config import ConfigDict, ConfigValue
 from pyomo.common.fileutils import import_file, PYOMO_ROOT_DIR
-from pyomo.contrib.appsi.solvers.gurobi import Gurobi
 from pyomo.contrib.gdpopt.create_oa_subproblems import (
     add_util_block,
     add_disjunct_list,
@@ -767,6 +766,9 @@ def test_time_limit(self):
             results.solver.termination_condition, TerminationCondition.maxTimeLimit
         )
 
+    @unittest.skipUnless(
+        license_available, "No BARON license--8PP logical problem exceeds demo size"
+    )
     def test_LOA_8PP_logical_default_init(self):
         """Test logic-based outer approximation with 8PP."""
         exfile = import_file(join(exdir, 'eight_process', 'eight_proc_logical.py'))
@@ -870,6 +872,9 @@ def test_LOA_8PP_maxBinary(self):
         )
         ct.check_8PP_solution(self, eight_process, results)
 
+    @unittest.skipUnless(
+        license_available, "No BARON license--8PP logical problem exceeds demo size"
+    )
     def test_LOA_8PP_logical_maxBinary(self):
         """Test logic-based OA with max_binary initialization."""
         exfile = import_file(join(exdir, 'eight_process', 'eight_proc_logical.py'))
@@ -1050,7 +1055,11 @@ def assert_correct_disjuncts_active(
 
         self.assertTrue(fabs(value(eight_process.profit.expr) - 68) <= 1e-2)
 
-    @unittest.skipUnless(Gurobi().available(), "APPSI Gurobi solver is not available")
+    @unittest.skipUnless(
+        SolverFactory('appsi_gurobi').available(exception_flag=False)
+        and SolverFactory('appsi_gurobi').license_is_valid(),
+        "Legacy APPSI Gurobi solver is not available",
+    )
     def test_auto_persistent_solver(self):
         exfile = import_file(join(exdir, 'eight_process', 'eight_proc_model.py'))
         m = exfile.build_eight_process_flowsheet()
@@ -1126,6 +1135,9 @@ def test_RIC_8PP_default_init(self):
         )
         ct.check_8PP_solution(self, eight_process, results)
 
+    @unittest.skipUnless(
+        license_available, "No BARON license--8PP logical problem exceeds demo size"
+    )
     def test_RIC_8PP_logical_default_init(self):
         """Test logic-based outer approximation with 8PP."""
         exfile = import_file(join(exdir, 'eight_process', 'eight_proc_logical.py'))
diff --git a/pyomo/contrib/gdpopt/util.py b/pyomo/contrib/gdpopt/util.py
index 2cb70f0ea60..babe0245d57 100644
--- a/pyomo/contrib/gdpopt/util.py
+++ b/pyomo/contrib/gdpopt/util.py
@@ -553,6 +553,13 @@ def _add_bigm_constraint_to_transformed_model(m, constraint, block):
     # making a Reference to the ComponentData so that it will look like an
     # indexed component for now. If I redesign bigm at some point, then this
     # could be prettier.
-    bigm._transform_constraint(Reference(constraint), parent_disjunct, None, [], [])
+    bigm._transform_constraint(
+        Reference(constraint),
+        parent_disjunct,
+        None,
+        [],
+        [],
+        1 - parent_disjunct.binary_indicator_var,
+    )
     # Now get rid of it because this is a class attribute!
     del bigm._config
diff --git a/pyomo/gdp/plugins/bigm.py b/pyomo/gdp/plugins/bigm.py
index 3f450dbbd4f..d715d913db8 100644
--- a/pyomo/gdp/plugins/bigm.py
+++ b/pyomo/gdp/plugins/bigm.py
@@ -213,21 +213,15 @@ def _apply_to_impl(self, instance, **kwds):
         bigM = self._config.bigM
         for t in preprocessed_targets:
             if t.ctype is Disjunction:
-                self._transform_disjunctionData(
-                    t,
-                    t.index(),
-                    bigM,
-                    parent_disjunct=gdp_tree.parent(t),
-                    root_disjunct=gdp_tree.root_disjunct(t),
-                )
+                self._transform_disjunctionData(t, t.index(), bigM, gdp_tree)
 
         # issue warnings about anything that was in the bigM args dict that we
         # didn't use
         _warn_for_unused_bigM_args(bigM, self.used_args, logger)
 
-    def _transform_disjunctionData(
-        self, obj, index, bigM, parent_disjunct=None, root_disjunct=None
-    ):
+    def _transform_disjunctionData(self, obj, index, bigM, gdp_tree):
+        parent_disjunct = gdp_tree.parent(obj)
+        root_disjunct = gdp_tree.root_disjunct(obj)
         (transBlock, xorConstraint) = self._setup_transform_disjunctionData(
             obj, root_disjunct
         )
@@ -236,7 +230,7 @@ def _transform_disjunctionData(
         or_expr = 0
         for disjunct in obj.disjuncts:
             or_expr += disjunct.binary_indicator_var
-            self._transform_disjunct(disjunct, bigM, transBlock)
+            self._transform_disjunct(disjunct, bigM, transBlock, gdp_tree)
 
         if obj.xor:
             xorConstraint[index] = or_expr == 1
@@ -249,7 +243,7 @@ def _transform_disjunctionData(
         # and deactivate for the writers
         obj.deactivate()
 
-    def _transform_disjunct(self, obj, bigM, transBlock):
+    def _transform_disjunct(self, obj, bigM, transBlock, gdp_tree):
         # We're not using the preprocessed list here, so this could be
         # inactive. We've already done the error checking in preprocessing, so
         # we just skip it here.
@@ -261,6 +255,12 @@ def _transform_disjunct(self, obj, bigM, transBlock):
 
         relaxationBlock = self._get_disjunct_transformation_block(obj, transBlock)
 
+        indicator_expression = 0
+        node = obj
+        while node is not None:
+            indicator_expression += 1 - node.binary_indicator_var
+            node = gdp_tree.parent_disjunct(node)
+
         # This is crazy, but if the disjunction has been previously
         # relaxed, the disjunct *could* be deactivated.  This is a big
         # deal for Hull, as it uses the component_objects /
@@ -270,13 +270,21 @@ def _transform_disjunct(self, obj, bigM, transBlock):
         # comparing the two relaxations.
         #
         # Transform each component within this disjunct
-        self._transform_block_components(obj, obj, bigM, arg_list, suffix_list)
+        self._transform_block_components(
+            obj, obj, bigM, arg_list, suffix_list, indicator_expression
+        )
 
         # deactivate disjunct to keep the writers happy
         obj._deactivate_without_fixing_indicator()
 
     def _transform_constraint(
-        self, obj, disjunct, bigMargs, arg_list, disjunct_suffix_list
+        self,
+        obj,
+        disjunct,
+        bigMargs,
+        arg_list,
+        disjunct_suffix_list,
+        indicator_expression,
     ):
         # add constraint to the transformation block, we'll transform it there.
         transBlock = disjunct._transformation_block()
@@ -348,7 +356,13 @@ def _transform_constraint(
             bigm_src[c] = (lower, upper)
 
             self._add_constraint_expressions(
-                c, i, M, disjunct.binary_indicator_var, newConstraint, constraint_map
+                c,
+                i,
+                M,
+                disjunct.binary_indicator_var,
+                newConstraint,
+                constraint_map,
+                indicator_expression=indicator_expression,
             )
 
             # deactivate because we relaxed
diff --git a/pyomo/gdp/plugins/bigm_mixin.py b/pyomo/gdp/plugins/bigm_mixin.py
index 510b36b5102..1c3fcb2c64a 100644
--- a/pyomo/gdp/plugins/bigm_mixin.py
+++ b/pyomo/gdp/plugins/bigm_mixin.py
@@ -232,7 +232,14 @@ def _estimate_M(self, expr, constraint):
         return tuple(M)
 
     def _add_constraint_expressions(
-        self, c, i, M, indicator_var, newConstraint, constraint_map
+        self,
+        c,
+        i,
+        M,
+        indicator_var,
+        newConstraint,
+        constraint_map,
+        indicator_expression=None,
     ):
         # Since we are both combining components from multiple blocks and using
         # local names, we need to make sure that the first index for
@@ -244,6 +251,8 @@ def _add_constraint_expressions(
         # over the constraint indices, but I don't think it matters a lot.)
         unique = len(newConstraint)
         name = c.local_name + "_%s" % unique
+        if indicator_expression is None:
+            indicator_expression = 1 - indicator_var
 
         if c.lower is not None:
             if M[0] is None:
@@ -251,7 +260,7 @@ def _add_constraint_expressions(
                     "Cannot relax disjunctive constraint '%s' "
                     "because M is not defined." % name
                 )
-            M_expr = M[0] * (1 - indicator_var)
+            M_expr = M[0] * indicator_expression
             newConstraint.add((name, i, 'lb'), c.lower <= c.body - M_expr)
             constraint_map.transformed_constraints[c].append(
                 newConstraint[name, i, 'lb']
@@ -263,7 +272,7 @@ def _add_constraint_expressions(
                     "Cannot relax disjunctive constraint '%s' "
                     "because M is not defined." % name
                 )
-            M_expr = M[1] * (1 - indicator_var)
+            M_expr = M[1] * indicator_expression
             newConstraint.add((name, i, 'ub'), c.body - M_expr <= c.upper)
             constraint_map.transformed_constraints[c].append(
                 newConstraint[name, i, 'ub']
diff --git a/pyomo/gdp/tests/test_bigm.py b/pyomo/gdp/tests/test_bigm.py
index c6ac49f6d36..3174a95292e 100644
--- a/pyomo/gdp/tests/test_bigm.py
+++ b/pyomo/gdp/tests/test_bigm.py
@@ -19,8 +19,11 @@
     Set,
     Constraint,
     ComponentMap,
+    LogicalConstraint,
+    Objective,
     SolverFactory,
     Suffix,
+    TerminationCondition,
     ConcreteModel,
     Var,
     Any,
@@ -1879,12 +1882,11 @@ def test_m_value_mappings(self):
     # many of the transformed constraints look like this, so can call this
     # function to test them.
     def check_bigM_constraint(self, cons, variable, M, indicator_var):
-        repn = generate_standard_repn(cons.body)
-        self.assertTrue(repn.is_linear())
-        self.assertEqual(repn.constant, -M)
-        self.assertEqual(len(repn.linear_vars), 2)
-        ct.check_linear_coef(self, repn, variable, 1)
-        ct.check_linear_coef(self, repn, indicator_var, M)
+        assertExpressionsEqual(
+            self,
+            cons.body,
+            variable - float(M) * (1 - indicator_var.get_associated_binary()),
+        )
 
     def check_inner_xor_constraint(self, inner_disjunction, outer_disjunct, bigm):
         inner_xor = inner_disjunction.algebraic_constraint
@@ -1949,6 +1951,10 @@ def test_transformed_constraints(self):
                                                     .binary_indicator_var,
                                                 )
                                             ),
+                                            1,
+                                            EXPR.MonomialTermExpression(
+                                                (-1, m.disjunct[1].binary_indicator_var)
+                                            ),
                                         ]
                                     ),
                                 )
@@ -1958,37 +1964,69 @@ def test_transformed_constraints(self):
                 ]
             ),
         )
-        self.assertIsNone(cons1ub.lower)
-        self.assertEqual(cons1ub.upper, 0)
-        self.check_bigM_constraint(
-            cons1ub, m.z, 10, m.disjunct[1].innerdisjunct[0].indicator_var
+        assertExpressionsEqual(
+            self,
+            cons1ub.expr,
+            m.z
+            - 10.0
+            * (
+                1
+                - m.disjunct[1].innerdisjunct[0].binary_indicator_var
+                + 1
+                - m.disjunct[1].binary_indicator_var
+            )
+            <= 0.0,
         )
 
         cons2 = bigm.get_transformed_constraints(m.disjunct[1].innerdisjunct[1].c)
         self.assertEqual(len(cons2), 1)
         cons2lb = cons2[0]
-        self.assertEqual(cons2lb.lower, 5)
-        self.assertIsNone(cons2lb.upper)
-        self.check_bigM_constraint(
-            cons2lb, m.z, -5, m.disjunct[1].innerdisjunct[1].indicator_var
+        assertExpressionsEqual(
+            self,
+            cons2lb.expr,
+            5.0
+            <= m.z
+            - (-5.0)
+            * (
+                1
+                - m.disjunct[1].innerdisjunct[1].binary_indicator_var
+                + 1
+                - m.disjunct[1].binary_indicator_var
+            ),
         )
 
         cons3 = bigm.get_transformed_constraints(m.simpledisjunct.innerdisjunct0.c)
         self.assertEqual(len(cons3), 1)
         cons3ub = cons3[0]
-        self.assertEqual(cons3ub.upper, 2)
-        self.assertIsNone(cons3ub.lower)
-        self.check_bigM_constraint(
-            cons3ub, m.x, 7, m.simpledisjunct.innerdisjunct0.indicator_var
+        assertExpressionsEqual(
+            self,
+            cons3ub.expr,
+            m.x
+            - 7.0
+            * (
+                1
+                - m.simpledisjunct.innerdisjunct0.binary_indicator_var
+                + 1
+                - m.simpledisjunct.binary_indicator_var
+            )
+            <= 2.0,
         )
 
         cons4 = bigm.get_transformed_constraints(m.simpledisjunct.innerdisjunct1.c)
         self.assertEqual(len(cons4), 1)
         cons4lb = cons4[0]
-        self.assertEqual(cons4lb.lower, 4)
-        self.assertIsNone(cons4lb.upper)
-        self.check_bigM_constraint(
-            cons4lb, m.x, -13, m.simpledisjunct.innerdisjunct1.indicator_var
+        assertExpressionsEqual(
+            self,
+            cons4lb.expr,
+            m.x
+            - (-13.0)
+            * (
+                1
+                - m.simpledisjunct.innerdisjunct1.binary_indicator_var
+                + 1
+                - m.simpledisjunct.binary_indicator_var
+            )
+            >= 4.0,
         )
 
         # Here we check that the xor constraint from
@@ -2088,35 +2126,6 @@ def innerIndexed(d, i):
                 m._pyomo_gdp_bigm_reformulation.relaxedDisjuncts,
             )
 
-    def check_first_disjunct_constraint(self, disj1c, x, ind_var):
-        self.assertEqual(len(disj1c), 1)
-        cons = disj1c[0]
-        self.assertIsNone(cons.lower)
-        self.assertEqual(cons.upper, 1)
-        repn = generate_standard_repn(cons.body)
-        self.assertTrue(repn.is_quadratic())
-        self.assertEqual(len(repn.linear_vars), 1)
-        self.assertEqual(len(repn.quadratic_vars), 4)
-        ct.check_linear_coef(self, repn, ind_var, 143)
-        self.assertEqual(repn.constant, -143)
-        for i in range(1, 5):
-            ct.check_squared_term_coef(self, repn, x[i], 1)
-
-    def check_second_disjunct_constraint(self, disj2c, x, ind_var):
-        self.assertEqual(len(disj2c), 1)
-        cons = disj2c[0]
-        self.assertIsNone(cons.lower)
-        self.assertEqual(cons.upper, 1)
-        repn = generate_standard_repn(cons.body)
-        self.assertTrue(repn.is_quadratic())
-        self.assertEqual(len(repn.linear_vars), 5)
-        self.assertEqual(len(repn.quadratic_vars), 4)
-        self.assertEqual(repn.constant, -63)  # M = 99, so this is 36 - 99
-        ct.check_linear_coef(self, repn, ind_var, 99)
-        for i in range(1, 5):
-            ct.check_squared_term_coef(self, repn, x[i], 1)
-            ct.check_linear_coef(self, repn, x[i], -6)
-
     def simplify_cons(self, cons, leq):
         visitor = LinearRepnVisitor({}, {}, {}, None)
         repn = visitor.walk_expression(cons.body)
@@ -2142,30 +2151,76 @@ def check_hierarchical_nested_model(self, m, bigm):
 
         # outer disjunction constraints
         disj1c = bigm.get_transformed_constraints(m.disj1.c)
-        self.check_first_disjunct_constraint(disj1c, m.x, m.disj1.binary_indicator_var)
+        self.assertEqual(len(disj1c), 1)
+        cons = disj1c[0]
+        assertExpressionsEqual(
+            self,
+            cons.expr,
+            m.x[1] ** 2
+            + m.x[2] ** 2
+            + m.x[3] ** 2
+            + m.x[4] ** 2
+            - 143.0 * (1 - m.disj1.binary_indicator_var)
+            <= 1.0,
+        )
 
         disj2c = bigm.get_transformed_constraints(m.disjunct_block.disj2.c)
-        self.check_second_disjunct_constraint(
-            disj2c, m.x, m.disjunct_block.disj2.binary_indicator_var
+        self.assertEqual(len(disj2c), 1)
+        cons = disj2c[0]
+        assertExpressionsEqual(
+            self,
+            cons.expr,
+            (3 - m.x[1]) ** 2
+            + (3 - m.x[2]) ** 2
+            + (3 - m.x[3]) ** 2
+            + (3 - m.x[4]) ** 2
+            - 99.0 * (1 - m.disjunct_block.disj2.binary_indicator_var)
+            <= 1.0,
         )
 
         # inner disjunction constraints
         innerd1c = bigm.get_transformed_constraints(
             m.disjunct_block.disj2.disjunction_disjuncts[0].constraint[1]
         )
-        self.check_first_disjunct_constraint(
-            innerd1c,
-            m.x,
-            m.disjunct_block.disj2.disjunction_disjuncts[0].binary_indicator_var,
+        self.assertEqual(len(innerd1c), 1)
+        cons = innerd1c[0]
+        assertExpressionsEqual(
+            self,
+            cons.expr,
+            m.x[1] ** 2
+            + m.x[2] ** 2
+            + m.x[3] ** 2
+            + m.x[4] ** 2
+            - 143.0
+            * (
+                1
+                - m.disjunct_block.disj2.disjunction_disjuncts[0].binary_indicator_var
+                + 1
+                - m.disjunct_block.disj2.binary_indicator_var
+            )
+            <= 1.0,
         )
 
         innerd2c = bigm.get_transformed_constraints(
             m.disjunct_block.disj2.disjunction_disjuncts[1].constraint[1]
         )
-        self.check_second_disjunct_constraint(
-            innerd2c,
-            m.x,
-            m.disjunct_block.disj2.disjunction_disjuncts[1].binary_indicator_var,
+        self.assertEqual(len(innerd2c), 1)
+        cons = innerd2c[0]
+        assertExpressionsEqual(
+            self,
+            cons.expr,
+            (3 - m.x[1]) ** 2
+            + (3 - m.x[2]) ** 2
+            + (3 - m.x[3]) ** 2
+            + (3 - m.x[4]) ** 2
+            - 99.0
+            * (
+                1
+                - m.disjunct_block.disj2.disjunction_disjuncts[1].binary_indicator_var
+                + 1
+                - m.disjunct_block.disj2.binary_indicator_var
+            )
+            <= 1.0,
         )
 
     def test_hierarchical_badly_ordered_targets(self):
@@ -2193,6 +2248,46 @@ def test_decl_order_opposite_instantiation_order(self):
     def test_do_not_assume_nested_indicators_local(self):
         ct.check_do_not_assume_nested_indicators_local(self, 'gdp.bigm')
 
+    @unittest.skipUnless(gurobi_available, "Gurobi is not available")
+    def test_constraints_not_enforced_when_an_ancestor_indicator_is_False(self):
+        m = ConcreteModel()
+        m.x = Var(bounds=(0, 30))
+
+        m.left = Disjunct()
+        m.left.left = Disjunct()
+        m.left.left.c = Constraint(expr=m.x >= 10)
+        m.left.right = Disjunct()
+        m.left.right.c = Constraint(expr=m.x >= 9)
+        m.left.disjunction = Disjunction(expr=[m.left.left, m.left.right])
+        m.right = Disjunct()
+        m.right.left = Disjunct()
+        m.right.left.c = Constraint(expr=m.x >= 11)
+        m.right.right = Disjunct()
+        m.right.right.c = Constraint(expr=m.x >= 8)
+        m.right.disjunction = Disjunction(expr=[m.right.left, m.right.right])
+        m.disjunction = Disjunction(expr=[m.left, m.right])
+
+        m.equiv_left = LogicalConstraint(
+            expr=m.left.left.indicator_var.equivalent_to(m.right.left.indicator_var)
+        )
+        m.equiv_right = LogicalConstraint(
+            expr=m.left.right.indicator_var.equivalent_to(m.right.right.indicator_var)
+        )
+
+        m.obj = Objective(expr=m.x)
+
+        TransformationFactory('gdp.bigm').apply_to(m)
+        results = SolverFactory('gurobi').solve(m)
+        self.assertEqual(
+            results.solver.termination_condition, TerminationCondition.optimal
+        )
+        self.assertTrue(value(m.right.indicator_var))
+        self.assertFalse(value(m.left.indicator_var))
+        self.assertTrue(value(m.right.right.indicator_var))
+        self.assertFalse(value(m.right.left.indicator_var))
+        self.assertTrue(value(m.left.right.indicator_var))
+        self.assertAlmostEqual(value(m.x), 8)
+
 
 class IndexedDisjunction(unittest.TestCase):
     # this tests that if the targets are a subset of the
diff --git a/pyomo/gdp/util.py b/pyomo/gdp/util.py
index fe11975954d..a8c6393f0b3 100644
--- a/pyomo/gdp/util.py
+++ b/pyomo/gdp/util.py
@@ -144,13 +144,13 @@ def parent(self, u):
         Arg:
             u : A node in the tree
         """
+        if u in self._parent:
+            return self._parent[u]
         if u not in self._vertices:
             raise ValueError(
                 "'%s' is not a vertex in the GDP tree. Cannot "
                 "retrieve its parent." % u
             )
-        if u in self._parent:
-            return self._parent[u]
         else:
             return None