first email introducing myselft

162 views
Skip to first unread message

Bruno Nicenboim

unread,
Feb 23, 2021, 1:08:18 PM2/23/21
to sympy
Hi,
This is my first email. I'm starting to use sympy, which I find fascinating, and I have some questions that couldn't be answered after going through the documentation and stackoverflow website.


  • level of familiarity with python (years of programming, previous projects)
I learned python maybe 15 years ago, but I have a weird level of familiarity with python. Now I mainly program with R, and I use python for things that cannot be done in R :)
  • mathematical education level (high school / ... / PhD?)
Mostly from my work in statistics and computational modeling
  • particular expertise (physics? biology? ...subtopics)
Cognitive science.
  • particular algorithmic interests
Statistics, probability distributions, computational models, signal processing
  • level of familiarity with symbolic math systems "computer algebra". Perhaps even course work in (say) algebra or books read.

  • your familiarity with SymPy (e.g. how have you used SymPy)?
I'm using it to approximate complex probability distributions
  • other possibly relevant information -- geographical location? native language?
I'm an assistant professor in Cognitive Science and AI in Tilburg University, more information about me can be found in https://bnicenboim.github.io/

Oscar Benjamin

unread,
Feb 23, 2021, 3:53:03 PM2/23/21
to sympy
On Tue, 23 Feb 2021 at 18:08, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>
> Hi,
> This is my first email. I'm starting to use sympy, which I find fascinating, and I have some questions that couldn't be answered after going through the documentation and stackoverflow website.

Hi Bruno and welcome!

Feel free to ask questions here. I didn't see any actual questions in
your email though...

Oscar

Bruno Nicenboim

unread,
Feb 23, 2021, 4:06:12 PM2/23/21
to sympy
Yeah, I just wanted to check that the mail goes through :)


I'm trying to impose a constraint on the variable that I want to solve for with sympy. (That is setting a domain that needs to be respected). And I want to use `solveset` (I've seen other answers that use `solve`).

Say I have the following `K_p` and I want to solve `K_p - x = 0` for `t`:

```
from sympy import *
psi, mu_1, mu_2, tau, t, x = symbols("psi, mu_1, mu_2, tau, t, x", positive = True, real = True)
K_p =(mu_1*mu_2*psi*t**2*tau**4 - 2*mu_1*mu_2*t*tau**2 - mu_1*psi*t*tau**2 + mu_1 - mu_2*psi*t*tau**2 + mu_2 + psi)/(mu_1*mu_2*t**2*tau**4 - mu_1*t*tau**2 - mu_2*t*tau**2 + 1)
```

How do I impose the following constraints:
`t < 1/(mu_1 * tau**2)` and  `t < 1/(mu_2 * tau**2)`

If I do this I have an extra solution that shouldn't be there:
```
s_x = solveset(K_p - x, t, domain=S.Reals)
```

I've tried this:
```
cset  = ConditionSet(t, t < 1/(mu_1 * tau**2)).intersect(ConditionSet(t, t < 1/(mu_1 * tau**2)))
s_x = solveset(K_p - x, t, domain=cset)
int = Interval(-oo, 1/(maximum(mu_1,mu_2)*tau**2))
s_x = solveset(K_p - x, t, domain=int)
```
But then sympy doesn't really solve it.

Notice that I manage to get which one of the two solutions with this dumb brute force approach:

```{python}
s_x = solveset(K_p - x, t, domain=S.Reals)
eq1 = s_x.args[0].args[1].args[0]
eq2 = s_x.args[0].args[1].args[1]

v_mu_1 = 2
v_mu_2 = 3
v_x = 50
v_psi = 10
v_tau =1
mint = min(1/(v_mu_1* v_tau**2), 1/(v_mu_2* v_tau**2))
eq1.subs([(mu_1, v_mu_1), (mu_2,v_mu_2), (tau, v_tau), (x, v_x), (psi, v_psi)]) < mint
eq2.subs([(mu_1, v_mu_1), (mu_2,v_mu_2), (tau, v_tau), (x, v_x), (psi, v_psi)]) < mint
```
Only one of the equations give me `True` as it should be. I would like  that the result of solveset would be that equation. Is this even possible?


 
Oscar

Oscar Benjamin

