Algebraic rearrangement and obtaining coefficients of variables

[Nicholas Kinar]

Hello,

Suppose that I have an equation such as the following:

A*x**2 + B*y**2 + C*x*y + D*x**2 = E*z**2 + F

In the above, {A, B, C, D, E, F} are coefficients, and {x, y, z} are variables. Multiplication is denoted by the * symbol, addition by the + symbol, and powers by the ** symbol.  Also, = is the equals sign.

Is there a way to use sympy to manipulate this equation into the following form?

(A+D)*x**2 + B*y**2  = E*z**2 - C*x*y + F

Terms on the LHS of the equation are arranged in some sort of canonical order, and the terms have been grouped. What I would like to do is write some code to select the variables wanted on the LHS of the equation (x**2 and y**2), and then let sympy re-arrange the expression. How would I obtain the coefficients of the variables on the LHS of the expression?

Nicholas

[Aaron S. Meurer]

I think you want the collect() function and the .as_independent method.  The easiest way to manipulate Equalities is to first convert them into regular expressions.  Then, you would call collect on the expression using x**2.  collect collects all powers of the argument by default, but you can use the exact=False flag to turn this off (see the docstring of collect).  Actually, collect is optional, and is only if you want (A + D)*x**2 instead of A*x**2 + D*x**2.  

Then, you can gather the independent and dependent parts using as_independent.  Here is an example of how you would do it with your expression:

In [1]: var('A B C D E F')

Out[1]: (A, B, C, D, E, F)

In [2]: expr = Eq(A*x**2 + B*y**2 + C*x*y + D*x**2, E*z**2 + F)

In [3]: expr

Out[3]: 

           2      2      2          2

C⋅x⋅y + A⋅x  + B⋅y  + D⋅x  = F + E⋅z 

In [4]: expr = expr.lhs - expr.rhs

In [5]: expr

Out[5]: 

                2      2      2      2

-F + C⋅x⋅y + A⋅x  + B⋅y  + D⋅x  - E⋅z 

In [6]: expr = collect(expr, x**2, exact=True)

In [7]: expr

Out[7]: 

      2                      2      2

-F + x ⋅(A + D) + C⋅x⋅y + B⋅y  - E⋅z 

In [8]: expr.as_independent(x**2, y**2)

Out[8]: 

⎛                2   2              2⎞

⎝-F + C⋅x⋅y - E⋅z , x ⋅(A + D) + B⋅y ⎠

In [9]: expr = Eq(*reversed(expr.as_independent(x**2, y**2)))

In [10]: expr

Out[10]: 

 2              2                   2

x ⋅(A + D) + B⋅y  = -F + C⋅x⋅y - E⋅z 

It shouldn't be too hard to write a simple function that is passed an Equality and some variables and returns the result like in [10].

Aaron Meurer

--
You received this message because you are subscribed to the Google Groups "sympy" group.
To post to this group, send email to sy...@googlegroups.com.
To unsubscribe from this group, send email to sympy+un...@googlegroups.com.
For more options, visit this group at http://groups.google.com/group/sympy?hl=en.

[Nicholas Kinar]

Thanks, Aaron; that's exactly what I wanted and it is very exciting to
realize that this can be done. I've just started playing with sympy
and exploring its possibilities, so thank you very much for helping to
guide me and for your response!

Nicholas

[Nicholas Kinar]

I've been playing around with the algebraic rearrangement being
discussed in this post, and I am wondering if it is possible to do the
same operation with variables that are indexed. Let me perhaps
demonstrate what I would like to achieve.

In [1]: from sympy import *
In [2]: from sympy.tensor import Indexed, Idx
In [3]: var('A B C D')
In [4]: i, j, n, m = symbols('i j n m', integer=True)
In [5]: M = Indexed('M')
In [6]: expr = Eq(A*M(i,j)**(n+1) + B*M(i+1,j)**(n+1), C*M(i+1,j)**(n+1)
+ D*M(i,j)**n)

In [7]: expr
Out[7]: A*M(i, j)**(1 + n) + B*M(1 + i, j)**(1 + n) == D*M(i, j)**n +
C*M(1 + i, j)**(1 + n) + D*M(i, j)**(1 + n)

So essentially what I've done is created an equation with a variable M.
The variable M has superscripts (either 1+n or n), and subscripts
(i,j). I would now like to do the following:

(1) Set i = 1 and j = 1. These are just example integers, but i and j
could be set to other numbers. So then the expr should become:

A*M(1, 1)**(1 + n) + B*M(1 + 1, 1)**(1 + n) == D*M(1, 1)**n + C*M(1 + 1,
1)**(1 + n) + D*M(1, 1)**(1 + n)
A*M(1, 1)**(1 + n) + B*M(2, 1)**(1 + n) == D*M(1, 1)**n + C*M(2, 1)**(1
+ n) + D*M(1, 1)**(1 + n)

(2) Group all terms with similar subscripts and move all terms with (1 +
n) superscripts to the LHS of the equation:

(A - D)*M(1, 1)**(1 + n) + (B - C)*M(2, 1)**(1 + n) == D*M(1, 1)**n

(3) Print this new expression out in a "nice" format suitable for reading.

If anyone could point me in the right direction, I would be extremely
grateful. Could these operations be done using sympy?

Nicholas

[Aaron S. Meurer]

Hi again.

I'm not sure if exponentiation is the correct way to do superscripts. I guess you could use it, as an abstraction, but remember that SymPy is going to want to treat M(i, j)**(1 + n) as M(i, j)*M(i, j)**n. Probably a better idea would be to just put the "superscript" as another element of the index, such as M(i, j, 1 + n). The only disadvantage is that It won't print the same. Maybe we should implement special handling for super/subscripts in Indexed with printing support (feel free to send in a patch to do this).

Also, the things I say below will work much better if you make n an index of M.

>
> So essentially what I've done is created an equation with a variable M. The variable M has superscripts (either 1+n or n), and subscripts (i,j). I would now like to do the following:
>
> (1) Set i = 1 and j = 1. These are just example integers, but i and j could be set to other numbers. So then the expr should become:

Use .subs().

>
> A*M(1, 1)**(1 + n) + B*M(1 + 1, 1)**(1 + n) == D*M(1, 1)**n + C*M(1 + 1, 1)**(1 + n) + D*M(1, 1)**(1 + n)
> A*M(1, 1)**(1 + n) + B*M(2, 1)**(1 + n) == D*M(1, 1)**n + C*M(2, 1)**(1 + n) + D*M(1, 1)**(1 + n)
>
> (2) Group all terms with similar subscripts and move all terms with (1 + n) superscripts to the LHS of the equation:

So the same thing I said before would work, except collect() and .as_independent aren't so smart to handle more complex things like this. One idea I had was to use .atoms(), but it seems that it doesn't work as expected for Indexed:

In [98]: expr = Eq(A*M(i,j)**(n+1) + B*M(i+1,j)**(n+1), C*M(i+1,j)**(n+1) + D*M(i,j)**n)

In [100]: expr.atoms()
Out[100]: set(1, A, B, C, D, M, i, j, n)

In [101]: expr.atoms(Indexed)
Out[101]: set(M)

