MatrixSymbol and diff function

241 views
Skip to first unread message

dibus2

unread,
Aug 22, 2013, 9:10:19 AM8/22/13
to sy...@googlegroups.com
Hi all,

I am trying to figure out why I cannot do the following anymore with sympy 0.7.3 while I could with sympy 0.7.2 :


A = MAtrixSymbol('A','d',1)
xpr = A[3,0] * A[0,0] * 6
xpr.diff(A[0,0])

returns 6*A30 in sympy 0.7.2 but returns can't differentiate with respect to A[0,0]

How can I do the same thing with sympy 0.7.3 ?

Thanks,

Cheers,

F.

Matthew Rocklin

unread,
Aug 22, 2013, 9:55:58 AM8/22/13
to sy...@googlegroups.com
We had a slight regression with elements of MatrixSymbols.  A substantially beneficial code refactoring left us with an expression object that sometimes fails.  We discussed solutions for this in a few PRs (listed below) but failed to resolve this issue before the 0.7.3 release (I was in crunch mode with other projects at the time.)   I'll take a look at this issue later this afternoon.


--
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 post to this group, send email to sy...@googlegroups.com.
Visit this group at http://groups.google.com/group/sympy.
For more options, visit https://groups.google.com/groups/opt_out.

Matthew Rocklin

unread,
Aug 22, 2013, 9:57:25 AM8/22/13
to sy...@googlegroups.com

Matthew Rocklin

unread,
Aug 22, 2013, 12:23:36 PM8/22/13
to sy...@googlegroups.com
Fortunately this issue does not descend from the larger issue and can be solved independently.  See https://github.com/sympy/sympy/pull/2402

dibus2

unread,
Aug 22, 2013, 12:58:58 PM8/22/13
to sy...@googlegroups.com
Alright I had a look at it and tried to run my program which uses it extensively and it does not complain anymore. Therefore I think that it works.

Thanks a lot.

Cheers,

F.

Matthew Rocklin

unread,
Aug 22, 2013, 12:59:40 PM8/22/13
to sy...@googlegroups.com
Happy to hear it. 

Clifford Wolf

unread,
Apr 21, 2014, 7:50:51 AM4/21/14
to sy...@googlegroups.com
Hello,

On Thursday, August 22, 2013 6:23:36 PM UTC+2, Matthew wrote:
Fortunately this issue does not descend from the larger issue and can be solved independently.  See https://github.com/sympy/sympy/pull/2402

I think setting MatrixElement._diff_wrt = True was a bad decision.


When running "expr.diff(wrt)", then the code in Derivative.__new__ that handles _diff_wrt works under the assumption that the terms in expr are either independent from wrt or contain an identical replica. Consider code like the following:

i,n = symbols('i,n')
w
= MatrixSymbol('w', n, 1)
result
= Sum(123*w[i, 0], (i, 0, n - 1)).diff(w[0,0])

Derivative.__new__ will replace w[0,0] with a dummy symbol and attempt to also replace it in the expression. But in the expression we only have w[i,0] and Expr.subs does not understand how w[i,0] relates to w[0,0], yielding something like:

result = Sum(123*w[i, 0], (i, 0, n - 1)).diff(_x1)

Therefore result is set to 0 instead of 123.

I can understand the need for supporting differentiation with respect to MatrixElements, but I think it can't be done with a simple Expr._diff_wrt property. Instead we would need something like MatrixElement._diff_subs_wrt(self, expr) that returns a tuple of the new expression and the dummy variable, or None if differentiation is not possible. Such a method could then check for terms in expr that might or might not overlap with wrt in a way that can not be handled correctly at the moment, or even rewrite the expression (for example unroll the case i==0 out of the sum) to make differentiation possible.

regards,
 - clifford

PS: I'm new to sympy (this is the first time I looked into the sympy source code), so I don't really have a much insight yet. If what I say makes no sense, it might be because I have no idea what I am talking about.. ;)

Aaron Meurer

unread,
Apr 24, 2014, 7:16:19 PM4/24/14
to sy...@googlegroups.com
If we override _eval_derivative on MatrixElement can we make it do the
right thing?

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 post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/sympy/b18af010-e570-4182-a5a4-389b159d5213%40googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