unread,
Feb 26, 2021, 3:01:05 PM2/26/21
to sympy
On Tue, 23 Feb 2021 at 21:06, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>
> On Tuesday, February 23, 2021 at 9:53:03 PM UTC+1 Oscar wrote:
>>
>> On Tue, 23 Feb 2021 at 18:08, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>> >
>> > Hi,
>> > This is my first email. I'm starting to use sympy, which I find fascinating, and I have some questions that couldn't be answered after going through the documentation and stackoverflow website.
>>
>> Hi Bruno and welcome!
>>
>> Feel free to ask questions here. I didn't see any actual questions in
>> your email though...
>>
> Yeah, I just wanted to check that the mail goes through :)
>
> I'm cross posting with stackoverflow where I didn't get any answer: https://stackoverflow.com/questions/66200069/restricting-the-domain-of-the-solution-of-solveset-with-sympy
>
> I'm trying to impose a constraint on the variable that I want to solve for with sympy. (That is setting a domain that needs to be respected). And I want to use `solveset` (I've seen other answers that use `solve`).
>
> Say I have the following `K_p` and I want to solve `K_p - x = 0` for `t`:
>
> ```
> from sympy import *
> psi, mu_1, mu_2, tau, t, x = symbols("psi, mu_1, mu_2, tau, t, x", positive = True, real = True)
> K_p =(mu_1*mu_2*psi*t**2*tau**4 - 2*mu_1*mu_2*t*tau**2 - mu_1*psi*t*tau**2 + mu_1 - mu_2*psi*t*tau**2 + mu_2 + psi)/(mu_1*mu_2*t**2*tau**4 - mu_1*t*tau**2 - mu_2*t*tau**2 + 1)
> ```
>
> How do I impose the following constraints:
> `t < 1/(mu_1 * tau**2)` and `t < 1/(mu_2 * tau**2)`
>
> If I do this I have an extra solution that shouldn't be there:
> ```
> s_x = solveset(K_p - x, t, domain=S.Reals)
> ```

I don't think your constraints are enough to determine what the
correct solution should be. For example if we have the values:
{mu_1:1, mu_2:2, tau:2, psi:1, x:2}

Then neither solution satisfies the constraint:

In [218]: rep = {mu_1:1, mu_2:2, tau:2, psi:1, x:S.Half}

In [219]: conditions = And(t < 1/(mu_1 * tau**2), t < 1/(mu_2 * tau**2))

In [220]: eq = K_p - x

In [221]: solve(eq.subs(rep))
Out[221]:
⎡11 √65 √65 11⎤
⎢── - ───, ─── + ──⎥
⎣16 16 16 16⎦

In [222]: r1, r2 = solve(eq.subs(rep))

In [223]: conditions.subs(rep).subs(t, r1)
Out[223]: False

In [224]: conditions.subs(rep).subs(t, r2)
Out[224]: False

--
Oscar

Bruno Nicenboim

unread,
Mar 1, 2021, 8:43:03 AM3/1/21
to sy...@googlegroups.com
On Fri, Feb 26, 2021 at 9:01 PM Oscar Benjamin <oscar.j....@gmail.com> wrote:
On Tue, 23 Feb 2021 at 21:06, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>
> On Tuesday, February 23, 2021 at 9:53:03 PM UTC+1 Oscar wrote:
>>
>> On Tue, 23 Feb 2021 at 18:08, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>> >
>> > Hi,
>> > This is my first email. I'm starting to use sympy, which I find fascinating, and I have some questions that couldn't be answered after going through the documentation and stackoverflow website.
>>
>> Hi Bruno and welcome!
>>
>> Feel free to ask questions here. I didn't see any actual questions in
>> your email though...
>>
> Yeah, I just wanted to check that the mail goes through :)
>
> I'm cross posting with stackoverflow where I didn't get any answer: https://stackoverflow.com/questions/66200069/restricting-the-domain-of-the-solution-of-solveset-with-sympy
>
> I'm trying to impose a constraint on the variable that I want to solve for with sympy. (That is setting a domain that needs to be respected). And I want to use `solveset` (I've seen other answers that use `solve`).
>
> Say I have the following `K_p` and I want to solve `K_p - x = 0` for `t`:
>
> ```
> from sympy import *
> psi, mu_1, mu_2, tau, t, x = symbols("psi, mu_1, mu_2, tau, t, x", positive = True, real = True)
> K_p =(mu_1*mu_2*psi*t**2*tau**4 - 2*mu_1*mu_2*t*tau**2 - mu_1*psi*t*tau**2 + mu_1 - mu_2*psi*t*tau**2 + mu_2 + psi)/(mu_1*mu_2*t**2*tau**4 - mu_1*t*tau**2 - mu_2*t*tau**2 + 1)
> ```
>
> How do I impose the following constraints:
> `t < 1/(mu_1 * tau**2)` and  `t < 1/(mu_2 * tau**2)`
>
> If I do this I have an extra solution that shouldn't be there:
> ```
> s_x = solveset(K_p - x, t, domain=S.Reals)
> ```

