diff --git a/sklearn/cluster/optics_.py b/sklearn/cluster/optics_.py
index 55f286404c9ee..c211b86a30eab 100755
--- a/sklearn/cluster/optics_.py
+++ b/sklearn/cluster/optics_.py
@@ -844,7 +844,10 @@ def _xi_cluster(reachability_plot, predecessor_plot, ordering, xi, min_samples,
                     # Find the first index from the right side which is almost
                     # at the same level as the beginning of the detected
                     # cluster.
-                    while (reachability_plot[c_end - 1] < D_max
+                    # Our implementation corrects a mistake in the original
+                    # paper, i.e., in Definition 11 4c, r(x) < r(sD) should be
+                    # r(x) > r(sD).
+                    while (reachability_plot[c_end - 1] > D_max
                            and c_end > U_start):
                         c_end -= 1
 
diff --git a/sklearn/cluster/tests/test_optics.py b/sklearn/cluster/tests/test_optics.py
index 4aa7899f2296c..7dac8e3e40d7d 100644
--- a/sklearn/cluster/tests/test_optics.py
+++ b/sklearn/cluster/tests/test_optics.py
@@ -98,7 +98,7 @@ def test_extract_xi():
     X, expected_labels = shuffle(X, expected_labels, random_state=rng)
 
     clust = OPTICS(min_samples=3, min_cluster_size=2,
-                   max_eps=np.inf, cluster_method='xi',
+                   max_eps=20, cluster_method='xi',
                    xi=0.4).fit(X)
     assert_array_equal(clust.labels_, expected_labels)
 
@@ -110,7 +110,7 @@ def test_extract_xi():
     pytest.skip('FIXME (#13739): This is not stable across platforms.')
 
     clust = OPTICS(min_samples=3, min_cluster_size=3,
-                   max_eps=np.inf, cluster_method='xi',
+                   max_eps=20, cluster_method='xi',
                    xi=0.1).fit(X)
     # this may fail if the predecessor correction is not at work!
     assert_array_equal(clust.labels_, expected_labels)
@@ -129,9 +129,10 @@ def test_extract_xi():
 
 
 def test_cluster_hierarchy_():
+    rng = np.random.RandomState(0)
     n_points_per_cluster = 100
     C1 = [0, 0] + 2 * rng.randn(n_points_per_cluster, 2)
-    C2 = [0, 0] + 10 * rng.randn(n_points_per_cluster, 2)
+    C2 = [0, 0] + 50 * rng.randn(n_points_per_cluster, 2)
     X = np.vstack((C1, C2))
     X = shuffle(X, random_state=0)