Matthew Rocklin

unread,
Apr 24, 2014, 8:49:15 PM4/24/14
to sy...@googlegroups.com
I haven't thought about this issue much.  Generally speaking though MatrixElement and derivatives and Matrix Expressions are absolutely in a not so perfect state.  Help here is welcome.


Clifford Wolf

unread,
Apr 25, 2014, 6:23:13 AM4/25/14
to sy...@googlegroups.com
On Friday, April 25, 2014 1:16:19 AM UTC+2, Aaron Meurer wrote:
If we override _eval_derivative on MatrixElement can we make it do the
right thing?

No. Because in the case that fails, MatrixElement._eval_derivative is never called. That's the whole point of the existing _diff_wrt interface: Derivative.__new__ will try to replace the instances of the MatrixElement we differentiate wrt in the expression so we can differentiate wrt an ordinary symbol.

I've now written a patch and created a pull request that adds a _eval_derivative_wrt interface that can be used to address that issue.


I preserved the current behavior of the other classes with _diff_wrt=True, but I think the is room for improvement there as well. For example should the following code really return "x"?

from sympy import *
x
= Symbol('x')
f
= Function('f')
diff
(x*f(x), f(x))

With the new interface Function._eval_derivative_wrt could detect that even after substituting f(x) in x*f(x), the resulting term x*xi_1 still contains members of f(x).free_symbols and thus make diff() return "Derivative(x*f(x), f(x))" instead.

Please comment on the pull request if you disagree with the approach I chose.

Clifford Wolf

unread,
Apr 25, 2014, 6:37:47 AM4/25/14
to sy...@googlegroups.com
On Friday, April 25, 2014 2:49:15 AM UTC+2, Matthew wrote:
I haven't thought about this issue much.  Generally speaking though MatrixElement and derivatives and Matrix Expressions are absolutely in a not so perfect state.  Help here is welcome.

I'll see what I can do.

Imo the next thing to fix is this: type(MatrixSymbol.args[0]) is str. It is the name of the symbol. I think there is even a docstring somewhere that states that all members of .args must be derived from Basic, and it actually breaks a couple of functions if this is not the case. For example the following code throws an AttributeError because str has no .match:

from sympy import *
from sympy.abc import n, m
X
= MatrixSymbol('X', n, m)
X
.find(n)

At the moment MatrixSymbol.is_Symbol=False, which kind of makes sense because it has arguments (the dimensions) that can be arbitrary expressions. But it also means that there is code that tries to clone it by performing X.func(*X.args), so the symbol name must be contained in args somehow.

I would recommend to sympify the symbol name, i.e. effectively turning MatrixSymbol into an operator that turns an ordinary symbol into a matrix symbol. This has the downside that mixing Symbol('X') and MatrixSymbol('X', n, m) can have strange effects, but it is the best solution I can think of.

Any other ideas how to fix this problem?



Matthew Rocklin

unread,
Apr 25, 2014, 11:04:47 AM4/25/14
to sy...@googlegroups.com
Yeah, Aaron has been on me to remove the string from MatrixSymbol.  My opinion is that we should allow non-basics in args.  This will help us avoid recreating all of Python.  I generally lose this argument.

Where I win though is with the claim that Symbol is way too magical and big to be used in Matrix Expressions.  I'm going to claim that Matrix expressions are much simpler in design than the core.  Part of this is the use of new assumptions.  Injecting a Symbol, a huge source of magic, right in the center of matrix expressions makes my skin crawl a tiny bit.  A good compromise is to make a Symbol object that is only a Symbol object, and not a scalar value with old assumptions info attached.


--
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 post to this group, send email to sy...@googlegroups.com.
Visit this group at http://groups.google.com/group/sympy.

Clifford Wolf

unread,
Apr 25, 2014, 12:18:58 PM4/25/14
to sy...@googlegroups.com
On Friday, April 25, 2014 5:04:47 PM UTC+2, Matthew wrote:
Yeah, Aaron has been on me to remove the string from MatrixSymbol.  My opinion is that we should allow non-basics in args.  This will help us avoid recreating all of Python.  I generally lose this argument.

