more unified (and usable) solve output

[Chris Smith]

The output from solve for algebraic equations is not usable without testing output because of the variable types of output. 

CASE 1: If no variables are given we get the following from the input:

* univariate equation -> list of values

>>> solve(x**2 - 4)

[-2, 2]

* multivariate equation -> list of one or more dictionaries

>>> solve(x**2 - y)

[{y: x**2}]

CASE 2: If variable(s) are given we get

* list of solutions

>>> solve(x**2-y,x)

[-sqrt(y), sqrt(y)]

* tuples of solutions

>>> solve(x**2-y,x,y)
[(-sqrt(y),y), (sqrt(y),y)]  --> [{x: -sqrt(y)}, {x: sqrt(y)}]

>>> solve([x**2-y],x)

[(-sqrt(y),), (sqrt(y),)]  --> [{x: -sqrt(y)}, {x: sqrt(y)}]

* dictionary for single solution

>>> solve([x-y],x)
{x: y}  --> [{x: y}]

Although the dict=True or set=True will give a standard output, can we at least unify the case for when more than one variable or more than one expression is given so we always get a list of one or more dictionaries?

This would make the output from solve always be a list. It would be a list of values for a univariate expression given without the variable specified or for a expression given with a single variable of interest, the first examples for each of the CASEs above.  (And, again, those could be given as list(s) of dictionaries with the `dict=True` flag.)

/c

[Aaron Meurer]

Changing the output type could break code that solves a specific
equation. I am doubtful whether any users actually understand the
output type behavior of solve without the dict=True flag. So
personally I think we should clean it up. We already recommend using
dict=True to get consistent output types, and this would only affect
users who aren't doing that.

Aaron Meurer

On Thu, Mar 3, 2022 at 8:28 AM Chris Smith <smi...@gmail.com> wrote:
>
> The output from solve for algebraic equations is not usable without testing output because of the variable types of output.
>
> CASE 1: If no variables are given we get the following from the input:
> * univariate equation -> list of values
> >>> solve(x**2 - 4)
> [-2, 2]
> * multivariate equation -> list of one or more dictionaries
> >>> solve(x**2 - y)
> [{y: x**2}]
>
> CASE 2: If variable(s) are given we get

> * tuples of solutions

> >>> solve([x**2-y],x)
> [(-sqrt(y),), (sqrt(y),)]

> * dictionary for single solution
> >>> solve([x-y],x)
> {x: y}
>

> Although the dict=True or set=True will give a standard output, can we at least unify the case for when variables are given so we always get a list of one or more dictionaries? So the above would be `[{x: -sqrt(y)}, {x: sqrt(y)}]` and `[{x: y}]`, respectively. This would then make `solve` always give a list of a) values for a univariate expression, b) a list of one or more dictionaries for every other case. (Case (a) will give a list of dictionaries if `dict=True`.)
>
> /c
>
> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/c391d370-5d64-4b94-ae15-a56901bbb219n%40googlegroups.com.

[Oscar Benjamin]

On Thu, 3 Mar 2022 at 20:28, Aaron Meurer <asme...@gmail.com> wrote:
> On Thu, Mar 3, 2022 at 8:28 AM Chris Smith <smi...@gmail.com> wrote:
> >
> > Although the dict=True or set=True will give a standard output, can we at least unify the case for when variables are given so we always get a list of one or more dictionaries? So the above would be `[{x: -sqrt(y)}, {x: sqrt(y)}]` and `[{x: y}]`, respectively. This would then make `solve` always give a list of a) values for a univariate expression, b) a list of one or more dictionaries for every other case. (Case (a) will give a list of dictionaries if `dict=True`.)
>

> Changing the output type could break code that solves a specific
> equation. I am doubtful whether any users actually understand the
> output type behavior of solve without the dict=True flag. So
> personally I think we should clean it up. We already recommend using
> dict=True to get consistent output types, and this would only affect
> users who aren't doing that.

It would be much better if solve always worked like that but obviously
that's not a backwards compatible change and I'm not sure what would
break. I don't think that the type of the output should depend on the
number of solutions although it could depend on the type of the
arguments (e.g a single equation vs a list of equations). At the
moment the type of the output can depend on the equations themselves
which is wildly awkward:

In [3]: solve([x-1], [x])
Out[3]: {x: 1}

In [4]: solve([x**2-1], [x])
Out[4]: [(-1,), (1,)]

The ideal output for solve is absolutely a list of dicts though and it
would be good to make that the default at least when given a list of
equations. I guarantee that a lot of notebooks or other little bits of
code will depend on each of these random special cases though.

--
Oscar

[Jason Moore]