If the second one had given set([M(i, j) ,M(1 + i, j)]), then you could have had something to work with (Chris or Oyvind, do you know why it does this?).

Of course, the way I gave again above will still work, if you know what the M's are in the expression:

In [104]: expr = expr.lhs - expr.rhs

In [105]: expr.subs({i:1}).as_independent(M)
Out[105]:
⎛ n 1 + n 1 + n 1 + n⎞
⎝0, - D⋅M(1, j) + A⋅M(1, j) + B⋅M(2, j) - C⋅M(2, j) ⎠

In [106]: expr.subs({i:1}).as_independent(M(1, j))
Out[106]:
⎛ 1 + n 1 + n n 1 + n⎞
⎝B⋅M(2, j) - C⋅M(2, j) , - D⋅M(1, j) + A⋅M(1, j) ⎠

The problem is knowing that M(1, j) and M(2, j) are the indexed objects in the expression (that is where I would have liked to use atoms).

>
> (A - D)*M(1, 1)**(1 + n) + (B - C)*M(2, 1)**(1 + n) == D*M(1, 1)**n
>
> (3) Print this new expression out in a "nice" format suitable for reading.

You can use the same collect() idea that I gave above, but again, you have to give each full indexed object at this point in time.

>
> If anyone could point me in the right direction, I would be extremely grateful. Could these operations be done using sympy?
>
> Nicholas

Indexed() is kind of new, so most functions don't know about it. But I think if we fix the atoms bug above, what you want to do should be possible.

Aaron Meurer

[Nicholas Kinar]

> If the second one had given set([M(i, j) ,M(1 + i, j)]), then you could have had something to work with (Chris or Oyvind, do you know why it does this?).
>
> Of course, the way I gave again above will still work, if you know what the M's are in the expression:
>
> In [104]: expr = expr.lhs - expr.rhs
>
> In [105]: expr.subs({i:1}).as_independent(M)
> Out[105]:
> ⎛ n 1 + n 1 + n 1 + n⎞
> ⎝0, - D⋅M(1, j) + A⋅M(1, j) + B⋅M(2, j) - C⋅M(2, j) ⎠
>
> In [106]: expr.subs({i:1}).as_independent(M(1, j))
> Out[106]:
> ⎛ 1 + n 1 + n n 1 + n⎞
> ⎝B⋅M(2, j) - C⋅M(2, j) , - D⋅M(1, j) + A⋅M(1, j) ⎠
>
> The problem is knowing that M(1, j) and M(2, j) are the indexed objects in the expression (that is where I would have liked to use atoms).
>

It would have been neat to use atoms as well, and I'm also curious as to
whether this would be possible.

> You can use the same collect() idea that I gave above, but again, you have to give each full indexed object at this point in time.
>
>

Hmm, too bad. What I would like to do with sympy is generate equations
with many terms, and then call collect without knowing the full indexed
object.

> Indexed() is kind of new, so most functions don't know about it. But I think if we fix the atoms bug above, what you want to do should be possible.
>
> Aaron Meurer
>
>

Thanks Aaron; this is greatly appreciated! I only hope that I can get
proficient enough with using sympy that I'll be able to send in some
patches.

Nicholas

[Aaron S. Meurer]

.atoms() is the standard way to get all the similar kinds of objects in an expression, so it is a shame that it isn't working in this case. I opened issue (2058) for this (see http://code.google.com/p/sympy/issues/detail?id=2058).

Feel free to lend your help in fixing it. Something just as simple as finding why the issue occurs even without sending in a patch to fix it can help (though don't be afraid to send in a patch, too!).

Aaron Meurer

[Nicholas Kinar]

On 10-09-16 10:33 PM, Aaron S. Meurer wrote:
>
> .atoms() is the standard way to get all the similar kinds of objects in an expression, so it is a shame that it isn't working in this case. I opened issue (2058) for this (see http://code.google.com/p/sympy/issues/detail?id=2058).
>
> Feel free to lend your help in fixing it. Something just as simple as finding why the issue occurs even without sending in a patch to fix it can help (though don't be afraid to send in a patch, too!).
>
> Aaron Meurer
>
>

That's great; thanks again Aaron!

Nicholas

[Øyvind Jensen]

Exactly, Indexed is very fresh, and will probably change a bit before
the 0.7 release. See issue 2059 for instance. But as smichr figured
out in issue 2058, the .atoms() method works as it should. The
confusion comes from the fact that when an Indexed object recieves
indices, the returned object is of type IndexedElement. So you need to
use .atoms(IndexedElement)

It is very cool that you have started to use the tensor module,
Nicholas! And even cooler that you communicate the issues you have, so
that we can fix it before the next release :-)

Øyvind

> Aaron Meurer
>

[Nicholas Kinar]

> Exactly, Indexed is very fresh, and will probably change a bit before
> the 0.7 release. See issue 2059 for instance. But as smichr figured
> out in issue 2058, the .atoms() method works as it should. The
> confusion comes from the fact that when an Indexed object recieves
> indices, the returned object is of type IndexedElement. So you need to
> use .atoms(IndexedElement)
>
> It is very cool that you have started to use the tensor module,
> Nicholas! And even cooler that you communicate the issues you have, so
> that we can fix it before the next release :-)
>
> Øyvind
>
>

Thanks, Oyvind; I will give it a try and let the list know how I am
doing. I find that it is quite exciting to be able to do these very
interesting (and powerful) things with sympy.

Nicholas

[Nicholas Kinar]

Well, I've written some code in an attempt to re-arrange another example
expression, but for some reason this is not working. What am I doing
wrong here, and how might I collect terms with (n+1) powers and bring
them to the RHS of the equation? Here's what I have done:

# begin code
from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

M = Indexed('M')

var('A B C D E')
i,j,n = symbols('i j n', integer = True)
expr = Eq( A*M(i,j)**n + B*M(i+1,j)**(n+1),
C*M(i,j+1)**n + D*M(i,j)**(n+1) +
E*M(i+1,j)**(n+1) )

expr1 = expr.lhs - expr.rhs

# substitute for the expressions
# expr2 = expr1.subs([(i,1), (j,1)])
expr2 = expr1.subs( {i:1, j:1} )

# collect all (n+1) powers and bring them to the RHS of the equation
expr2.atoms( IndexedElement )
expr3 = collect( expr2, M**(n+1), exact=True )
print expr3

# end code

Here is the output of the print statement:

A*M(1, 1)**n - C*M(1, 2)**n + B*M(2, 1)**(1 + n) - D*M(1, 1)**(1 + n) -
E*M(2, 1)**(1 + n)

In the above code, collect() does not seem to group the B*M(2, 1)**(1 +
n) and E*M(2, 1)**(1 + n) terms.

Once again, thank you very much for continuing to point me in the right
direction.

Nicholas

[Øyvind Jensen]

> expr3 = collect( expr2, M**(n+1), exact=True )

This doesn't work because it is M(2, 1)**(n+1) that is present in expr2,
and not M**(n+1). So you must use:

expr3 = collect( expr2, M(2, 1)**(n+1), exact=True)

Øyvind

[Nicholas Kinar]