I don't think your constraints are enough to determine what the
correct solution should be. For example if we have the values:
{mu_1:1, mu_2:2, tau:2, psi:1, x:2}

Then neither solution satisfies the constraint:

yes, it's true. I realize that a missing piece of information is that psi < min(mu_1,mu_2,x). Is there a way to include this to get only one solution from solveset?

And is changing the domain of solveset the right approach?




In [218]: rep = {mu_1:1, mu_2:2, tau:2, psi:1, x:S.Half}

In [219]: conditions = And(t < 1/(mu_1 * tau**2), t < 1/(mu_2 * tau**2))

In [220]: eq = K_p - x

In [221]: solve(eq.subs(rep))
Out[221]:
⎡11   √65  √65   11⎤
⎢── - ───, ─── + ──⎥
⎣16    16   16   16⎦

In [222]: r1, r2 = solve(eq.subs(rep))

In [223]: conditions.subs(rep).subs(t, r1)
Out[223]: False

In [224]: conditions.subs(rep).subs(t, r2)
Out[224]: False

--
Oscar

--
You received this message because you are subscribed to a topic in the Google Groups "sympy" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sympy/oIKJ0nJHZC8/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sympy+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/CAHVvXxTzN5zMMUkNjXh%3Dbf5s4_mxVc%2BGzQpMCu1CmCQkyr5tVA%40mail.gmail.com.


--
Bests,

Bruno Nicenboim

Oscar Benjamin

unread,
Mar 1, 2021, 10:58:56 AM3/1/21
to sympy
On Mon, 1 Mar 2021 at 13:43, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>
> On Fri, Feb 26, 2021 at 9:01 PM Oscar Benjamin <oscar.j....@gmail.com> wrote:
>>
>> I don't think your constraints are enough to determine what the
>> correct solution should be. For example if we have the values:
>> {mu_1:1, mu_2:2, tau:2, psi:1, x:2}
>>
>> Then neither solution satisfies the constraint:
>
> yes, it's true. I realize that a missing piece of information is that psi < min(mu_1,mu_2,x). Is there a way to include this to get only one solution from solveset?
>
> And is changing the domain of solveset the right approach?

In principle yes, but I wouldn't be surprised if `solveset` is unable
to handle a case like this.

There isn't really a sympy function that handles this kind of case.
The solve function can handle a mix of equations and inequalities but
it doesn't treat them in the way that you want. The solveset function
is only for a single equation so you would have to encode the
inequality constraints in the domain somehow (which is not so easy for
the last constraint you have shown).

There has been discussion recently on this last about adding a solver
for systems of inequalities.

--
Oscar

Bruno Nicenboim

unread,
Mar 1, 2021, 3:01:47 PM3/1/21
to sy...@googlegroups.com
ok, thanks!

But for example here:
```
K_pp = tau**2*(2*mu_1**2*mu_2**2*t**2*tau**4 + mu_1**2 + mu_2**2 + t*(-2*mu_1**2*mu_2*tau**2 - 2*mu_1*mu_2**2*tau**2))/((mu_1*t*tau**2 - 1)**2*(mu_2*t*tau**2 - 1)**2)
```
This is another version without the psi, so the two constraints for t should be enough. But this is still giving me two solutions

