Fix inconsistency of != with is_zero() and matrix symmetry

src/sage/symbolic/expression.pyx

@@ -3287,7 +3287,6 @@ cdef class Expression(Expression_abc):
         ::
 
             sage: assert(not x == 1)
-            sage: assert(not x != 1)
             sage: forget()
             sage: assume(x>y)
             sage: assert(not x==y)
@@ -3374,6 +3373,21 @@ cdef class Expression(Expression_abc):
             sage: expr = reduce(lambda u, v: 1/u -v, [1/pi] + list(continued_fraction(pi)[:20]))
             sage: expr.is_zero()
             False
+
+        Check that :issue:`41125` is fixed::
+
+            sage: y = SR.var("y")
+            sage: bool(y != 0)
+            True
+            sage: y = SR.var("y", domain="real")
+            sage: bool(y != 0)
+            True
+            sage: z = SR.var("z", domain="complex")
+            sage: bool(z != 0)
+            True
+            sage: z = SR.var("z", domain="integer")
+            sage: bool(z != 0)
+            True
         """
         if self.is_relational():
             # constants are wrappers around Sage objects, compare directly
@@ -3388,17 +3402,8 @@ cdef class Expression(Expression_abc):
                 return pynac_result == relational_true
 
             if pynac_result == relational_true:
-                if self.operator() == operator.ne:
-                    # this hack is necessary to catch the case where the
-                    # operator is != but is False because of assumptions made
-                    m = self._maxima_()
-                    s = m.parent()._eval_line('is (notequal(%s,%s))' % (repr(m.lhs()),repr(m.rhs())))
-                    if s == 'false':
-                        return False
-                    else:
-                        return True
-                else:
-                    return True
+                #In fact, it will return notimplemented for the unequal cases unknown to be true
+                return True
 
             # If assumptions are involved, falsification is more complicated...
             need_assumptions = False
@@ -3434,9 +3439,9 @@ cdef class Expression(Expression_abc):
             # associated with different semantics, different
             # precision, etc., that can lead to subtle bugs.  Also, a
             # lot of basic Sage objects can't be put into maxima.
-            from sage.symbolic.relation import check_relation_maxima
+            from sage.symbolic.relation import check_relation_maxima_neq_as_not_eq
             if self.variables():
-                return check_relation_maxima(self)
+                return check_relation_maxima_neq_as_not_eq(self)
             else:
                 return False
 

src/sage/symbolic/relation.py

@@ -485,32 +485,42 @@
         [k == 1/2*I*sqrt(3) - 1/2, k == -1/2*I*sqrt(3) - 1/2]
         sage: assumptions()
         [k is noninteger]
+
+    Check that boolean values are handled correctly::
+
+        sage: check_relation_maxima(2 == 2)
+        True
+        sage: check_relation_maxima(2 == 3)
+        False
     """
-    m = relation._maxima_()
+    # Handle boolean values directly (e.g., when comparing Python integers)
+    if isinstance(relation, bool):
+        return relation
 
-    # Handle some basic cases first
-    if repr(m) in ['0=0']:
-        return True
-    elif repr(m) in ['0#0', '1#1']:
-        return False
+    from sage.interfaces.maxima_lib import maxima
+
+    # Use _maxima_init_() to get proper string representation with _SAGE_VAR_ prefixes
+    # This ensures assumptions are properly recognized by Maxima
+    lhs_str = relation.lhs()._maxima_init_()
+    rhs_str = relation.rhs()._maxima_init_()
 
     if relation.operator() == operator.eq:  # operator is equality
         try:
-            s = m.parent()._eval_line('is (equal(%s,%s))' % (repr(m.lhs()),
-                                                             repr(m.rhs())))
+            s = maxima._eval_line('is (equal(%s,%s))' % (lhs_str, rhs_str))
         except TypeError:
             raise ValueError("unable to evaluate the predicate '%s'" % repr(relation))
 
     elif relation.operator() == operator.ne:  # operator is not equal
         try:
-            s = m.parent()._eval_line('is (notequal(%s,%s))' % (repr(m.lhs()),
-                                                                repr(m.rhs())))
+            s = maxima._eval_line('is (notequal(%s,%s))' % (lhs_str, rhs_str))
         except TypeError:
             raise ValueError("unable to evaluate the predicate '%s'" % repr(relation))
 
     else:  # operator is < or > or <= or >=, which Maxima handles fine
         try:
-            s = m.parent()._eval_line('is (%s)' % repr(m))
+            # For inequalities, use the full relation string
+            relation_str = relation._maxima_init_()
+            s = maxima._eval_line('is (%s)' % relation_str)
         except TypeError:
             raise ValueError("unable to evaluate the predicate '%s'" % repr(relation))
 
@@ -526,28 +536,49 @@
     if difference.is_trivial_zero():
         return True
 
-    # Try to apply some simplifications to see if left - right == 0.
-    #
-    # TODO: If simplify_log() is ever removed from simplify_full(), we
-    # can replace all of these individual simplifications with a
-    # single call to simplify_full(). That would work in cases where
-    # two simplifications are needed consecutively; the current
-    # approach does not.
-    #
-    simp_list = [difference.simplify_factorial(),
-                 difference.simplify_rational(),
-                 difference.simplify_rectform(),
-                 difference.simplify_trig()]
-    for f in simp_list:
-        try:
-            if f().is_trivial_zero():
-                return True
-                break
-        except Exception:
-            pass
+    # Try simplify_full() to see if left - right == 0.
+    # Note: simplify_full() does not include simplify_log(), which is
+    # unsafe for complex variables, so this is safe to call here.
+    try:
+        if difference.simplify_full().is_trivial_zero():
+            return True
+    except Exception:
+        pass
     return False
 
 
+def check_relation_maxima_neq_as_not_eq(relation):
+    """
+    A variant of :func:`check_relation_maxima` that treats `x != y`
+    as `not (x == y)` for consistency with Python's boolean semantics.
+
+    For inequality relations (!=), this function checks the corresponding
+    equality and returns its logical negation, ensuring that
+    ``bool(x != y) == not bool(x == y)``.
+
+    EXAMPLES::
+
+        sage: from sage.symbolic.relation import check_relation_maxima_neq_as_not_eq
+        sage: x = var('x')
+        sage: check_relation_maxima_neq_as_not_eq(x != x)
+        False
+        sage: check_relation_maxima_neq_as_not_eq(x == x)
+        True
+        sage: check_relation_maxima_neq_as_not_eq(x != 1)
+        True
+        sage: check_relation_maxima_neq_as_not_eq(x == 1)
+        False
+    """
+    # For inequality (!=), check equality and return the opposite.
+    # This ensures bool(x != y) == not bool(x == y) for semantic consistency.
+    if relation.operator() == operator.ne:
+        eq_relation = (relation.lhs() == relation.rhs())
+        return not check_relation_maxima(eq_relation)
+
+    # For all other relations, delegate to check_relation_maxima
+    return check_relation_maxima(relation)
+
+
 def string_to_list_of_solutions(s):
     r"""
     Used internally by the symbolic solve command to convert the output
@@ -1627,7 +1658,7 @@
 
     We can solve with respect to a bigger modulus if it consists only of small prime factors::
 
         sage: [d] = solve_mod([5*x + y == 3, 2*x - 3*y == 9], 3*5*7*11*19*23*29, solution_dict = True)
         sage: d[x]
         12915279
         sage: d[y]
@@ -1674,7 +1705,7 @@
         sage: solve_mod([2*x^2 + x*y, -x*y+2*y^2+x-2*y, -2*x^2+2*x*y-y^2-x-y], 12)
         [(0, 0), (4, 4), (0, 3), (4, 7)]
         sage: eqs = [-y^2+z^2, -x^2+y^2-3*z^2-z-1, -y*z-z^2-x-y+2, -x^2-12*z^2-y+z]
         sage: solve_mod(eqs, 11)
         [(8, 5, 6)]
 
     Confirm that modulus 1 now behaves as it should::
@@ -1793,7 +1824,7 @@
     Confirm we can reproduce the first few terms of :oeis:`A187719`::
 
         sage: from sage.symbolic.relation import _solve_mod_prime_power
         sage: [sorted(_solve_mod_prime_power([x^2==41], 10, i, [x]))[0][0] for i in [1..13]]
         [1, 21, 71, 1179, 2429, 47571, 1296179, 8703821, 26452429, 526452429,
         13241296179, 19473547571, 2263241296179]
     """