>
>> expr3 = collect( expr2, M**(n+1), exact=True )
>>
> This doesn't work because it is M(2, 1)**(n+1) that is present in expr2,
> and not M**(n+1). So you must use:
>
> expr3 = collect( expr2, M(2, 1)**(n+1), exact=True)
>

> �yvind
>
>

Thank you very much for your response, �yvind! That works well for me,
so once again, thank you!

Nicholas

[Nicholas Kinar]

Playing around with these concepts again, I've written another example
script:

# begin code

from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

p = Indexed('p')

i,j,n = symbols('i j n', integer = True)

exprA = Eq(p(1,1)**(n+1) + p(2,1)**(n+1), p(3,3)**(n+1) + p(1,3)**(n+1)
+ p(1,1)**n)
exprA.atoms(IndexedElement)

exprB = exprA.lhs - exprA.rhs
exprC = collect(exprB, p(1,1)**(n+1), exact = True)

vars_to_collect = [p(1,1)**(n+1), p(3,3)**(n+1)]

#error occurs here
exprD = Eq( *reversed( exprC.as_independent( vars_to_collect ) ) )
print exprD

#end code

However, this script generates the following error:

Traceback (most recent call last):
File "test.py", line 18, in <module>
exprD = Eq( *reversed( exprC.as_independent( vars_to_collect ) ) )
File "/usr/local/lib/python2.6/dist-packages/sympy/core/expr.py",
line 323, in as_independent
if term.has(*deps):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 923, in has
if p.is_Atom and not isinstance(p, Wild):
AttributeError: 'list' object has no attribute 'is_Atom'

What would I have to change to be able to ensure that the script will
run without errors?

Nicholas

[smichr]

But that only groups together one variety of M**(n+1). If you want
them all together you have to give them a common coefficient that can
be used for collection:

>>> u = Dummy('u')

These are the quantities that may have the n+1 power
>>> expr2.atoms(IndexedElement)
set([M(1, 2), M(1, 1), M(2, 1)])

So we replace them with u*M**(n+1) and colelct on u
>>> collect(expr2.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)
_u*(B*M(2, 1)**(1 + n) - D*M(1, 1)**(1 + n) - E*M(2, 1)**(1 + n)) +

A*M(1, 1)**n - C*M(1, 2)**n

Then split the expression into u-independent and u-dependent pieces
with as_independent and feed that back to Eq with the *prefix and
remove the u by replacing with 1.
>>> Eq(*_.as_independent(u)).subs(u, 1)
A*M(1, 1)**n - C*M(1, 2)**n == B*M(2, 1)**(1 + n) - D*M(1, 1)**(1 + n)

- E*M(2, 1)**(1 + n)

All the n+1 powers are on the right now.

[Nicholas Kinar]

That's exactly what I wanted; thank you very much smichr!

Nicholas

[Nicholas Kinar]

>
>>> u = Dummy('u')

I can't seem to use the u = Dummy('u') variable, since I now receive the
following error:

Traceback (most recent call last):

File "test.py", line 12, in <module>

u = Dummy('u')

NameError: name 'Dummy' is not defined

Moreover, I think that something like this would be extremely useful as
a function incorporated into the sympy distribution:

[Ondrej Certik]

On Sat, Sep 18, 2010 at 8:54 PM, Nicholas Kinar <n.k...@usask.ca> wrote:
>>
>>>> u = Dummy('u')
>
>
> I can't seem to use the u = Dummy('u') variable, since I now receive the
> following error:
>
> Traceback (most recent call last):
>  File "test.py", line 12, in <module>
>    u = Dummy('u')
> NameError: name 'Dummy' is not defined

I think that you need to use:

In [1]: a = Symbol("x", dummy=True)

>
> Moreover, I think that something like this would be extremely useful as a
> function incorporated into the sympy distribution:
>>>>
>>>> collect(expr2.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)

If you have time, send us a patch for this please.

Btw --- you should not be using "from sympy import *", but rather use
explicit imports, then it would be clear where things are coming from
and also Python can guide you more when some imports change.

Ondrej

[Øyvind Jensen]

> Moreover, I think that something like this would be extremely useful as
> a function incorporated into the sympy distribution:
> >>> collect(expr2.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)

Another approach would be to use Wild symbols when calling collect().
Like this:

>>> from sympy.tensor import Indexed

>>> M = Indexed('M')

>>> expr = x*M(k,m)**(n+1) + y*M(k,m)**(n+1) + M(k,m)**n + z*M(k,m)**n
>>> w = Wild('w1')
>>> collect(expr, w**(n+1))
M(k, m)**n + z*M(k, m)**n + (x + y)*M(k, m)**(1 + n)
>>> collect(expr, w**n)
(1 + z)*M(k, m)**n + x*M(k, m)**(1 + n) + y*M(k, m)**(1 + n)

However, this does not work. Does anyone know if there is a particular
reason for that, or is it a bug?

Øyvind

[Nicholas Kinar]

Thanks Ondrej and �yvind; that's greatly appreciated!

Ondrej:

As you suggest, I've attempted to use

a = Symbol("x", dummy=True)

In the example code given below, I am still using "from sympy import *" since I am not certain
exactly what I would have to import. However, I do agree that it is far better to use explicit imports.

�yvind: I like the idea of the wild symbols when calling collect; that's an interesting way of thinking, and thank
you for posting it.

Here's some example code where I've attempted to use the suggestions being posted in this thread:

# begin code
from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

M = Indexed('M')

var('A B C D E')
i,j,n = symbols('i j n', integer = True)

i = 1
j = 1

expr = Eq( A*M(i,j)**n + B*M(i+1,j)**(n+1),
C*M(i,j+1)**n + D*M(i,j)**(n+1) +
E*M(i+1,j)**(n+1) )
expr1 = expr.lhs - expr.rhs

expr1.atoms(IndexedElement)

a = Symbol("x", dummy=True)

print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)
# end code

However, this does not seem to work, since the following error is produced:

Traceback (most recent call last):

File "test-1.py", line 17, in<module>
print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)
NameError: name '_' is not defined

[smichr]

> a = Symbol("x", dummy=True)
> print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)
> # end code
>

If you try C.Dummy('u') that will work for you; Symbol('u',
dummy=True) is equivalent.

> However, this does not seem to work, since the following error is produced:
>
> Traceback (most recent call last):
>    File "test-1.py", line 17, in<module>
>      print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)
> NameError: name '_' is not defined

The underscore in an interactive session refers to the last output. In
a script it doesn't have that meaning; it's just a variable name which
you haven't defined yet.

You are also missing some of what happened in the example I gave.
Atoms was used to find the bases of the powers of interest and then
those, raised to the n+1 power, where replaced with themselves*u. What
you wrote above should appear something like this in a script:

...
u = C.Dummy('u')
IE_atoms = expr1.atoms(IndexedElement)
print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in IE_atoms]),
u)


And Øyvind, I don't think it's a bug that collect doesn't respect
Wild. It's just not implemented. Someone suggested that subs might do
substitutions using Wilds and I put that work into t2. And this thread
provided a test case for it, but doing expr.subs(W**(n+1), u*W**(n+1),
wild=1) failed because Idx (as it is being traversed by subs) doesn't
have methods that normal Basic objects do, e.g. as_base_exp(), so it
failed.

[Aaron S. Meurer]