Well, I've only been reading sympy code for a day or so now, but I have to agree with Aaron: I also don't like the ideas of non-basics in args..

Where I win though is with the claim that Symbol is way too magical and big to be used in Matrix Expressions.  I'm going to claim that Matrix expressions are much simpler in design than the core.  Part of this is the use of new assumptions.  Injecting a Symbol, a huge source of magic, right in the center of matrix expressions makes my skin crawl a tiny bit.  A good compromise is to make a Symbol object that is only a Symbol object, and not a scalar value with old assumptions info attached.

You mean something like a MatrixSymbolName class with is_Symbol = True but otherwise just a container for the one string?

Are there any plans or existing infrastructure for storing assumptions on matrices, such as is_positive_definite? When I suggested to store a Symbol in MatrixSymbol.args[0] I also thought about maybe using this Symbol to store assumptions regarding the matrix.. What are your thoughts on that?

Matthew Rocklin

unread,
Apr 25, 2014, 12:33:26 PM4/25/14
to sy...@googlegroups.com
Matrix expressions use new assumptions so assumptions are stored outside of the expression.  There are good reasons for this.  See distinction here http://matthewrocklin.com/blog/work/2013/02/05/Assuming/ .  See example here http://scicomp.stackexchange.com/questions/74/symbolic-software-packages-for-matrix-expressions

Regarding Symbol in MatrixSymbol what we really need is a Basic that just holds a name/string.  The current Symbol object does that plus a whole lot of other things that we don't need.  Probably we just need a String(Basic) class or something similar.


--
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 post to this group, send email to sy...@googlegroups.com.
Visit this group at http://groups.google.com/group/sympy.

Clifford Wolf

unread,
Apr 25, 2014, 1:34:11 PM4/25/14
to sy...@googlegroups.com
On Friday, April 25, 2014 6:33:26 PM UTC+2, Matthew wrote:
Regarding Symbol in MatrixSymbol what we really need is a Basic that just holds a name/string.  The current Symbol object does that plus a whole lot of other things that we don't need.  Probably we just need a String(Basic) class or something similar.

Just thinking out loud: What's about storing the name directly in MatrixSymbol.name and overloading MatrixSymbol.func to [return a closure that] inject[s] the name into a new instance.

This would make sense because afaics the whole idea behind the .args interface is make it easy to implement generic transformations, and I don't see why such transformations should ever want to change the name of a MatrixSymbol. In this sense the dimensions are the only 'arguments' to a MatrixSymbol. (The symbol name is afaics not an 'argument' of Symbol objects either.)

Matthew Rocklin

unread,
Apr 25, 2014, 4:18:43 PM4/25/14
to sy...@googlegroups.com
It is convenient to know that all information about a term is in it's type and it's .args.  This is true for almost all of sympy and it makes it very convenient to interact with from other systems.  Any occasion where this doesn't occur (e.g. Symbol, Integer, AppliedPredicate) ends up being a huge headache when inter-operating with sympy.  Essentially whenever I hook something else into sympy I tell it about type and .args.  Then I go and special case every exception. 

I interoperate with sympy a fair amount from external projects so I tend to resist this this sort of solution pretty strongly.

Matrix Expressions is currently pure in this sense.  It serves as a nice example of a very easy-to-deal-with system.


--
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 post to this group, send email to sy...@googlegroups.com.
Visit this group at http://groups.google.com/group/sympy.

Aaron Meurer

unread,
Apr 25, 2014, 5:19:15 PM4/25/14
to sy...@googlegroups.com
On Fri, Apr 25, 2014 at 5:23 AM, Clifford Wolf <cliffor...@gmail.com> wrote:
> On Friday, April 25, 2014 1:16:19 AM UTC+2, Aaron Meurer wrote:
>>
>> If we override _eval_derivative on MatrixElement can we make it do the
>> right thing?
>
>
> No. Because in the case that fails, MatrixElement._eval_derivative is never
> called. That's the whole point of the existing _diff_wrt interface:
> Derivative.__new__ will try to replace the instances of the MatrixElement we
> differentiate wrt in the expression so we can differentiate wrt an ordinary
> symbol.
>
> I've now written a patch and created a pull request that adds a
> _eval_derivative_wrt interface that can be used to address that issue.
>
> https://github.com/sympy/sympy/pull/7437
>
> I preserved the current behavior of the other classes with _diff_wrt=True,
> but I think the is room for improvement there as well. For example should
> the following code really return "x"?

