Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Suppose I have a vector in $\mathbb{R}^n$ but $n$ is not known in advance. I want to be able to write functions which operate on the components of that vector, and then I'd like to be able to take derivatives with respect to the components. As an example, consider the relation: $$\frac{\partial}{\partial x_j} \sum_i x_i$$

Under the assumption that the $x_i$ are independent, I want a call to Simplify[] to return 1. Similarly, calling Simplify[] on $\frac{\partial x_i}{\partial x_j}$ should give KroneckerDelta[i,j]. It's not clear how I should represent generic coordinates like this. I've seen this, but I'm not sure it provides an answer. As the linked post suggests, I could do this for a fixed $n$, but that's not situation I'm working on, especially since I want to see the generic form for any $n$.

For reference, it seems that sympy let's you do something close to this.

from sympy.tensor import IndexedBase, Idx
x = IndexBase('x')
i, j = map(Idx, ['i', 'j'])
x[i]

x[i].diff
# <bound method Indexed.diff of x[i]>

x[i].diff(x[j])
# ValueError: Can't differentiate wrt the variable: x[j], 1

But despite being able to represent the variables abstractly, I can't seem to differentiate them.

Here is another example. Suppose you wanted to calculate the derivative of the entropy with respect to one of the components of the input distribution (again, assuming all the variables are independent). The final form is the same no matter what $n$-simplex the distribution lives on, so you'd like to be able to do this for any dimension.

$$ \frac{\partial H}{\partial p_j} = - \frac{\partial}{\partial p_j} \sum_i p_i \log p_i = - (\log p_j + 1)$$

Problems like this come up in optimization, when you need to provide the gradient and Hessian to numerical algorithms.

Update: Here are two other related posts:

how to differentiate formally?

How to customize derivative behavior via upvalues?

share|improve this question
    
I don't think this is possible with the built-in symbolic tensor manipulation functions. You'd have to look for packages. It will probably help if you complement your post by describing a realistic use case to give a better idea of what you are looking for. –  Szabolcs Feb 7 at 22:06

2 Answers 2

up vote 3 down vote accepted

I think this can be hacked more or less case by case with UpValues, I think this is one of the most flexible aspects of Mathematica.

For instance if you just want partial derivatives to interact with sums you can just define your sum function MySum (or you can maybe unprotect Sum, not sure if this is possible) and define UpValues

MySum /: D[MySum[s_, i_], x[j_]] := D[s /. i -> j, x[j]]

This gives the desired results for

D[MySum[x[i], i], x[j]]

(1)

-D[MySum[x[i] Log[x[i]], i], x[j]]

(-1 - Log[x[j]])

The x here can also be modified to a more general pattern

share|improve this answer
1  
This is great, but its not quite as flexible as I had hoped. Suppose the summation had a nested summation (e.g., s_ has a MySum inside it) that involved x[k]. This hack wouldn't work as desired then. –  Tom Feb 18 at 6:51
    
With Nested Sums it does indeed become more difficult, it becomes a problem of splitting your sums where some number of indices are the same, I think it can still be done. I am gonna have a think about this and hopefully post something in the next few days –  Eduardo Serna Feb 19 at 19:05
    
Be sure to check out the related posts I added to the question. I was able to get something working with them. –  Tom Feb 19 at 22:31
    
Yeah I was actually doing that now, it seems that the answer there is the right one if you want to work with deltas and not sums at least. This problem was useful to something I was doing anyway –  Eduardo Serna Feb 19 at 22:58
    
Putting the UpValue on x might help. –  Daniel Mahler Sep 18 at 18:30

You can also try unprotecting and changing D.

For example:

D[var_[is___], var_[js___]] := Times @@ MapThread[KroneckerDelta, {{is}, {js}}]

Allows

D[x[i, a], x[j, b]]
(* KroneckerDelta[a, b] KroneckerDelta[i, j] *)

This is by no means a complete solution, since it doesn't implement the product rule, for example:

D[x[i, a] x[j, b], x[k, c]]
(* 0 *)

but I think it's a decent start.

share|improve this answer
    
I don't think this is the way to go because of the number of case by case redefinitions of D you have to is much larger than if you focus on redefining the interaction between the D and MySum –  Eduardo Serna Feb 19 at 19:08

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.