That would be cool to have implemented, though. Or maybe some kind of lambda thing, so you could just say

collect(expr, lambda i: i.is_Indexed, lambda=True)

and it would collect all Indexed terms.

Also, as far as making the Indexed separation work, maybe some kind of lambda support in as_independent would work too, like

expr.as_independent(lambda i: i.is_Indexed and i[3] == n + 1, lambda=True)

I am assuming is_Indexed is an actual property of Basic, and that you are using M(i, j, n) instead of M(i, j)**n. (By the way, don't try these Nicholas. They don't work yet; they are just ideas).

Any ideas?

And of course, if Atoms worked correctly with Indexed (as per issue 2058), then you could just run a loop around collect() and as_independent().

Aaron Meurer

[Nicholas Kinar]

Thanks for your response, smichr; that's greatly appreciated.

>
>> However, this does not seem to work, since the following error is produced:
>>
>> Traceback (most recent call last):
>> File "test-1.py", line 17, in<module>
>> print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in _]), u)
>> NameError: name '_' is not defined
>>
> The underscore in an interactive session refers to the last output. In
> a script it doesn't have that meaning; it's just a variable name which
> you haven't defined yet.
>

Ah yes, I haven't been using iPython or something similar.

> You are also missing some of what happened in the example I gave.
> Atoms was used to find the bases of the powers of interest and then
> those, raised to the n+1 power, where replaced with themselves*u. What
> you wrote above should appear something like this in a script:
>
> ...
> u = C.Dummy('u')
> IE_atoms = expr1.atoms(IndexedElement)
> print collect(expr.subs([(a**(n+1), u*a**(n+1)) for a in IE_atoms]),
> u)
>
>

I've tried changing some sections of my script, but I still don't know
what I am doing wrong:

# begin code
from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

M = Indexed('M')
var('A B C D E')
i,j,n = symbols('i j n', integer = True)
i = 1
j = 1
expr = Eq( A*M(i,j)**n + B*M(i+1,j)**(n+1),
C*M(i,j+1)**n + D*M(i,j)**(n+1) +
E*M(i+1,j)**(n+1) )
expr1 = expr.lhs - expr.rhs

u = Symbol('u', dummy=True)
IE_atoms = expr1.atoms(IndexedElement)

print collect(expr1.subs([(a**(n+1), u*a**(n+1)) for a in IE_atoms]), u)
# begin code
# end code

After a few seconds, this script generates a very long output, ending
with the following error:

RuntimeError: maximum recursion depth exceeded while calling a Python object

So what would I have to change to make this script work as intended?

Nicholas

[Aaron S. Meurer]

So given Chris's tip in issue 2058 (which flew over my head the first time for some reason), I think I might have this figured out, except for one thing. Given M(i, j, n), how to I get the third argument, n? I tried M(i, j, n)[2] and M(i, j, n).args, but neither works.

Aaron Meurer

[Aaron S. Meurer]

I think you have found a bug with subs. Can you report it in the issues, with the traceback?

Aaron Meurer

[Nicholas Kinar]

> That would be cool to have implemented, though. Or maybe some kind of lambda thing, so you could just say
>
> collect(expr, lambda i: i.is_Indexed, lambda=True)
>
> and it would collect all Indexed terms.
>
> Also, as far as making the Indexed separation work, maybe some kind of lambda support in as_independent would work too, like
>
> expr.as_independent(lambda i: i.is_Indexed and i[3] == n + 1, lambda=True)
>
> I am assuming is_Indexed is an actual property of Basic, and that you are using M(i, j, n) instead of M(i, j)**n. (By the way, don't try these Nicholas. They don't work yet; they are just ideas).
>
> Any ideas?
>
> And of course, if Atoms worked correctly with Indexed (as per issue 2058), then you could just run a loop around collect() and as_independent().
>

I kind of like this idea as well. Perhaps it would be also possible for
the user to indicate on which side of an expression the terms should be
collected.

Nicholas

[Aaron S. Meurer]

Generally it's easiest to convert an Equality into a regular expression and then convert it back (like I showed in my first email to this thread). If you only want to manipulate one side of an Equality, you could do Eq(manipulation(expr.lhs), expr.rhs), i.e., only do it to whatever side and then rebuild the Equality.

Basically, Equality support in SymPy right now is shabby at best, because only a few functions really handle them correctly, so your best bet is to not use them in intermediate calculations.

Aaron Meurer

[Nicholas Kinar]

> I think you have found a bug with subs. Can you report it in the issues, with the traceback?
>
> Aaron Meurer
>

Done, the issue has been reported here
(http://code.google.com/p/sympy/issues/detail?id=2060). I expect that
it will be closed soon, but at least it gets everyone thinking ;-)

Nicholas

[Nicholas Kinar]

> Generally it's easiest to convert an Equality into a regular expression and then convert it back (like I showed in my first email to this thread). If you only want to manipulate one side of an Equality, you could do Eq(manipulation(expr.lhs), expr.rhs), i.e., only do it to whatever side and then rebuild the Equality.
>
> Basically, Equality support in SymPy right now is shabby at best, because only a few functions really handle them correctly, so your best bet is to not use them in intermediate calculations.
>
> Aaron Meurer
>
>

Agreed, thanks Aaron!

Nicholas

[Aaron S. Meurer]

Well, now I am having trouble reproducing the problem I had. I was getting:

In [23]: M(i, j, n).args
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)

/Users/aaronmeurer/Documents/Python/sympy/sympy/<ipython console> in <module>()

/Users/aaronmeurer/Documents/Python/sympy/sympy/sympy/tensor/indexed.pyc in __call__(self, *indices, **kw_args)
121
122 def __call__(self, *indices, **kw_args):
--> 123 return IndexedElement(self, *indices, **kw_args)
124
125 @property

/Users/aaronmeurer/Documents/Python/sympy/sympy/sympy/tensor/indexed.pyc in __new__(cls, stem, *args, **kw_args)
154
155 # FIXME: 2.4 compatibility

--> 156 args = map(_ensure_Idx, args)
157 # args = tuple([ a if isinstance(a, Idx) else Idx(a) for a in args ])

158 return Expr.__new__(cls, stem, *args, **kw_args)

/Users/aaronmeurer/Documents/Python/sympy/sympy/sympy/tensor/indexed.pyc in _ensure_Idx(arg)
135 return arg
136 else:
--> 137 return Idx(arg)
138
139 class IndexedElement(Expr):

/Users/aaronmeurer/Documents/Python/sympy/sympy/sympy/tensor/indexed.pyc in __new__(cls, label, range, **kw_args)
235
236 if not label.is_integer:
--> 237 raise TypeError("Idx object requires an integer label")
238
239 elif isinstance(range, (tuple, list, Tuple)):

TypeError: Idx object requires an integer label

But now I just get

In [21]: M(i, j, n).args
Out[21]: (M, i, j, n)

By the way, shouldn't it be just (i, j, n), like Function does?

In [20]: f(x, y).args
Out[20]: (x, y)

But anyway, now that it works, I think I can get somewhere.

Aaron Meurer

[Aaron S. Meurer]

OK, ignore that. I was using i as a loop variable (using i as a Symbol is a bad idea because of this, imho).