```
cset  = ConditionSet(t, t < 1/(mu_1 * tau**2)).intersect(ConditionSet(t, t < 1/(mu_1 * tau**2)))
s_x = solveset(K_p - x, t, domain=cset)
```
As follows.
```
Intersection(ConditionSet(t, t < 1/(mu_1*tau**2)), {(-2*mu_1*mu_2 + mu_1*x + mu_2*x)/(2*mu_1*mu_2*tau**2*x) - sqrt(4*mu_1**2*mu_2**2 + mu_1**2*x**2 - 2*mu_1*mu_2*x**2 + mu_2**2*x**2)/(2*mu_1*mu_2*tau**2*x), (-2*mu_1*mu_2 + mu_1*x + mu_2*x)/(2*mu_1*mu_2*tau**2*x) + sqrt(4*mu_1**2*mu_2**2 + mu_1**2*x**2 - 2*mu_1*mu_2*x**2 + mu_2**2*x**2)/(2*mu_1*mu_2*tau**2*x)}) \ Intersection(ConditionSet(t, t < 1/(mu_1*tau**2)), {1/(mu_2*tau**2)})
```

I understand correctly that this is the right approach, but solveset is just unable to find the right solution given the specific domain?

There has been discussion recently on this last about adding a solver
for systems of inequalities.

--
Oscar

--
You received this message because you are subscribed to a topic in the Google Groups "sympy" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sympy/oIKJ0nJHZC8/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sympy+un...@googlegroups.com.

Oscar Benjamin

unread,
Mar 1, 2021, 3:31:45 PM3/1/21
to sympy
On Mon, 1 Mar 2021 at 20:01, Bruno Nicenboim <bruno.n...@gmail.com> wrote:
>
> ```
> cset = ConditionSet(t, t < 1/(mu_1 * tau**2)).intersect(ConditionSet(t, t < 1/(mu_1 * tau**2)))
> s_x = solveset(K_p - x, t, domain=cset)
> ```
> As follows.
> ```
> Intersection(ConditionSet(t, t < 1/(mu_1*tau**2)), {(-2*mu_1*mu_2 + mu_1*x + mu_2*x)/(2*mu_1*mu_2*tau**2*x) - sqrt(4*mu_1**2*mu_2**2 + mu_1**2*x**2 - 2*mu_1*mu_2*x**2 + mu_2**2*x**2)/(2*mu_1*mu_2*tau**2*x), (-2*mu_1*mu_2 + mu_1*x + mu_2*x)/(2*mu_1*mu_2*tau**2*x) + sqrt(4*mu_1**2*mu_2**2 + mu_1**2*x**2 - 2*mu_1*mu_2*x**2 + mu_2**2*x**2)/(2*mu_1*mu_2*tau**2*x)}) \ Intersection(ConditionSet(t, t < 1/(mu_1*tau**2)), {1/(mu_2*tau**2)})
> ```
>
> I understand correctly that this is the right approach, but solveset is just unable to find the right solution given the specific domain?

The output from solveset is I think correct in a formal sense. It is
just possible that solveset could do a better job of simplifying it to
get the result into a usable form. Conversely the input you are
providing to solveset is correct in a formal sense but is also not in
a form that solveset can easily use. That's essentially because you're
trying to use it for something that it has not really been designed
for.

The primary purpose of the domain argument to solveset is really to be
explicit about whether you want to solve over the reals or the complex
numbers for example. When you pass a complicated domain like this all
solveset will do is:
1. Check if the domain is a subset of the reals and if so compute
solutions over the reals
2. Return the intersection of the solution set with the domain
The intersection code will simplify the solution set where possible
but it only tries to handle simple cases. While it is possible to pass
an interval with symbolic conditions or a conditionset as the domain
that is not really intended as an alternative to having a solver that
can handle systems of inequalities. It might work in simple cases but
I wouldn't expect it to work in complicated cases.

I would approach this by solving the equation without any constraints
and then checking the conditions afterwards. That is at best all that
solveset does in this case but actually the way solveset approaches
this is more indirect and less likely to work because the conditions
are encoded in the domain. Checking the conditions yourself means you
can do it more systematically and you don't have the awkwardness of
trying to extract the expressions from the sets.


--
Oscar

Bruno Nicenboim

unread,
Mar 2, 2021, 12:11:18 AM3/2/21
to sy...@googlegroups.com
Understood! Thanks!

--
Oscar

--
You received this message because you are subscribed to a topic in the Google Groups "sympy" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sympy/oIKJ0nJHZC8/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sympy+un...@googlegroups.com.
Reply all
Reply to author
Forward
0 new messages