I think changing this will break tons of code in the wild. Isnt it best make a new "solve_new" and then leave solve be (maybe with a deprecation warning. You could call it `solve_equations` or something.

Jason

moorepants.info
+01 530-601-9791

--
You received this message because you are subscribed to the Google Groups "sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHVvXxROqs0CpGmSX2%3Dk8xJgTos9FPpP1x0GYP0aHCdiYErDnA%40mail.gmail.com.

[gu...@uwosh.edu]

I think the name should be something like "solve2" and provide a couple year deprecation warning after which "solve" points to "solve2".

Jonathan

[Chris Smith]

We could just make a wrapper named `solved = lambda f,*s,**k: solve(f, *s, **k, dict=True)`

/c

[Oscar Benjamin]

On Fri, 4 Mar 2022 at 14:57, Chris Smith <smi...@gmail.com> wrote:
>
> We could just make a wrapper named `solved = lambda f,*s,**k: solve(f, *s, **k, dict=True)`

If there is to be a new API then there are a *lot* of other things
that I would want to change about it compared to solve so it
definitely wouldn't be just a thin wrapper around solve. It would be
nice if nonlinsolve could be the new API but it doesn't have good
return types either.

--
Oscar

[Aaron Meurer]

It would really be a shame to give up the name solve(). Maybe we can
put this on a list of big changes that we can do in 2.0? Although it
would be good to have a working replacement ready before we do any
renaming.

It might be helpful to enumerate some of those things. Here are a few
I've noticed

- In general, the output of solve() should not depend on the specific
mathematical type of the input. This means for instance, linear and
nonlinear equations should be handled in the exact same way. Right
now, even dict=True can't save you:

>>> solve([x**2 - 2, x > 1], x)
Eq(x, sqrt(2))
>>> solve([x**2 - 2, x > 1], x, dict=True)
Eq(x, sqrt(2))

One can imagine wanting to use x > 1 just to limit the solutions
(which can't be done just with the old assumptions), but this has
caused solve() to treat it as a completely different type of equation.

- One issue I see is that it's not always clear what solve()
will/should do in some cases, especially undetermined cases. It seems
like "solving" and "variable elimination" are too often mixed
together. It would be better if they were completely separated, since
it's not clear to me why you would want one when solve can't produce
the other. We should very carefully define what it means to "solve" an
equation (note this is not so trivial once you include systems,
multiple variables, infinite solutions, and inequalities into the
mix).

- Similarly, solve() should "refuse the temptation to guess". In an
API like diff() or integrate(), you are allowed to pass a single
expression without a variable, but only in the case where the
expression has a single variable. In all other cases, it raises an
exception and forces the user to provide a variable. But consider
solve

>>> solve(x*y - y)
[{x: 1}, {y: 0}]

Here solve(x*y - y) could mean one of three things: "solve for x",
"solve for y", or "solve for x and y simultaneously". solve() has
decided to do the latter. But what's worse, consider a similar
example:

>>> solve(sin(x)*x - y)
[{y: x*sin(x)}]

Now solving for x and y simultaneously is impossible, as is solving
for x, so solve() has just picked the one that it knows how to do and
solved for y.

Both of the above should be an error. The user should be required to
specify which variable(s) they want to solve for in ambiguous cases
like this.

- It would be nice to have an API that can work uniformly even in
cases where there are infinitely many solutions (both parametric and
inequality). The traditional solution to this has been solveset(),
which definitely fixes the "uniform API" problem, but the issue is
that the outputs of it aren't exactly easy to manipulate
programmatically. There's also no guarantees at all about the exact
format of a solveset() output. It might be a ConditionSet or an
ImageSet or some Union.

One idea is to make things programmatic with some sort of Solution
object that provides different things like solution parameters in
consistent attributes. This would also allow moving some of the more
complex keyword-arguments from solve() itself into methods of the
Solution object.

- The API should also work uniformly for systems as well as single
equations. Perhaps the cleanest way would be for solve() to always
treat a single equation as a system of one equation. That way the
return type is always uniform (the downside is the common case might
involve an extra layer of indirection, which could be annoying).

There's also this old issue "solve is a mess" which discusses a lot of
these things. https://github.com/sympy/sympy/issues/6659

Aaron Meurer

>
> --
> Oscar
>
> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHVvXxRY4xwA-ST6BDjFLfiOwWykj2eLW_iFhzsfEbmAK%3DSTVQ%40mail.gmail.com.

[Chris Smith]

Aaron mentioned:

>>> solve(sin(x)*x - y)
[{y: x*sin(x)}]

Now solving for x and y simultaneously is impossible, as is solving
for x, so solve() has just picked the one that it knows how to do and
solved for y.

-------------------------------

I find this useful, especially for a system of equations. I just want some sort of closed form solution. If you force the user to specify (and there aren't valid solution paths for some subsets of variables) then the user is forced to iteratte through potential subsets of variables. The issue [here](https://github.com/sympy/sympy/issues/5849) talks about such a system that has caused issues in this regard.

/c

[Chris Smith]

One sort of "Solution" structure that could be used is something like this:

Solution({a: [-1, 1], x: [2 + y], y: [4, 5]})

This structure could be iterated with backsubstitution done in JIT manner: [{a:-1, x:6, y:4}, ...

/c

[Oscar Benjamin]

On Fri, 4 Mar 2022 at 22:02, Aaron Meurer <asme...@gmail.com> wrote:
>
> It would really be a shame to give up the name solve(). Maybe we can
> put this on a list of big changes that we can do in 2.0? Although it
> would be good to have a working replacement ready before we do any
> renaming.

Let's just work on a replacement first and then we can consider what
to do with legacy solve later. Remember there was a previous attempt
to do exactly this and resulted in linsolve, nonlinsolve and solveset.
Unfortunately none of those can replace solve and there are failings
in their APIs that mean we still just want to make a new solve
function that isn't compatible with any of the existing ones.

> On Fri, Mar 4, 2022 at 8:27 AM Oscar Benjamin
> <oscar.j....@gmail.com> wrote:
> >
> > On Fri, 4 Mar 2022 at 14:57, Chris Smith <smi...@gmail.com> wrote:
> > >
> > > We could just make a wrapper named `solved = lambda f,*s,**k: solve(f, *s, **k, dict=True)`
> >
> > If there is to be a new API then there are a *lot* of other things
> > that I would want to change about it compared to solve so it
> > definitely wouldn't be just a thin wrapper around solve. It would be
> > nice if nonlinsolve could be the new API but it doesn't have good
> > return types either.
>
> It might be helpful to enumerate some of those things. Here are a few
> I've noticed
>
> - In general, the output of solve() should not depend on the specific
> mathematical type of the input. This means for instance, linear and
> nonlinear equations should be handled in the exact same way. Right
> now, even dict=True can't save you:
>
> >>> solve([x**2 - 2, x > 1], x)
> Eq(x, sqrt(2))
> >>> solve([x**2 - 2, x > 1], x, dict=True)
> Eq(x, sqrt(2))
>
> One can imagine wanting to use x > 1 just to limit the solutions
> (which can't be done just with the old assumptions), but this has
> caused solve() to treat it as a completely different type of equation.

Most people who use an inequality with solve do not expect the current
behaviour. Indeed they usually want to filter a finite set of
solutions. However if inequalities can be included in the inputs then
there needs to be a way of representing unsolved inequalities in the
output format. In fact representing unsolved inequalities is needed
for other cases such as a solving a quadratic with symbolic
coefficients over the reals where the symbolic condition b^2 - 4ac > 0
needs to be satisfied.

> - One issue I see is that it's not always clear what solve()
> will/should do in some cases, especially undetermined cases. It seems
> like "solving" and "variable elimination" are too often mixed
> together. It would be better if they were completely separated, since
> it's not clear to me why you would want one when solve can't produce
> the other. We should very carefully define what it means to "solve" an
> equation (note this is not so trivial once you include systems,
> multiple variables, infinite solutions, and inequalities into the
> mix).

The first point that needs to be understood is that the output of
solve is for use in *substitutions*. The design decision that you want
something that can be used for substitutions places significant
constraints on how mathematically correct its output can be and what
kinds of problems can be solved. This is what solve is for though.
Mathematica has separate functions solve and reduce where solve
returns substitution rules and reduce simplifies a system of boolean
statements about integer/real/complex symbols. The output of reduce is
not directly usable for substitution but it can be always formally
correct. The design goal for solve is to be as correct as possible
while still returning an explicit list of substitution rules.

Understanding that the output is for substitutions we can see why:

1. The list of dicts output from solve(..., dict=True) is the best
output format from all of the existing solve functions because each
solution can be used directly with subs:

In [7]: equation = x**2 + x + 1

In [8]: sol1, sol2 = solve(equation, x, dict=True)

In [9]: print(sol1, sol2)
{x: -1/2 - sqrt(3)*I/2} {x: -1/2 + sqrt(3)*I/2}

In [10]: equation.subs(sol1).expand()
Out[10]: 0

2. The new solvers all failed to have usable output:

The output of solveset is completely unstructured so attempting to use
it for any further processing will fail very often:

In [18]: print(solveset(sin(x)-1, x))
ImageSet(Lambda(_n, 2*_n*pi + pi/2), Integers)

In [19]: print(solveset(x**2 - y, x, Reals))
Intersection({-sqrt(y), sqrt(y)}, Reals)

It's very unclear how to write some code that can extract substitution
rules from the Set object returned by solveset.

The linsolve function is a bit better in that the return type is
consistently a FiniteSet containing a single tuple. You can use it for
substitutions like this:

In [21]: (sol,) = linsolve(eq, [x, y])

In [22]: sol
Out[22]: (1/3, 1/3)

In [23]: rep = dict(zip([x, y], sol))

In [24]: rep
Out[24]: {x: 1/3, y: 1/3}

In [25]: eq[0].subs(rep)
Out[25]: 0

This is still an awkward format though because the FiniteSet isn't
indexable so you have to extract with (sol,) = or with list(sol)[0].
You also need to wrap it up with dict(zip(...)) before you can do any
substitution.

The FiniteSet does seem particularly pointless when it's guaranteed to
contain only a single tuple (why not just return the tuple?).
Supposedly the reason for doing this is to have a "consistent format"
across different solvers. However it isn't consistent because solveset
returns proper sets like ImageSet for infinite sets but linsolve still
returns a FiniteSet in that case:

In [26]: linsolve([x+y], [x, y])
Out[26]: {(-y, y)}

That's not a proper set and it can't be used properly with
intersections or other set operations. At the same time it's not a
directly usable output like a list of dicts so it's the worst of both
worlds: an awkward format that is also not mathematically correct.

With nonlinsolve at least returning a FiniteSet makes sense because
there can be a non-unique but still finite solution set:

In [27]: nonlinsolve([x**2 - 1, y], [x, y])
Out[27]: {(-1, 0), (1, 0)}

However solveset then does something really bizarre for infinite solutions:

In [28]: nonlinsolve([sin(x), y], [x, y])
Out[28]: {({2⋅n⋅π + π │ n ∊ ℤ}, 0), ({2⋅n⋅π │ n ∊ ℤ}, 0)}

Here you have a set of tuples but where an expression (0) is given for
y in each tuple we have a set ({2⋅n⋅π + π │ n ∊ ℤ}) as the
"expression" for x. Again this is not a mathematically correct set and
at the same time this is not a convenient format for further
processing and it is not clear how to use this to do substitutions.

> One idea is to make things programmatic with some sort of Solution
> object that provides different things like solution parameters in
> consistent attributes. This would also allow moving some of the more
> complex keyword-arguments from solve() itself into methods of the
> Solution object.

I could imagine something like having a solve function that just
returns a list of dicts but then a function like solve_full that
returns a more formally correct representation of a set of solutions.
The basic goal of what solveset tries to achieve is to represent the
different cases but the problem is that it is unstructured and does
not easily resolve into explicit expressions for the output. More
generally it should always be possible to represent the set of
solutions as a union of sets where each set is described more or less
like a Python list-comprehension and also is included only
conditionally:

S = {f(x1, x2, ...) for (x1, x2, ...) in (S1, S2, ...) where P(x1, x2,
...)} if condition

Examples of what this could look like:

>>> solve_full([x**2 - 1], [x])
[({x: 1}, {}, True, True),
({x: -1}, {}, True, True)]

>>> solve_full([sin(x)], [x])
[({x: n*p}, {n: Integers}, True, True)]

>>> solve_full([a*x**2 + b*x + c], [x], {x:Reals}])
[({x: (-b + sqrt(b**2-4*a*c))/(2*a)}, {}, True, (b**2-4*a*c > 0) & (a>0)),
({x: (-b - sqrt(b**2-4*a*c))/(2*a)}, {}, True, (b**2-4*a*c > 0) & (a>0))]

>>> solve_full([x + y, x>2], [x, y])
[({x: t, y: -t}, {t: Reals}, t>2, True)]

>>> solve_full([a*x - b], [x])
[({x: b/a}, {}, True, Ne(a, 0))]

Of course the representation above could be wrapped in a solution
object of some kind. It's also possible to inject much of the above
into the expressions themselves using ExprCondPair or new types of
symbols for the parameters n and t above.

Then the solve function can be a wrapper around this that just distils
it down to the list-of-dicts substitution rules with the understanding
that solve_full is needed to see full solution information such as the
conditional existence of solutions above.

A key question is whether solve is allowed to introduce new symbols.
At the moment solve doesn't do that but that makes it impossible to
correctly represent the solution set of basic equations like sin(x)=0
in an explicit form. If new symbols are to be returned then there has
to be some way to know what set those symbols are supposed to
represent (e.g. n is an integer) and to get those symbols somehow in
the output from the solving function.

Another question is how to handle solutions that exist for "most"
values of a parameter e.g. a solution for x that is valid provided a
does not equal 0. I think it's right that solve should return the
"generic" solution, however it should return a sufficient condition
describing that genericity like `Ne(a,0)`. None of the current solve
functions can do this. Also sometimes the system is generically
inconsistent e.g.:

In [31]: solve([x+y+a, x+y+1], [x, y])
Out[31]: []

Here there are no solutions unless a=1. There should be some way to
know that solve has assumed that a does not equal 1 and likewise
someway to tell solve that you want it to get the solution when a is
equal to 1 e.g. this doesn't work with current solve:

In [32]: solve([x+y+a, x+y+1, a-1], [x, y])
Out[32]: []

It might make more sense to have an assumptions=Eq(a, 1) argument
where various inequalities or other kinds of statements can be
assumed.

Another basic issue with solve is that it always tries to use roots to
get explicit radical expressions but RootOf is better when there are
no symbolic coefficients and in fact it is often better to have an
implicit description for the roots of a polynomial anyway e.g.:

>>> solve_full([exp(5*a*x) - exp(a*x) + 1], [x])
[({x: (log(t) + I*n*pi)/a}, {n: Integer, t: RootSet(t**5 - t + 1,
t)}, True, True)]

With implicit descriptions of roots we can handle polynomials of
arbitrary degree without the Abel-Ruffini limit. It's also much more
efficient in the case of large degree polynomials just to express the
solution in terms of a parameter t representing any of the roots. Of
course there should be some way to convert to an explicit list of the
roots using e.g. RootOf but a RootSet object could be correct even for
polynomials with symbolic coefficients and unknown total degree (e.g.
if the leading coefficient is not known to be nonzero as in a*t**5 - t
+ 1 or something like t**n - 1).

> - The API should also work uniformly for systems as well as single
> equations. Perhaps the cleanest way would be for solve() to always
> treat a single equation as a system of one equation. That way the
> return type is always uniform (the downside is the common case might
> involve an extra layer of indirection, which could be annoying).

The list of dicts form works just as well for a system as for a single
equation if you remember that the goal is substitution.

> There's also this old issue "solve is a mess" which discusses a lot of
> these things. https://github.com/sympy/sympy/issues/6659

Also:
https://github.com/sympy/sympy/issues/16861

--
Oscar

[Chris Smith]

For the case of `eqs = [x+y+a, x+y+1]` you get the solution for `a = 1` if you give it the freedom to tell you that: `solve(eqs) -> {a: 1, x: -y - 1}`. This gives you the conditions for a solution: a value for `a` and a relationship between `x` and `y`. 

As for extra parameters, do we ever need anything other than an integer and a Range to indicate the values for the integers? e.g. positive integer multiples of `pi` can be represented `Piecewise((i*pi, And(Eq(i, floor(i)), i > 0)))`.

/c

[Oscar Benjamin]

On Wed, 9 Mar 2022 at 19:29, Chris Smith <smi...@gmail.com> wrote:
>
> For the case of `eqs = [x+y+a, x+y+1]` you get the solution for `a = 1` if you give it the freedom to tell you that: `solve(eqs) -> {a: 1, x: -y - 1}`. This gives you the conditions for a solution: a value for `a` and a relationship between `x` and `y`.

In very simple cases that approach might seem to work but it's not
hard to show when it doesn't if you go beyond trivial examples:

In [22]: a, b, c, d, e, f, x, y = symbols('a:f, x, y')

In [23]: eqs = [a*x + b*y + c, a*x + b*y + d]

In [24]: solve(eqs, [x, y])
Out[24]: []

In [26]: print(solve(eqs))
[{a: -(b*y + d)/x, c: d}]

Clearly I didn't want to solve for a and c here: that's why solve has
an argument to say what the unknowns are.

--
Oscar

[Aaron Meurer]

I think the problem is that solve() is designed as an interactive
function. If you are sitting down in an interactive session and using
solve(), then its output makes sense. It has different return types,
but that's fine if you are just passing it something and looking
directly at the solution. It sometimes gives partial results, but this
is also fine in this setting because you can iterate. In these
situations, solve() is actually kind of OK, because it just tries to
do something helpful, although it's still not great, because with
solve(), it can be hard to know how to tell it to do exactly what you
do want.

But this whole idea falls apart once you want to work
programmatically, because it's harder to programmatically react to
different return types or different kinds of solution output,
especially if you don't even know what the possible options are that
you might get.

This is sort of similar to simplify() vs. targetted simplification
functions. If you use simplify(), you could get really anything. It's
SymPy's prerogative to do whatever it wants. This is fine if you are
working interactively because you can check yourself whether you like
what it is doing. But in situations where you aren't manually
eyeballing the result of every function call, it's a bad idea, and
you're better off trying to think about what the exact simplifications
you need are, and calling the appropriate functions for those.

This disconnect exists in a lot of places in SymPy. We might try to
document it better to make it clearer to users what the best practices
are for interactive use vs. the best practices for
programmatic/library use.

So I'm starting to wonder if the real fix here isn't so much to "fix
solve" (although solve() should definitely be improved and cleaned
up), but rather to treat solve() as the "interactive only" function
for equation solving, just as simplify() is the "interactive only"
function for simplification. Users who know what they want should be
encouraged to call more targeted solution routines, which only work on
a specific type of equation, but also always give a fixed known output
type with certain guarantees. solve() itself should be no more than a
wrapper on top of these with some heuristics to guess what you want,
just as simplify() just calls other simplification functions.

Aaron Meurer

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHVvXxRqw_1da01a%2BUvkeGkF5xA3ge%3DRuPYWv5UZDE-ez_reyg%40mail.gmail.com.

[Oscar Benjamin]

On Wed, 23 Mar 2022 at 19:42, Aaron Meurer <asme...@gmail.com> wrote:
>
> So I'm starting to wonder if the real fix here isn't so much to "fix
> solve" (although solve() should definitely be improved and cleaned
> up), but rather to treat solve() as the "interactive only" function
> for equation solving, just as simplify() is the "interactive only"
> function for simplification.

I don't think we can save solve this way. Even if I wanted to make an
"interactive" solve function I wouldn't want it to be so inconsistent.
Most of the return types are awkward even for interactive use.
Remember the number one thing someone wants to do with the output:
substitute it into something.

There are so many other problems with solve internally and externally
that whichever way you look at it the whole thing needs to be
rewritten from scratch. It's impossible to do that while preserving
the existing behaviour which is itself impossible to even define (the
docstring doesn't even try!).

The only actionable takeaway I see from the idea that solve could be
considered "interactive only" is that maybe backwards compatibility
could be ignored and it could be okay to change the output to be a
list of dicts always.

I agree with the idea that solve should be implemented as a wrapper
around more carefully defined routines but ultimately there still
always needs to be a "find explicit solutions to these arbitrary
equations if possible" routine and it should be possible to use that
in a programmatic way.

--
Oscar

[Chris Smith]

Except for booleans, consider single and system output:

If it is a single equation the output depends on whether symbols are given or not:

* list of expressions if a single variable is given (or the expression is univariate and no symbol is given)

* if multiple variables are given and the system is solved by linear undetermined coefficients or is split into real and imaginary parts you get a dictionary, otherwise you get a list of tuples of values

* a list of dictionaries if no variables are given (since otherwise you wouldn't know for which symbol the solution corresponds)

If it is a system of equations (defined as 1 or more equations in an iterable) then the solution is:

* a single dictionary if the simplified system was solved as a linear system

* if multiple variables are given and the system is solved as a linear system you get a dictionary (and extraneous symbols are ignored), otherwise a list of tuples of values

* a list of dictionaries if no variables were specified (or else you wouldn't know for which symbol the solution corresponds)

Although the docstring does not lay out the return cases as I have just done, it does try to show the situation under which the different types are returned. (I wrote them and I was mainly the messenger.) 

Any user that is concerned about the return value has had the option to select `set` or `dict`ionary output for several releases now. There should be no reason that it can't be used reliably in a programmatical way.

As has been pointed out, some of the issues arise from `solve` trying to guess what the user might have wanted. This is against the admonition to "refuse the temptation to guess". Others have pointed out the mixing of idioms of "eliminating" vs "solving".

/c

[Aaron Meurer]

Yes, it's true that there are still a lot of problems with solve().
But the point is that if you are going to go for a "one function"
approach, the best you can hope for is "you get what you get". There
aren't a dozen flags to simplify() that let you tweak what it does
here and what it does there (and in fact I see we should probably
fight this temptation. There's already "inverse" and "doit" flags,
which IMO shouldn't be there). If you use simplify(), and you don't
like what it does, you should instead use specific simplification
functions. Those functions each work in very specialized ways and do
let you precisely control what they do.

Similarly, if we want a single solve() function that handles solving a
single equation, systems of equations, overdetermined and
underdetermined equations (and systems), inequalities, diophantine
equations, equations with symbolic constants, solutions with
parameters, restricting solutions with assumptions, elimination, etc.
all at once, then I think it's reasonable to just let it do whatever
it does. If you don't like what operation it chose to do or how it
chose to represent the answer, then you should instead call the
specific function that does exactly what you want. That's a much
better API than trying to precisely control all these things at once
with keyword arguments, and trying to unify all those very different
things into a single return type.

For example, in my example above, a systems of inequalities solver
would treat [x**2 - 2, x > 1] as a system of inequalities (or error
because the first term is not an inequality). A single equation solver
would treat it as an equation with an assumptions restriction. Both
are valid things but importantly they are distinct things.

Another example is solving the equation |z| = 1 (abs(x) - 1), where z
is a complex number. Do you want re(z)**2 + im(z)**2 = 1? Or do you
want something like z = exp(I*theta), where theta is a new parameter
created by the solver (with "theta in [0, 2*pi)" specified somehow)?
Or do you want the set {x + I*y | x**2 + y**2 = 1}? Or something else?
All are legitimate things that you might want. And only one of them (z
= exp(I*theta)) actually makes sense with subs(), so I'm not even sure
if we can really say "the output of solve() should be passable to
subs()". Inequality solutions also aren't generally meaningful as
subs() inputs either.

Aaron Meurer

>
> --
> Oscar
>
> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHVvXxSM2JOLxthB7qr1sAUHTbVMj5TY18iwn0eAcRJJGPmW-g%40mail.gmail.com.

[Chris Smith]

From a historical perspective, seeing one inequality and treating all equations as inequalites came at a time when `Eq(expr)` was an allowed input meaning `Eq(expr,0)` so the interpretation of `(expr, inequality)` as a system of relationals was a meaningful interpretation at the time. (I think it makes more sense, now, to interpret the relationals in a mixed system as constraints on the solution. Or maybe that isn't the job of solve and the user should apply relational conditions on their own.)

/c

[Aaron Meurer]

On Wed, Mar 23, 2022 at 5:52 PM Chris Smith <smi...@gmail.com> wrote:
>
> Except for booleans, consider single and system output:
>
> If it is a single equation the output depends on whether symbols are given or not:
>
> * list of expressions if a single variable is given (or the expression is univariate and no symbol is given)
> * if multiple variables are given and the system is solved by linear undetermined coefficients or is split into real and imaginary parts you get a dictionary, otherwise you get a list of tuples of values
> * a list of dictionaries if no variables are given (since otherwise you wouldn't know for which symbol the solution corresponds)
>
> If it is a system of equations (defined as 1 or more equations in an iterable) then the solution is:
>
> * a single dictionary if the simplified system was solved as a linear system
> * if multiple variables are given and the system is solved as a linear system you get a dictionary (and extraneous symbols are ignored), otherwise a list of tuples of values
> * a list of dictionaries if no variables were specified (or else you wouldn't know for which symbol the solution corresponds)
>
> Although the docstring does not lay out the return cases as I have just done, it does try to show the situation under which the different types are returned. (I wrote them and I was mainly the messenger.)
>
> Any user that is concerned about the return value has had the option to select `set` or `dict`ionary output for several releases now. There should be no reason that it can't be used reliably in a programmatical way.
>
> As has been pointed out, some of the issues arise from `solve` trying to guess what the user might have wanted. This is against the admonition to "refuse the temptation to guess". Others have pointed out the mixing of idioms of "eliminating" vs "solving".

Another problem is how to return solutions that have a parameter,
where the parameter is a new variable. For example, the exp(I*theta)
solution of the abs(z) = 1 example I gave above. This also occurs with
most diophantine equations. At best it could just return the
parameterized solution, and hope the user knows how to extract it. But
even that only really works if the range of the parameter is something
you can represent using the old assumptions (e.g., you could add
Dummy('theta', real=True) but "theta in [0, 2pi)" wouldn't be
representable). Ideally the parameters would be returned separately
somehow, including their range of values as a set.

This was actually one of the original motivations for solveset
(consider solve(sin(x))), but ironically, even though such things are
representable from solveset, extracting the parameter from a set ends
up being harder than it would be if solve() just returned Dummy('n',
integer=True)*pi.

Aaron Meurer

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/8e0c6289-8e46-41e7-b755-5216eb9a2dc9n%40googlegroups.com.

[Chris Smith]

Is there an example of when a parameter other than an integer would be necessary?

/c

[Chris Smith]

Oscar gave an example in https://github.com/sympy/sympy/issues/16590

solve(x**2+y**2-1, [x, y]) -> [(cos(C₁), sin(C₁))]

This is a good example of where the purpose of solve should be clarified. I have always approached it as trying to solve for some (or all) of the given symbols explicitly in terms of others. I am not sure what I would call the above example - parametrization?

/c

[Chris Smith]

Here is the current logic of the input and output for `solve`

if equation is a boolean expression or a list contains a boolean
    boolean value or expression
        >>> solve((x**2-4, x>1))
        Eq(x, 2)
        >>> solve((x**2-4, x>2))
        False

else expression or list of expression(s) is passed
    and there is no solution
        return []
    else
        equation is single expression
            0 symbols
                univariate - list of values
                else - list of dictionaries
            1 symbol
                list of values
            > 1 symbol
                ordered symbols (e.g. list) -> list of tuples
                else -> list of dictionaries
        list of 1 or more expressions for equation with
            0 symbols
                linear in all variables - dict
                else - list of dictionaries
            1 symbol
                single solution - dict
                else - list of tuples
            > 1 symbol
                linear system and #eqs = #symbols -> dict
                ordered symbols (e.g., list) -> list of tuples
                else -> list of dictionaries

Except for booleans, using `dict=True` should return a list of dictionaries.

There is some inconsistency when using functions instead of symbols that I will try to fix in https://github.com/sympy/sympy/pull/23699.

/c

[Chris Smith]

In my most recent PR concerning `solve` Oscar has suggested (and I concur) that the varieties of output for `solve` not be changed -- though I am not sure that we are doing so for the same. Here is a draft for the documentation's "Explanations" section:

<p>

The output of the `solve` function can seem very unwieldy since it may appear to arbitrarily return one of 6 different types of output (in addition to raising errors). The reasons for this are historical and are biased toward human interaction rather than programatic use. The type of output will depend on the type of equation(s) (and how they are entered) and the number of symbols that are provided (and how they are provided).

Since this output is controllable by using the `dict=True` flag to give either a boolean expression or a list of dictionaries, the following is an explanation for those that might prefer the default output but wish to understand better why it is being given to make better sense of the output of an interactive session.

1. an empty list indicates that no solution was found

    >>> solve(sqrt(x) + 1)
    []

2. a list of values corresponding the solutions for a single variable which should be obvious in context because a) the equation was univariate or b) a single symbols was specified as being of interest.

    >>> solve(x**2 - 4)
    [-2, 2]

    >>> solve(x - y - 1, x)
    [y + 1]

3. a single dictionary with keys being symbols and values being the solutions for those symbols - this is the result when one or more equations passed as a list has only a single, unique solution (e.g. the system is linear or else is monotonic).

    >>> solve([x + y - 2, x - y + 2], x, y)
    {x: 0, y: 2}

4. a list of tuples, each giving a set of values for the symbols in the order they were given - this is the result when more than one symbol was given in a well defined order with 1) a single expression or 2) a list of one more more expressions was given and one or more of them was not linear/monotonic.

    >>> solve(x - 1, x, y)
    [(1, y)]
    >>> solve([x**2], x)
    [(0,)]
    >>> solve([x**2 - 1], x)
    [(-1,), (1,)]
    >>> solve([x**2 - 1], x, y)
    [(-1, y), (1, y)]

5. a list of dictionaries, each containing a solution for the variables appearing as keys in the dictionary - this is used when the expression (or an expression in the list) passed for ``f`` is neither linear nor monotonic and the symbols were either not provided or else provided in an unordered container like a list.

    >>> solve([exp(x) - 1, x**2])
    [{x: 0}]
    >>> solve([x + y - 2, x**2 - y + 2])
    [{x: -1, y: 3}, {x: 0, y: 2}]
    >>> solve([x + y - 2, x**2 - y + 2])
    [{x: -1, y: 3}, {x: 0, y: 2}]

6. a boolean expression - this will happen when an relational expression other than an Equality is given as an expression to solve. A single equality, or more complicated relational expression might be returned. Simple constraints to limit output to integers or positive values can be accomplished by setting the assumptions on the variables for which a solution is being sought, e.g. ``x = Symbol('x', positive=True)``, before defining the expression and passing it to `solve`.

    >>> solve([x**2 - 4, Ne(x,-2)])
    Eq(x, 2)

For those that used a version of SymPy older than [TBD], a complicating factor was that `solve` automatically tried to detect when an equation might be solved by matching coefficients on powers of an unspecified symbol. That would happen when a single expression was passed with all but one of the free symbols given as symbols for which a solution was desired. The output would be a single dictionary (if the symbols appeared in a linear fashion) or else a list of tuples if the system could be solved, otherwise an algebraic solution was sought for one or more of the symbols.

    >>> solve(a*x + b - (3*x + 4), a, b)
    {a: 3, b: 4}
    >>> solve(a*x + b**2 - (3*x + 4), a, b)
    [(3, -2), (3, 2)]

Now, only the algebraic solution for one or more of the specified symbols is sought. Determination of values for undetermined coefficients is obtained by explicitly passing the constraints on the coefficients or by calling `solve_undetermined_coeffs` directly.

    >>> solve_undetermined_coeffs(a*x + b - (3*x + 4), [a, b], x)
    {a: 3, b: 4}
    >>> solve([a - 3, b**2 - 4])
    [(3, -2), (3, 2)]
</p>

The problem with a global `solve` flag like is used for `evaluate` to force use the `dict=True` is that this would break internal code that depended on the current ouput. Instead, locally, at the beginning of an interactive session, the following could be defined:

    from sympy import solve as __solve
    def solve(a, *b, **c):
        from sympy.utilities.iterables import iterable
        if b:
            if len(b) == 1 and iterable(b[0]):
                return __solve(a, b[0], **c, dict=True)
            return __solve(a, *b, **c, dict=True)
        return __solve(a, **c, dict=True)

The first step in removing the automatic identification of undetetermined coefficient systems in implemented in PR #23699: a flag can disable this detection so an algebraic result is returned.

/c