But now I get this strange error:

In [48]: M(i, j, n).args[3]
Out[48]: n

In [49]: M(i, j, n).args[3] == n
Out[49]: False

Aaron Meurer

[Øyvind Jensen]

sø., 19.09.2010 kl. 12.01 -0600, skrev Aaron S. Meurer:
> OK, ignore that. I was using i as a loop variable (using i as a Symbol is a bad idea because of this, imho).
>
> But now I get this strange error:
>
> In [48]: M(i, j, n).args[3]
> Out[48]: n
>
> In [49]: M(i, j, n).args[3] == n
> Out[49]: False

That is due to issue 2059.

Øyvind

[Øyvind Jensen]

>
> And Øyvind, I don't think it's a bug that collect doesn't respect
> Wild. It's just not implemented.

Yeah, that was kind of what I meant :-)

> Someone suggested that subs might do
> substitutions using Wilds and I put that work into t2. And this thread
> provided a test case for it, but doing expr.subs(W**(n+1), u*W**(n+1),
> wild=1) failed because Idx (as it is being traversed by subs) doesn't
> have methods that normal Basic objects do, e.g. as_base_exp(), so it
> failed.
>

But Idx is a subclass of Basic, and there are tests of .subs() for both
Idx, Indexed and IndexedElement, so it does work for simple cases.
However, it appears that .as_base_exp() is defined in Expr, not Basic.
If you implement features of subs() that depend on .as_base_exp(), maybe
those things should go into a Expr.subs()?

Øyvind

[Øyvind Jensen]

>
> > And Øyvind, I don't think it's a bug that collect doesn't respect
> > Wild. It's just not implemented. Someone suggested that subs might do
> > substitutions using Wilds and I put that work into t2. And this thread
> > provided a test case for it, but doing expr.subs(W**(n+1), u*W**(n+1),
> > wild=1) failed because Idx (as it is being traversed by subs) doesn't
> > have methods that normal Basic objects do, e.g. as_base_exp(), so it
> > failed.
>
> That would be cool to have implemented, though. Or maybe some kind of lambda thing, so you could just say
>
> collect(expr, lambda i: i.is_Indexed, lambda=True)
>
> and it would collect all Indexed terms.
>

The lambda function is very good idea. It would be extremely flexible,
and you could also use it with wildcards:

pattern = Wild('w1')**(n + 1)
collect(expr, lambda i: i.match(pattern))

But, I think it will make the learning curve much steeper if this is the
primary way to do these things. At least I had a hard time digesting
the functional aspects of python, and it took a gsoc to finally get it
in ;-)

Øyvind

[Aaron S. Meurer]

Unfortunately, you usually have to sacrifice learning curve for flexibility, especially in SymPy where we are limited to the Python language semantics. Of course, if it's just the lambda that is confusing, you can make it more explicit (at the cost of using more lines):

def is_index(i):
return i.is_Indexed

collect(expr, is_index, lambda=True)

Another alternative is to just use the wilds like you suggested, but to make Wild smarter so that it can encompass the lambda idea, like

w = Wild('w', lambda i: i.is_Indexed)
collect(expr, w**n)

I'm not so sure how easy that one would be to implement, but I think it would be preferable to the first one, actually.

Aaron Meurer

[Øyvind Jensen]

sø., 19.09.2010 kl. 14.43 -0600, skrev Aaron S. Meurer:
> Unfortunately, you usually have to sacrifice learning curve for flexibility, especially in SymPy where we are limited to the Python language semantics. Of course, if it's just the lambda that is confusing, you can make it more explicit (at the cost of using more lines):
>
> def is_index(i):
> return i.is_Indexed
>
> collect(expr, is_index, lambda=True)
>
> Another alternative is to just use the wilds like you suggested, but to make Wild smarter so that it can encompass the lambda idea, like
>
> w = Wild('w', lambda i: i.is_Indexed)
> collect(expr, w**n)

That is an excellent idea!

Actually, from the source code it looks like something similar to that
could be already implemented under the name of "properties". But
unfortunately the documentation is extremely sparse in that corner of
the code.

In the mean time I found an easy way to do the Wild() trick in collect,
and have uploaded to github, the branch is fix_collect. I am trying to
use the pull request feature, without any luck so far...

Øyvind

[Nicholas Kinar]

> In the mean time I found an easy way to do the Wild() trick in collect,
> and have uploaded to github, the branch is fix_collect. I am trying to
> use the pull request feature, without any luck so far...
>

> �yvind

That's very exciting, �yvind - it will be a neat experience to try your
patch. I've taken a look at github and I can't see the fix_collect
branch, so I guess that something will have to happen before the code
shows up. I take it that the Wild() trick is the same as the one that
you demonstrated in a previous post?

Nicholas

[Øyvind Jensen]

sø., 19.09.2010 kl. 17.19 -0600, skrev Nicholas Kinar:
> > In the mean time I found an easy way to do the Wild() trick in collect,
> > and have uploaded to github, the branch is fix_collect. I am trying to
> > use the pull request feature, without any luck so far...
> >

> > Ųyvind
>
> That's very exciting, Ųyvind - it will be a neat experience to try your

> patch. I've taken a look at github and I can't see the fix_collect
> branch, so I guess that something will have to happen before the code
> shows up. I take it that the Wild() trick is the same as the one that
> you demonstrated in a previous post?

Yes, it is that trick. Here is the branch:

http://github.com/jegerjensen/sympy/tree/fix_collect

If you want to, you can take part in the review process here:

http://github.com/sympy/sympy/pulls

Øyvind

>
> Nicholas
>
>

[Nicholas Kinar]

Done - I've tested these patches, and I've also signed up for github
access to make a comment. Thanks to all of you for the help that I've
received. Although I'm an old C/C++ programmer, I'm learning Python
just for the purpose of using sympy. With the help of this great
community, I've also found that the learning process is much easier.

Many thanks,

Nicholas

[Nicholas Kinar]

Working again on this example program, I am now trying to use
as_independent() on expr1 after collect() has been called on expr.

# begin code
from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

M = Indexed('M')
var('A B C D E')
i,j,n = symbols('i j n', integer = True)
i = 1
j = 1

expr = Eq( A*M(i,j)**n + A*M(i,j)**(n+1) + B*M(i+1,j)**(n+1),

C*M(i,j+1)**n + D*M(i,j)**(n+1) +
E*M(i+1,j)**(n+1) )

expr = expr.lhs - expr.rhs

w = Wild('w')
expr1 = collect(expr, w**(n+1))
print expr1.as_independent( w**(n+1) )
# end code

However, as_independent() does not seem to work with Wild variables,

since the following error is produced:

Traceback (most recent call last):

File "test-1.py", line 25, in <module>
print expr1.as_independent( w**(n+1) )
File "/usr/local/lib/python2.6/dist-packages/sympy/core/expr.py",
line 323, in as_independent
if term.has(*deps):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 928, in has
if isinstance(e, Basic) and e.has(p):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 928, in has
if isinstance(e, Basic) and e.has(p):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 928, in has
if isinstance(e, Basic) and e.has(p):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 925, in has
if p.matches(self) is not None:
File "/usr/local/lib/python2.6/dist-packages/sympy/core/power.py",
line 550, in matches
b, e = expr.as_base_exp()
AttributeError: 'Idx' object has no attribute 'as_base_exp'

