Opened 7 years ago

Last modified 5 weeks ago

#21071 needs_work defect

substitution in denominator is skipped

Reported by: dkrenn Owned by:
Priority: major Milestone:
Component: symbolics Keywords:
Cc: tscrim, dkrenn Merged in:
Authors: Ralf Stephan Reviewers:
Report Upstream: N/A Work issues:
Branch: u/rws/substitution_in_denominator_is_skipped (Commits, GitHub, GitLab) Commit: 5e39b7b65f7f3ae6e778adacf9514f137c930900
Dependencies: #22102 Stopgaps:

Status badges

Description (last modified by rws)

The following comes very unexpected

sage: ((1+x^2)/x^2).subs({x^2: 42})

The problem is the internal representation as

sage: ((1+x^2)/x^2).operands()
[x^2 + 1, x^(-2)]

Reported upstream as

Change History (18)

comment:1 Changed 6 years ago by rws

Description: modified (diff)
Report Upstream: N/AReported upstream. Developers acknowledge bug.
Summary: substitution in denomintor is skippedsubstitution in denominator is skipped

comment:2 Changed 6 years ago by rws

Note that there is the walkaround of substituting twice, i.e. as given followed by a substitution of the inverse pattern with the inverse replacement, e.g., ((1+x^2)/(x^2)).subs({x^2: 42}).subs({x^-2: 1/42})

comment:3 Changed 6 years ago by rws

Branch: u/rws/substitution_in_denominator_is_skipped

comment:4 Changed 6 years ago by rws

Authors: Ralf Stephan
Cc: tscrim dkrenn added
Commit: 67f46511fb82b1c96182c212cb6112bceda240e0
Milestone: sage-7.3sage-7.4
Report Upstream: Reported upstream. Developers acknowledge bug.N/A

Remaining fail, hopefully the author(s) can see where they worked around this bug so that their workaround can be removed:

File "src/sage/rings/asymptotic/", line 1717, in sage.rings.asymptotic.asymptotics_multivariate_generating_functions.FractionWithFactoredDenominator.?
Failed example:
    asy # long time
    (4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 47/216*sqrt(3)/(sqrt(pi)*sqrt(r)),
     1, 4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 47/216*sqrt(3)/(sqrt(pi)*sqrt(r)))
    (4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 587/216*sqrt(3)/(sqrt(pi)*sqrt(r)),
     4/3*sqrt(3)*sqrt(r)/sqrt(pi) + 587/216*sqrt(3)/(sqrt(pi)*sqrt(r)))
File "src/sage/rings/asymptotic/", line 1720, in sage.rings.asymptotic.asymptotics_multivariate_generating_functions.FractionWithFactoredDenominator.?
Failed example:
    F.relative_error(asy[0], alpha, [1, 2, 4, 8], asy[1]) # long time
    [((3, 3, 2), 0.9812164307, [1.515572606], [-0.54458543...]),
     ((6, 6, 4), 1.576181132, [1.992989399], [-0.26444185...]),
     ((12, 12, 8), 2.485286378, [2.712196351], [-0.091301338...]),
     ((24, 24, 16), 3.700576827, [3.760447895], [-0.016178847...])]
    [((3, 3, 2), 0.9812164307, [3.958585166], [-3.034364939]),
     ((6, 6, 4), 1.576181132, [3.720460147], [-1.360426775]),
     ((12, 12, 8), 2.485286378, [3.933702631], [-0.5827965202]),
     ((24, 24, 16), 3.700576827, [4.624183269], [-0.2495844526])]

New commits:

67f465121071: substitution in denominator is skipped

comment:5 Changed 6 years ago by git

Commit: 67f46511fb82b1c96182c212cb6112bceda240e0deeea0ac1ab28e246944059bbd576c554a6f4579

Branch pushed to git repo; I updated commit sha1. New commits:

deeea0a21071: refine algorithm, fixes fail

comment:6 Changed 6 years ago by rws

Status: newneeds_review

comment:7 Changed 6 years ago by git

Commit: deeea0ac1ab28e246944059bbd576c554a6f45795e39b7b65f7f3ae6e778adacf9514f137c930900

Branch pushed to git repo; I updated commit sha1. New commits:

5e39b7bMerge branch 'develop' into t/21071/substitution_in_denominator_is_skipped

comment:8 Changed 6 years ago by tscrim

