try an idea

Author: fchapoton

src/sage/combinat/posets/posets.py

@@ -5232,7 +5232,7 @@ def factor(self):
         """
         Factor the poset as a Cartesian product of smaller posets.
 
-        This only works for connected posets for the moment.
+        This only works for connected posets.
 
         The decomposition of a connected poset as a Cartesian product
         of posets (prime in the sense that they cannot be written as
@@ -5309,22 +5309,23 @@ def factor(self):
         is_product, dic = G.is_cartesian_product(relabeling=True)
         if not is_product:
             return [self]
-        dic = {key: tuple(flatten(dic[key])) for key in dic}
+        dic = {key: tuple(flatten(val)) for key, val in dic.items()}
 
         prod_dg = dg.relabel(dic, inplace=False)
         v0 = next(iter(dic.values()))
         n = len(v0)
         factors_range = range(n)
 
         def edge_color(va, vb):
-            for i in range(n):
-                if va[i] != vb[i]:
-                    return i
-
-        neighbors_table = {x: [[y for y in prod_dg.neighbor_iterator(x)
-                                if edge_color(x, y) == i]
-                               for i in factors_range]
-                           for x in prod_dg}
+            return next(i for i, (vai, vbi) in enumerate(zip(va, vb))
+                        if vai != vbi)
+
+        neighbors_table = {}
+        for x in prod_dg:
+            z = [[] for _ in range(n)]
+            for y in prod_dg.neighbor_iterator(x):
+                z[edge_color(x, y)].append(y)
+            neighbors_table[x] = z
 
         fusion_edges = []
         for i0, i1 in Subsets(factors_range, 2):