Is there a way around this?

Nicholas

[Nicholas Kinar]

I am interested in re-arranging equations in this fashion since I would
like to create a function to generate systems for the solution of
finite-difference time domain (FDTD) schemes, which are very useful in
computational physics. I plan to create a function such as
generate_fdtd(expr) that accepts an expression and returns an expression
in the required form. Would such a function be useful in sympy, and
could it be plugged into the git repository as a patch?

Nicholas

[smichr]

On Sep 19, 10:20 pm, "Aaron S. Meurer" <asmeu...@gmail.com> wrote:
> So given Chris's tip in issue 2058 (which flew over my head the first time for some reason), I think I might have this figured out, except for one thing.  Given M(i, j, n), how to I get the third argument, n?  I tried M(i, j, n)[2] and M(i, j, n).args, but neither works.
>

>>> m = M(i,j,n)
>>> m.args
(M, i, j, n)
>>> m.args[-1]
n

[Aaron S. Meurer]

Yeah, see one of my previous emails on this thread. I was doing something stupid, and so it didn't work when I tried it. The problem now is with

>>> m.args[-1] == n
False

Aaron Meurer

[Øyvind Jensen]

ma., 20.09.2010 kl. 22.06 -0600, skrev Aaron S. Meurer:
> Yeah, see one of my previous emails on this thread. I was doing something stupid, and so it didn't work when I tried it. The problem now is with
>
> >>> m.args[-1] == n
> False
>

A workaround is to use :

>>> m.args[-1] == Idx(n)
True

There is actually also a test for this behaviour. But it was a
horrible design decision, hence issue 2059.

Øyvind

[Aaron S. Meurer]

OK, so here is how to do it (with all the work arounds). First off, again I HIGHLY RECOMMEND that you use M(i, j, n) instead of M(i, j)**n. As a simple example of how M(i, j)**n can lead to unintended consequences, consider what will happen if you set n = 0. You will get M(i, j)**0 == 1!

So, doing that, here is how you would do what you want with the expression given earlier. I skipped the collect() stage, but you already seem to know how to get that working (plus it will work better with Øyvind's upcoming patch).

In [68]: expr = Eq(A*M(i, j, n + 1) + B*M(i + 1, j, n + 1), C*M(i + 1, j, n+1) + D*M(i, j, n))

In [69]: expr = expr.lhs - expr.rhs

In [70]: expr = expr.subs({i:1, j:1})

In [71]: expr
Out[71]: A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n) - D⋅M(1, 1, n)

In [73]: filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement))
Out[73]: [M(1, 1, n)]

In [74]: expr.as_independent(*filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement)))
Out[74]: (A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n), -D⋅M(1, 1, n))

In [75]: Eq(*reversed(expr.as_independent(*filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement)))))
Out[75]: -D⋅M(1, 1, n) = A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n)

Now, it would also be great if as_independent() worked easier, without such a workaround. Maybe we should implement the Wild thing there too. Actually, I was surprised how easy it was to implement the Wild trick in collect(). It basically entailed modifying one line of code a little bit to have a .match() call in it. Probably as_independent() (and others!) would be just as easy to fix.

Aaron Meurer

[smichr]

>
> That would be cool to have implemented, though.  Or maybe some kind of lambda thing, so you could just say
>
> collect(expr, lambda i: i.is_Indexed, lambda=True)
>

Part of the problem I have run into in the past is a non-kludge way to
tell if you have a lambda argument. I think I ended up seeing if
'<lambda>' in str(arg) *...fiddling in interactive...*

OK, lambdas have a func_name method and the name is '<lambda>' so any
routine should be able to test an argument with:

>>> arg = lambda x: x
>>> if getattr(arg, 'func_name', None) == '<lambda>':
... print 'a lambda'
...
a lambda

Would it be worth putting an is_lambda function (not method since
lambda isn't a Basic object) in core that does the above?

[Nicholas Kinar]

> OK, so here is how to do it (with all the work arounds). First off, again I HIGHLY RECOMMEND that you use M(i, j, n) instead of M(i, j)**n. As a simple example of how M(i, j)**n can lead to unintended consequences, consider what will happen if you set n = 0. You will get M(i, j)**0 == 1!
>
> So, doing that, here is how you would do what you want with the expression given earlier. I skipped the collect() stage, but you already seem to know how to get that working (plus it will work better with Øyvind's upcoming patch).
>
> In [68]: expr = Eq(A*M(i, j, n + 1) + B*M(i + 1, j, n + 1), C*M(i + 1, j, n+1) + D*M(i, j, n))
>
> In [69]: expr = expr.lhs - expr.rhs
>
> In [70]: expr = expr.subs({i:1, j:1})
>
> In [71]: expr
> Out[71]: A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n) - D⋅M(1, 1, n)
>
> In [73]: filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement))
> Out[73]: [M(1, 1, n)]
>
> In [74]: expr.as_independent(*filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement)))
> Out[74]: (A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n), -D⋅M(1, 1, n))
>
> In [75]: Eq(*reversed(expr.as_independent(*filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement)))))
> Out[75]: -D⋅M(1, 1, n) = A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n)

Thanks, Aaron, that's very nice to know.

> Now, it would also be great if as_independent() worked easier, without such a workaround. Maybe we should implement the Wild thing there too. Actually, I was surprised how easy it was to implement the Wild trick in collect(). It basically entailed modifying one line of code a little bit to have a .match() call in it. Probably as_independent() (and others!) would be just as easy to fix.

Perhaps it would be possible to do this? I am certain that this
additional functionality would be extremely useful.

Nicholas

[smichr]

<copied from issue 2059>
Just a follow up about using subs(wild=True). If Idx is derived from
Expr rather than Basic, then the wild subs shown below (in t2) works:

>>> eq=A*M(i, j)**n + B*M(1 + i, j)**(1 + n)- D*M(i, j)**(1 + n)
>>> Eq(*reversed(collect(eq.subs(w**(1+n), w**(1+n)*u, wild=True), u).as_independent( u))).subs(u, 1) B*M(1 + i, j)**(1 + n) - D*M(i, j)**(1 + n) == A*M(i, j)**n
>>>

I still think the lambda for collection, as_independent (and
coefficient, atoms, has, subs, etc...) would be useful. Also,
something that would be nice in the wild=True subs would be if there
was a way to refer (as in a grep pattern) to the thing that was found.
What could be used for that -- the no-name Dummy as in C.Dummy(''):
>>> C.Dummy('')
_
>>> _.name
''

[Nicholas Kinar]

I agree - thanks smichr!

[Aaron S. Meurer]

That's why I included the lambda=True argument to collect() in that example.

But I think the Wild trick is better than the lambda idea anyway, and I think whatever lambda stuff is implemented should be implemented in Wild and matches().

Aaron Meurer

[Nicholas Kinar]

On 10-09-21 01:18 AM, Aaron S. Meurer wrote:
> OK, so here is how to do it (with all the work arounds). First off, again I HIGHLY RECOMMEND that you use M(i, j, n) instead of M(i, j)**n. As a simple example of how M(i, j)**n can lead to unintended consequences, consider what will happen if you set n = 0. You will get M(i, j)**0 == 1!
>