Yes it should. See the docstring of Derivative.

Aaron Meurer

>
> from sympy import *
> x = Symbol('x')
> f = Function('f')
> diff(x*f(x), f(x))
>
> With the new interface Function._eval_derivative_wrt could detect that even
> after substituting f(x) in x*f(x), the resulting term x*xi_1 still contains
> members of f(x).free_symbols and thus make diff() return "Derivative(x*f(x),
> f(x))" instead.
>
> Please comment on the pull request if you disagree with the approach I
> chose.
>
> --
> 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 post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/sympy/177e7ddd-655e-46bf-a64d-a2aae98d5aae%40googlegroups.com.

Aaron Meurer

unread,
Apr 25, 2014, 5:22:36 PM4/25/14
to sy...@googlegroups.com
Funny you should mention that, because that was actually the idea that
we came up with a while back that was considered the best (well at
least I considered it the best). See
https://groups.google.com/d/msg/sympy/pyzfmCq_thI/UwFLv_RDX2IJ. You're
starting to dive into a controversial topic here, which deals in the
core of the design of SymPy, so fair warning.

The next best solution I think was to create a non-Expr Symbol and use
that, which would take care of Matthew's worries about assumptions and
so on. That should not be that difficult to do, I think.

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 post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/sympy/26b50466-8966-49e0-94d1-075011a6a063%40googlegroups.com.

Matthew Rocklin

unread,
Apr 25, 2014, 5:30:17 PM4/25/14
to sy...@googlegroups.com
Yeah, just to be clear, while I'll argue pretty strongly against injecting Symbol I should say that it's a perfectly rational thing to do.  

Clifford I'm actually pretty impressed that you touched this fairly deep issue.  Have you been developing on SymPy long?


Clifford Wolf

unread,
Apr 25, 2014, 10:10:10 PM4/25/14
to sy...@googlegroups.com
On Friday, April 25, 2014 11:30:17 PM UTC+2, Matthew wrote:
Clifford I'm actually pretty impressed that you touched this fairly deep issue.  Have you been developing on SymPy long?

haha. I actually don't think of myself as developing on SymPy yet. I started _using_ sympy on monday (i.e. 5 days ago) and just tried to get the things working that I never managed to get right with Maxima. I think I filed Issue #7421 only a few hours after starting using sympy. Maybe I was just lucky to try differentiating matrix expression on my first day.. ;)

Anyways: It was quite easy to find my way around in the SymPy source code and figure out why things went wrong in this test case. Awesome!

So I think I might use SymPy a lot in the future, in which case I am very interested in getting Derivatives involving Matrix expressions right and am also willing to invest some time..

Clifford Wolf

unread,
Apr 27, 2014, 8:26:10 AM4/27/14
to sy...@googlegroups.com
On Friday, April 25, 2014 11:19:15 PM UTC+2, Aaron Meurer wrote:
> I preserved the current behavior of the other classes with _diff_wrt=True,
> but I think the is room for improvement there as well. For example should
> the following code really return "x"?

Yes it should. See the docstring of Derivative.

That docstring documents that it is this way, but it does not give any reason why it _should_ be this way. For example, the docstring states that

>>> (2*sqrt(1 - sin(x)**2)).diff(cos(x))
0

but it does not really state why this would be a good idea. Its more like a warning, telling my to avoid such expressions. Now consider the following code for Function._eval_derivative_wrt:

def _eval_derivative_wrt(self, expr, new_name):
    new_self
= Dummy(new_name)
    new_expr
= expr.subs(self, new_self)
   
if new_expr.free_symbols.isdisjoint(self.free_symbols):
       
