Sum of partial derivatives of a product over a symbolic number of variables

Question

I would like SymPy to evaluate an expression like the following:

How would I define the symbols and the expression so that SymPy could handle it nicely? I would like to keep N as just a symbol, i.e. not make an actual finite list of x's. I have tried various combinations of IndexedBase and Sum /Product, but didn't get it working right.

Behemoth

It is not difficult to compute that analytically. Why even use SymPy? – Duyne

By the way, the desired result is twice the sum of the x's. – Duyne

@RodrigodeAzevedo The expression I gave would be part a larger one, which I cut out in the interest of keeping the question about SymPy, rather than my specific problem. I have done it all analytically, but I would like to make sure it is correct. – Behemoth

Answer 1

Ideally it would be this:

x = IndexedBase("x")
i, j, N = symbols("i j N")
expr = Sum(Product(exp(-x[j]**2), (j, 1, N)).diff(x[i]), (i, 1, N))

So far this is unevaluated, expr is

Sum(Derivative(Product(exp(-x[j]**2), (j, 1, N)), x[i]), (i, 1, N)) 

The method doit can be used to evaluate it. Unfortunately the differentiation of a product doesn't quite work yet: expr.doit() returns

N*Derivative(Product(exp(-x[j]**2), (j, 1, N)), x[i])

Rewriting the product as the sum prior to differentiation helps:

expr = Sum(Product(exp(-x[j]**2), (j, 1, N)).rewrite(Sum).diff(x[i]), (i, 1, N))
expr.doit()

returns

Sum(Piecewise((-2*exp(Sum(log(exp(-x[j]**2)), (j, 1, N)))*x[i], (1 <= i) & (i <= N)), (0, True)), (i, 1, N))

which is the correct result of differentiation. Sadly we have that extraneous condition in Piecewise, and also log(exp(...)) that should have been simplified. SymPy doesn't infer that (1 <= i) & (i <= N) is True from the context of the outer sum, and it also hesitates to simplify log(exp thinking x[j] might be complex. So I resort to surgical procedure with Piecewise, replacing it by the first piece, and to forcefully expanding logs:

e = expr.doit()
p = next(iter(e.atoms(Piecewise)))
e = expand_log(e.xreplace({p: p.args[0][0]}), force=True)

Now e is

Sum(-2*exp(Sum(-x[j]**2, (j, 1, N)))*x[i], (i, 1, N))

Couldn't get exp(Sum(..)) to become a Product again, unfortunately.

Gnatcatcher

Also, factor(e) recognizes that the double sum can be factored as a product of two sums. – Gnatcatcher