Thanks again, Aaron; I agree that there are unwanted consequences when n
= 0.

> So, doing that, here is how you would do what you want with the expression given earlier. I skipped the collect() stage, but you already seem to know how to get that working (plus it will work better with Øyvind's upcoming patch).
>

Hmm, well - perhaps I'm not running collect properly either:

w = Wild('w')

expr = collect(expr, M(w,w,n+1) )
print expr

But as the traceback shows this didn't work:

File "test-2.py", line 16, in <module>
expr = collect(expr, M(w,w,n+1) )
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
123, in __call__
return IndexedElement(self, *indices, **kw_args)
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
156, in __new__
args = map(_ensure_Idx, args)
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
137, in _ensure_Idx
return Idx(arg)
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
237, in __new__

raise TypeError("Idx object requires an integer label")

TypeError: Idx object requires an integer label

> In [68]: expr = Eq(A*M(i, j, n + 1) + B*M(i + 1, j, n + 1), C*M(i + 1, j, n+1) + D*M(i, j, n))
>
> In [69]: expr = expr.lhs - expr.rhs
>
> In [70]: expr = expr.subs({i:1, j:1})
>
> In [71]: expr
> Out[71]: A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n) - D⋅M(1, 1, n)
>
> In [73]: filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement))
> Out[73]: [M(1, 1, n)]
>
> In [74]: expr.as_independent(*filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement)))
> Out[74]: (A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n), -D⋅M(1, 1, n))
>
> In [75]: Eq(*reversed(expr.as_independent(*filter(lambda t: t.args[-1] == Idx(n), expr.atoms(IndexedElement)))))
> Out[75]: -D⋅M(1, 1, n) = A⋅M(1, 1, 1 + n) + B⋅M(2, 1, 1 + n) - C⋅M(2, 1, 1 + n)
>
> Now, it would also be great if as_independent() worked easier, without such a workaround. Maybe we should implement the Wild thing there too. Actually, I was surprised how easy it was to implement the Wild trick in collect(). It basically entailed modifying one line of code a little bit to have a .match() call in it. Probably as_independent() (and others!) would be just as easy to fix.
>
> Aaron Meurer
>
>
>

This is an excellent idea, and I've tried to implement it. I don't know
if it should be tricky to implement, but as smichr mentions earlier in
this thread, Idx does not have 'as_base_exp' present. So of course,
when I try the following,

w = Wild('w')

expr = collect(expr, w**(n+1))
print expr.as_independent( w**(n+1) )

I receive an error in the traceback:

Traceback (most recent call last):

File "test-1.py", line 21, in <module>

print expr1.as_independent( w**(n+1) )
File "/usr/local/lib/python2.6/dist-packages/sympy/core/expr.py",
line 323, in as_independent
if term.has(*deps):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 928, in has
if isinstance(e, Basic) and e.has(p):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 928, in has
if isinstance(e, Basic) and e.has(p):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 928, in has
if isinstance(e, Basic) and e.has(p):
File "/usr/local/lib/python2.6/dist-packages/sympy/core/basic.py",
line 925, in has
if p.matches(self) is not None:
File "/usr/local/lib/python2.6/dist-packages/sympy/core/power.py",
line 550, in matches
b, e = expr.as_base_exp()
AttributeError: 'Idx' object has no attribute 'as_base_exp'

Now I believe that as_base_exp() and other related functions can be
found in expr.py, and I've tried copying these functions into the Idx
class (indexed.py) to serve as a template for what I would like to do.
I've then attempted to modify these functions, but based on these very
preliminary excursions I think that something else needs to be done to
fix as_independent().

Nicholas

[Aaron S. Meurer]

As Chris pointed out, the correct fix is to make Idx subclass from Expr.

Aaron Meurer

[Nicholas Kinar]

>> Now I believe that as_base_exp() and other related functions can be found in expr.py, and I've tried copying these functions into the Idx class (indexed.py) to serve as a template for what I would like to do. I've then attempted to modify these functions, but based on these very preliminary excursions I think that something else needs to be done to fix as_independent().
>>
>> Nicholas
>>
> As Chris pointed out, the correct fix is to make Idx subclass from Expr.
>
> Aaron Meurer
>
>

Thanks Aaron; so what I've done is changed the following line in
indexed.py from

class Idx(Base):

to the following:

class Idx(Expr):

However, I now receive another error:

Traceback (most recent call last):

File "test-2.py", line 16, in <module>
expr = collect(expr, M(w,w,n+1) )
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
123, in __call__
return IndexedElement(self, *indices, **kw_args)
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
156, in __new__
args = map(_ensure_Idx, args)
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
137, in _ensure_Idx
return Idx(arg)
File
"/usr/local/lib/python2.6/dist-packages/sympy/tensor/indexed.py", line
237, in __new__
raise TypeError("Idx object requires an integer label")
TypeError: Idx object requires an integer label

I received this same error when manually copying the functions between
the files, which isn't too extraordinary (given that we're subclassing
one class from the other).

Nicholas

[Nicholas Kinar]

>
> I received this same error when manually copying the functions between
> the files, which isn't too extraordinary (given that we're subclassing
> one class from the other).
>
> Nicholas
>
>

My mistake; I was running the wrong Python script! Now by subclassing
Idx from Expr, the following works quite nicely!

# begin code
from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

M = Indexed('M')
var('A B C D E')
i,j,n = symbols('i j n', integer = True)
i = 1
j = 1
expr = Eq( A*M(i,j)**n + A*M(i,j)**(n+1) + B*M(i+1,j)**(n+1),
C*M(i,j+1)**n + D*M(i,j)**(n+1) +
E*M(i+1,j)**(n+1) )

expr = expr.lhs - expr.rhs

# print expr

w = Wild('w')

expr1 = collect(expr, w**(n+1))
print Eq( *reversed( expr1.as_independent( w**(n+1) ) ) )
# end code

The output is as expected!

(A - D)*M(1, 1)**(1 + n) + (B - E)*M(2, 1)**(1 + n) == A*M(1, 1)**n -
C*M(1, 2)**n

I would really like to see this fixed in the trunk, and as you say, we
can probably do this for other things as well!

Nicholas

[smichr]

Another way of handling the lambda functionality would be to piggy
back off of match so that eq.match(lambda) would return True if the
lambda was true for eq. That way routines would gain matching and
lambda capabilities at the same time. Or, like for .has(), routines
already having wild matching would gain lambda matching.

[Nicholas Kinar]

> I am interested in re-arranging equations in this fashion since I
> would like to create a function to generate systems for the solution
> of finite-difference time domain (FDTD) schemes, which are very useful
> in computational physics. I plan to create a function such as
> generate_fdtd(expr) that accepts an expression and returns an
> expression in the required form. Would such a function be useful in
> sympy, and could it be plugged into the git repository as a patch?
>
> Nicholas
>

So given that Idx can subclass from Expr as suggested by smichr, and
Oyvind's fix gets pushed into the trunk, my proposed (extremely modest)
function would become:

def generate_fdtd(expr, nsolve):

w = Wild('w')