return (new_expr, new_self)
   
return None

This will detect such cases and not attempt to evaluate a derivative wrt cos(x) if the expression contains something like sin(x):

>>> from sympy import cos, sin, sqrt
>>> from sympy.abc import x
>>> (2*cos(x)).diff(cos(x))
2
>>> (2*sqrt(1 - sin(x)**2)).diff(cos(x))
Derivative(2*sqrt(-sin(x)**2 + 1), cos(x))

Afics the legit use cases for differentiating wrt a function are those where you actually get the right result. And my understanding of the documentation is that atm it is my responsibility as user to avoid the other cases. I think we should shift as much of this responsibility as possible from the user to SymPy.

PS: There may be cases where one would prefer a blind substitution, but I think that this should not be the default behavior of diff(), as it will return wrong results in some cases and therefore is not really a plain derivative anymore.

Aaron Meurer

unread,
Apr 27, 2014, 8:14:48 PM4/27/14
to sy...@googlegroups.com
On Sun, Apr 27, 2014 at 7:26 AM, Clifford Wolf <cliffor...@gmail.com> wrote:
> On Friday, April 25, 2014 11:19:15 PM UTC+2, Aaron Meurer wrote:
>>
>> > I preserved the current behavior of the other classes with
>> > _diff_wrt=True,
>> > but I think the is room for improvement there as well. For example
>> > should
>> > the following code really return "x"?
>>
>> Yes it should. See the docstring of Derivative.
>
>
> That docstring documents that it is this way, but it does not give any
> reason why it _should_ be this way. For example, the docstring states that
>
>>>> (2*sqrt(1 - sin(x)**2)).diff(cos(x))
> 0
>
> but it does not really state why this would be a good idea. Its more like a
> warning, telling my to avoid such expressions. Now consider the following
> code for Function._eval_derivative_wrt:

Oh yes, we need to remove that. See https://github.com/sympy/sympy/issues/7027.

>
> def _eval_derivative_wrt(self, expr, new_name):
> new_self = Dummy(new_name)
> new_expr = expr.subs(self, new_self)
> if new_expr.free_symbols.isdisjoint(self.free_symbols):
> return (new_expr, new_self)
> return None
>
> This will detect such cases and not attempt to evaluate a derivative wrt
> cos(x) if the expression contains something like sin(x):
>
>>>> from sympy import cos, sin, sqrt
>>>> from sympy.abc import x
>>>> (2*cos(x)).diff(cos(x))
> 2
>>>> (2*sqrt(1 - sin(x)**2)).diff(cos(x))
> Derivative(2*sqrt(-sin(x)**2 + 1), cos(x))
>
> Afics the legit use cases for differentiating wrt a function are those where
> you actually get the right result. And my understanding of the documentation
> is that atm it is my responsibility as user to avoid the other cases. I
> think we should shift as much of this responsibility as possible from the
> user to SymPy.
>
> PS: There may be cases where one would prefer a blind substitution, but I
> think that this should not be the default behavior of diff(), as it will
> return wrong results in some cases and therefore is not really a plain
> derivative anymore.

If you want, you can dig up the discussions about this. There were
quite a few, both on GitHub and the mailing list. We finally decided
that the only way to make sense of diff(expr, f(x)) is to think of
expr as a function of x and f(x), like F(x, f(x)), and that the only
way to think of *that* is to think of it as something like F(x,
y)|y=f(x). The primary motivation is Lagrangians, and in that context,
this is what is meant by that notation.

It's a little odd to think of this in "reverse", i.e., instead of
taking a function F(x, y) evaluating it at F(x, f(x)), you take an
expression in x and f(x) and think of it as a function in x and f(x),
but this all makes sense because there is a very precise way to define
this in SymPy, because each expression is a tree, and f(x) appears in
very precise places in that tree, and can be replaced with some other
symbol, and hence differentiated with respect to (with no thought to
other "trees" that are mathematically but not structurally
equivalent).

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 post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/sympy/f15e5572-06b8-452d-9201-00a7f1a660b6%40googlegroups.com.
Reply all
Reply to author
Forward
0 new messages