Fix subs failure involving integer-valued rational exponents

[EmmanuelCharpentier]

From Ask Sage question 54986:
subs fails with a memory error (or doesn't return) on some expressions.

The following works fine:

sage: _ = var('A, L, G, R, f, k, n, q, u, beta, gamma', domain="positive")
sage: a = I*R^2*f^3*k*q*A*u
sage: b = 2*pi*L*R^2*G*f^4*k^2*q - 2*pi*L*R^2*G*f^4*q - 2*pi*L*R^2*beta^2*G*q
sage: c = (2*I*pi*L*R^2*beta*gamma*q + 2*I*pi*L*R*(beta + q))*G*f^3
sage: d = 2*(pi*(beta^2 + 1)*L*R^2*q + pi*L*R*beta*gamma*q + pi*L*beta)*G*f^2
sage: e = (-2*I*pi*L*R^2*beta*gamma*q - 2*I*pi*(beta^2*q + beta)*L*R)*G*f
sage: expr = a / ((b + c + d + e)*n)
sage: R1 = (expr.real()^2 + expr.imag()^2).factor()
sage: R2 = (sqrt(expr.real()^2 + expr.imag()^2)^2).factor()
sage: R3 = ((sqrt(expr.real()^2 + expr.imag()^2).factor())^2).factor()
sage: all((R1 == R2, R1 == R3, R2 == R3))
True

These expressions have the same string representation:

sage: [len(str(ex)) for ex in (R1, R2, R3)]
[734, 734, 734]
sage: all((str(R1) == str(R2), str(R1) == str(R3), str(R2) == str(R3)))

Substituting f = 2*beta in R1 R2 R3 works:

sage: R1s = R1.subs(f = 2*beta)
sage: R2s = R2.subs(f = 2*beta)
sage: R3s = R3.subs(f = 2*beta)

However, the following never return:

sage: bool(R3s == R1s)  # hangs

sage: len(str(R1s))
520
sage: len(str(R2s))
520
sage: len(str(R3s))  # hangs

Attempting to interrupt the last call can result in a Sage crash;
the original poster reports a "Memory Error" exception.

As a comparison, substituting in R3 via SymPy works:

sage: import sympy
sage: bool(sympy.sympify(R3).subs(f, 2*beta)._sage_() == R1.subs(f = 2*beta))
True

but seems to fail via Mathematica somehow:

sage: mma = mathematica
sage: bool(mma.Replace(R3, mma.Rule(f, 2*beta)).sage() == R1.subs(f = 2*beta))
False

Priority set to critical, because it affects a very basic feature of Sage.

---

RESOLUTION: solved by the patch to pynac in #30446. There was a bug in evaluating an expression (perhaps sqrt(expr.real())^2) that has a rational exponent (such as 2/2) that simplifies to an integer. This ticket just adds a doctest to ensure that the problem is not allowed to return.

Depends on #30786
Depends on #30446

CC: @slel

Component: symbolics

Keywords: factor

Author: Dave Morris

Branch/Commit: c2b45c0

Reviewer: Dima Pasechnik

Issue created by migration from https://trac.sagemath.org/ticket/31137

[slel]

comment:1

Seems to have a lot to do with factor.

I also got a crash when interrupting the computation:

KeyboardInterrupt:
Exception ignored in: 'sage.libs.pynac.pynac.py_repr'
Traceback (most recent call last):
  File "sage/rings/integer.pyx", line 1013, in sage.rings.integer.Integer.__repr__ (build/cythonized/sage/rings/integer.
c:8427)
  File "sage/rings/integer.pyx", line 1134, in sage.rings.integer.Integer.str (build/cythonized/sage/rings/integer.c:901
1)
KeyboardInterrupt:
------------------------------------------------------------------------
...
------------------------------------------------------------------------
Unhandled SIGSEGV: A segmentation fault occurred.
This probably occurred because a *compiled* module has a bug
in it and is not properly wrapped with sig_on(), sig_off().
Python will now terminate.
------------------------------------------------------------------------
Segmentation fault

[slel]

comment:2

Seems expressions resulting from factor behave strangely.
Not sure what role subs plays here.

[slel]

Changed keywords from none to factor

[slel]

comment:3

What seemed to be hanging was just taking a long time.
Still, the result is surprising.

sage: len(str(R3s))  # long time -- 1,292,914,505
1292914505

[DaveWitteMorris]

Branch: public/31137

[DaveWitteMorris]

Dependencies: #30786

[DaveWitteMorris]

Author: Dave Morris

[DaveWitteMorris]

comment:5

The problem in this ticket seems to be solved by the patch to pynac in #30446. The PR just adds a doctest. It is based on #30786 in order to avoid a merge conflict.

---

New commits:

dc073df	fixes for trac #30446
4a5d053	doctests for trac 28620, 30304, 30786
5da750a	doctest for trac 31137

[DaveWitteMorris]

Commit: 5da750a

[dimpase]

Changed dependencies from #30786 to #30786, #30446

[dimpase]

comment:6

only the patch for expression.pyx is needed here, with the rest in #30446

[sagetrac-git]

Changed commit from 5da750a to c2b45c0

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. This was a forced push. Last 10 new commits:

f74f66c	sage --package update-latest: Accept package classes :standard:, :optional: etc., restrict to normal Python packages
182b3d2	sage -package fix-checksum: Handle package classes, ignore non-normal packages
9a57cf6	pynac update
563940e	doctest from #30446
51def3c	sane version name
4834dc8	remove upstreamed patch
214712d	tarball update
6ee3113	doctests for trac 28620, 30304, 30786
6161215	fixes for trac #30446
c2b45c0	doctest for trac 31137

[dimpase]

comment:8

LGTM

[dimpase]

Reviewer: Dima Pasechnik

[slel]

comment:9

Can we improve the ticket summary and ticket description?

[dimpase]

comment:10

it just adds an extra doctest, which used to fail before #30446

[slel]

comment:12

Thanks for the explanation.

[vbraun]

Changed branch from public/31137 to c2b45c0