I'm not completely sure I agree with this:

     sage: eq.substitute(a=x, x=1)
-    x + 1 == sin(1/x)
+    x + 1 == sin(1)

On the LHS, you simplify x then a (or simultaneously where one substitution does not affect the other), whereas on the RHS, you currently do a, then x (with the new x from the a). So I see this as an inconsistency. IMO, they should be done simultaneously, and the RHS should be sin(1/x).

comment:9 Changed 6 years ago by rws

I agree, but the inconsistency already exists in develop:

sage: sin(x/a).subs({1/x : 1, a : x, x : 1})
sage: sin(x/a).subs({a^-1 : x^-1, x^-1 : 1})

(EDIT: simplified)

Last edited 6 years ago by rws (previous) (diff)

comment:10 Changed 6 years ago by rws

However, the following shows that it is not due to sequential application but premature evaluation ('x/x' --> 1).

sage: sin(x/a).subs({a^-1:x^-1,x^-1:pi})
Version 0, edited 6 years ago by rws (next)

comment:11 Changed 6 years ago by nbruin

Quoting from, it looks like this is pretty fundamental:

  606 ex basic::subs(const exmap & m, unsigned options) const
  607 {
  608     size_t num = nops();
  609     if (num) {
  611         // Substitute in subexpressions
  612         for (size_t i=0; i<num; i++) {
  613             const ex & orig_op = op(i);
  614             const ex & subsed_op = orig_op.subs(m, options);
  615             if (!are_ex_trivially_equal(orig_op, subsed_op)) {
  617                 // Something changed, clone the object
  618                 basic *copy = duplicate();
  619                 copy->clearflag(status_flags::hash_calculated | status_flags::expanded);
  621                 // Substitute the changed operand
  622                 copy->let_op(i++) = subsed_op;
  624                 // Substitute the other operands
  625                 for (; i<num; i++)
  626                     copy->let_op(i) = op(i).subs(m, options);
  628                 // Perform substitutions on the new object as a whole
  629                 return copy->subs_one_level(m, options);
  630             }
  631         }
  632     }
  634     // Nothing changed or no subexpressions
  635     return subs_one_level(m, options);
  636 }

The problem is that substitutions are done on the subexpressions and *then* on the entire expression. So, we have:

sage: (x/a).subs({a:x,x:1})
sage: (x/a).subs({a:x,x:1,1/x:17})

I'm pretty sure this is because of line 629. This is not because of premature evaluation (which, as you show, can also be an issue)

comment:12 Changed 6 years ago by rws

Dependencies: #22102

I would like to solve this outside Pynac via an intermediate substitution step. It would depend on #22102.

comment:13 Changed 6 years ago by rws

There is an ambiguity. What is the expected result of

sage: f = piecewise([((0,2), x), ([-2,0], -x)]
sage: (x+f).subs(x==1)

piecewise(x|-->1 on (0, 2), x|-->-1 on [-2, 0]; x) + 1

Plot expects this to be 2.

Last edited 6 years ago by rws (previous) (diff)

comment:14 Changed 6 years ago by rws

I think subs should do nothing inside piecewise and only evaluate if the main variable is substituted. If one wishes to change the piece definition one should use another method that is still to be implemented, like e.g. substitute_piece. Otherwise I don't see how to untangle this issue.

comment:15 Changed 6 years ago by rws

Status: needs_reviewneeds_work

comment:16 Changed 6 years ago by tscrim

I think what it should do is return the constant function 3 as this would be in line with how subs works on general functions/expressions:

sage: f(x) = x^2
sage: f.subs(x=2)
x |--> 4
sage: f
x |--> x^2
sage: (f+2).subs(x=2)
x |--> 6
sage: (f(x)+2).subs(x=2)

comment:17 Changed 6 years ago by nbruin

Substitution in computer algebra systems doesn't tend to be semantically defined. It's an operation on syntax trees. So I think piecewise functions can be treated in two ways:

Either you regard them as "atomic" objects. Then substitution doesn't do anything within them. That's not unreasonable. A piecewise function cannot be evaluated at a symbolic input either ...

Or you regard them as a syntax tree (probably a list of (a,b,expression) tuples or so) and you do a substitution on that (hopefully raising an error when the substitution leads to something illegal).

comment:18 Changed 5 weeks ago by mkoeppe

Milestone: sage-7.4
Note: See TracTickets for help on using tickets.