expr0 = expr.lhs - expr.rhs
expr1 = collect(expr0, w**nsolve)
return Eq( *reversed( expr1.as_independent( w**nsolve ) ) )

Thanks to everyone who helped me out with this; it is greatly appreciated!

The function takes an expression expr and the next timestep nsolve,
which is usually given as n+1. All terms with n+1 superscripts are then
collected and taken to the LHS of the equation. The superscripts are
not powers in this type of specialized notation, but bearing in mind
what Aaron has said, it would also be nice to use the M(i,j,n) format.

Is there now a way to iterate through each term of the returned
expression and then find the

(1) coefficient of each term
(2) variable of each term ( such as M(1,1) ) for indexed variables?

Ideally I would like to use a loop, and then test for equality between a
known variable, such as

if term == M(1,1):

Is this possible?

Nicholas

[Nicholas Kinar]

+1, I like this as well, and it would be very neat to do.

Nicholas

[Aaron S. Meurer]

I think it would be better design wise if the lambda were entered in the Wild itself, so instead use something like:

w = Wild('w', lambda=lambda i: i.is_Whatever and <some other stuff maybe>)

It also gives more flexibility this way, because you can have just a part of the match expression use the lambda.

Aaron Meurer

[Nicholas Kinar]

>
> So given that Idx can subclass from Expr as suggested by smichr, and
> Oyvind's fix gets pushed into the trunk, my proposed (extremely
> modest) function would become:
>
> def generate_fdtd(expr, nsolve):
> w = Wild('w')
> expr0 = expr.lhs - expr.rhs
> expr1 = collect(expr0, w**nsolve)
> return Eq( *reversed( expr1.as_independent( w**nsolve ) ) )
>
> Thanks to everyone who helped me out with this; it is greatly
> appreciated!
>
> The function takes an expression expr and the next timestep nsolve,
> which is usually given as n+1. All terms with n+1 superscripts are
> then collected and taken to the LHS of the equation. The superscripts
> are not powers in this type of specialized notation, but bearing in
> mind what Aaron has said, it would also be nice to use the M(i,j,n)
> format.

However, this doesn't seem to work for more complicated FDTD stencils,
since terms are not distributed through the equation. So I added one
line and changed my proposed function to the following:

# function to generate the FDTD scheme

def generate_fdtd(expr, nsolve):
w = Wild('w')
expr0 = expr.lhs - expr.rhs

expr1 = expr0.expand()
expr2 = collect(expr1, w**nsolve)
return Eq( *reversed( expr2.as_independent( w**nsolve ) ) )

>
> Is there now a way to iterate through each term on the LHS of the

> returned expression and then find the
>
> (1) coefficient of each term
> (2) variable of each term ( such as M(1,1) ) for indexed variables?
>
> Ideally I would like to use a loop, and then test for equality between
> a known variable, such as
>
> if term == M(1,1):
>
> Is this possible?
>

Perhaps there would be a way to access each of the terms in the
expression in the same fashion as elements in a list?

[Aaron S. Meurer]

Try using .args, as in

expr.args

Aaron Meurer

[Nicholas Kinar]

>
> Try using .args, as in
>
> expr.args
>
> Aaron Meurer
>

Thanks, Aaron; that's the icing on the cake, and it is exactly what I
needed. Isn't it neat how expr.args can be cascaded?

expr1 = expr0.args[0].args[1].args[1]

Nicholas

[Aaron S. Meurer]

Yep. .args is basically the key to doing advanced stuff in SymPy. You can access any part of an expression you want by walking the .args tree, and you can use it to do any kind of manipulation you want. Basically, if you look at the source for many of these functions like collect() or as_independent(), they are doing nothing more than walking .args and moving things around there.

Happy coding!

Aaron Meurer

[Nicholas Kinar]

> Yep. .args is basically the key to doing advanced stuff in SymPy. You can access any part of an expression you want by walking the .args tree, and you can use it to do any kind of manipulation you want. Basically, if you look at the source for many of these functions like collect() or as_independent(), they are doing nothing more than walking .args and moving things around there.
>
> Happy coding!
>
> Aaron Meurer
>
>

To you as well, Aaron!

Many thanks,

Nicholas

[smichr]

On Sep 22, 8:18 am, "Aaron S. Meurer" <asmeu...@gmail.com> wrote:
> I think it would be better design wise if the lambda were entered in the Wild itself, so instead use something like:

I won't be working on this until I finish the t2 work and the solver
mods so someone else might want to look at this. Since I've heavily
modified lots of the core routines it would be great if it could be
done on top of t2 to make the merging process easier.

[Nicholas Kinar]

I've written another example program :-D

What I would like to do is separate each term into coefficients and
variables.

Is it possible to obtain the coefficients of a term using .coeff without
knowing the indexed variable? The following example program
demonstrates what I would like to achieve, but using Wild with .coeff
does not seem to work.

What would I have to change in my example program to:

(1) use .coeff with Wild
(2) obtain the indexed variable p(1, 3)**(1 + n), given term?

However, perhaps there is a better way to achieve what I am looking for?

Nicholas

# begin code
from sympy import *
from sympy.tensor import Indexed, Idx, IndexedElement

def main():
p = Indexed('p')
var('deltax deltay delta deltat A B')

i,j,n = symbols('i j n', integer = True)

LHS = ((deltax**(-4) + 0.5/(deltax**2*deltay**2))*p(3, 1)**(1 + n) +
(3/deltax**4 + 0.5/(deltax**2*deltay**2) +
A/(deltat**2*deltax**2) +
B/(2*deltat*deltax**2))*p(1, 1)**(1 + n) +
(-4/deltax**4 - 1.0/(deltax**2*deltay**2) -
A/(deltat**2*deltax**2)
- B/(2*deltat*deltax**2))*p(2, 1)**(1 + n) + 0.5*p(1,
3)**(1 + n)/(deltax**2*deltay**2)
+ 0.5*p(3, 3)**(1 + n)/(deltax**2*deltay**2) + 2.0*p(2,
2)**(1 + n)/(deltax**2*deltay**2)
- 1.0*p(1, 2)**(1 + n)/(deltax**2*deltay**2) - 1.0*p(2,
3)**(1 + n)/(deltax**2*deltay**2)
- 1.0*p(3, 2)**(1 + n)/(deltax**2*deltay**2)
)

print ''
print 'This is the expression displayed in printing order'
print ''
print LHS
print ''

print 'This is one of the terms of the expression'
print 'Note that .args does not seem to select parameters'
print 'of an instance in the same order as printing order'
print ''
term = LHS.args[3]
print term
print ''

print 'This does not seem to work when finding coefficients'
print ''
print term.args[0]
print term.args[1]
print term.args[2]
print term.args[3]

print ''
print 'But we can grab the coefficient if p(1, 3)**(1 + n) is known'
print 'and this works well'
print ''
term = term.coeff( p(1, 3)**(1 + n) )
print term

print ''
print 'Using Wild does not seem to work'
print ''

w = Wild('w')

print term.coeff( w )

print ''
print 'Now given term, what would I have to do to obtain p(1,
3)**(1 + n) '
print 'without knowing that this indexed variable is in the term?'

if __name__ == "__main__":
main()

#end code