Opened 2 years ago

Closed 21 months ago

# Fix subs failure involving integer-valued rational exponents

Reported by: Owned by: Emmanuel Charpentier minor sage-9.3 symbolics factor Samuel Lelièvre Dave Morris Dima Pasechnik N/A c2b45c0 c2b45c0b3181ff8803e116e4c92ecf4c4bd748d9 #30786, #30446

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.

### comment:1 Changed 2 years ago by Samuel Lelièvre

Cc: Samuel Lelièvre added modified (diff)

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
```

### comment:2 Changed 2 years ago by Samuel Lelièvre

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

### comment:3 Changed 2 years ago by Samuel Lelièvre

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
```

### comment:4 Changed 23 months ago by Dave Morris

Branch: → public/31137

### comment:5 Changed 23 months ago by Dave Morris

Authors: → Dave Morris → 5da750a7fca49b55bd131585f0c620490611834b → #30786 critical → minor new → needs_review

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`

### comment:6 Changed 23 months ago by Dima Pasechnik

Dependencies: #30786 → #30786, #30446 needs_review → needs_work

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

### comment:7 Changed 23 months ago by git

Commit: 5da750a7fca49b55bd131585f0c620490611834b → c2b45c0b3181ff8803e116e4c92ecf4c4bd748d9

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`

### comment:8 Changed 23 months ago by Dima Pasechnik

Reviewers: → Dima Pasechnik needs_work → positive_review

LGTM

### comment:9 Changed 23 months ago by Samuel Lelièvre

Can we improve the ticket summary and ticket description?

### comment:10 Changed 23 months ago by Dima Pasechnik

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

### comment:11 Changed 23 months ago by Dave Morris

Description: modified (diff)

### comment:12 Changed 23 months ago by Samuel Lelièvre

Description: modified (diff) subs fails in obscure circumstances → Fix subs failure involving integer-valued rational exponents

Thanks for the explanation.

### comment:13 Changed 21 months ago by Volker Braun

Branch: public/31137 → c2b45c0b3181ff8803e116e4c92ecf4c4bd748d9 → fixed positive_review → closed
Note: See TracTickets for help on using tickets.