# Simplification of double symbolic sums containing a DiscreteDelta without explicit summation range

I am trying to get Mathematica to automatically do simplifications like the following:

$$\sum\limits_{q}^{q\in qV}\sum\limits_{q'}^{q'\in q'V}{f(q)g(q')\delta(q-q')}=\sum_{q}^{q\in qV}{f(q)g(q)}.$$

Where the range of values that $q$ and $q'$ are denoted by $qV$ and $q'V$ respectively and are the same. With pen an paper it is really easy to do that kind of simplifications because I know that both variables $q$ and $q'$ have the same range and that the Discrete Delta function $\delta(q-q')$ kills one of the sums.

I have solved it, but it seems to me that I haven't done it in the most elegant/efficient way. My solution is:

Clear[qV, qpV, f,g, MyDiscreteDelta]
a = Sum[f[q] Sum[g[qp] MyDiscreteDelta[q, qp], {qp, qpV}], {q, qV}]
b = % /. {qpV -> {q, q1, q2, q3, q4}}
c = % /. {MyDiscreteDelta[x_, y_] -> HoldForm[If[x =!= y, 0, 1]]}
d = ReleaseHold[%]


$a$ defines the sum and I deliberately leave the sum range undefined as even though I know the real range of number I will be summing, it's too large to be put directly and more importantly, at this stage of the calculation I don't care about it, I am just interested in simplifying the result as much as possible, i.e., get rid of one of the sums by using the delta function. I am not using Mathematica defined function DiscreteDelta because a) it takes too long and b) when used symbolically it is left unevaluated.

First thing I need to do to get rid of one of the sums is to give it a range to sum through, in doing so, I am actually giving it a list of possible variables that it could take(assuming that in the general case I will be having more than two summations and I will have more than one delta function in many of the variables). That's done in $b$.

$c$ changes the empty definition of MyDiscreteDelta but I need to hold it as if I don't do so, it will all be evaluated to $0$. The last step, $d$ evaluates the expression and gives the results I want.

This works, but I was wondering if there is a simpler way of doing it.

-
Why not using the built-in KroneckerDelta ? –  b.gatessucks Jan 22 '13 at 12:56
I did a test, changing code MyDiscreetDelta[q,qp] code to code KroneckerDelta[q,qp] code in (a) and only evaluating (a) and (b). It is very slow. It does simplify to 1 when q==qp but on the other cases, it is just left unevaluated, the same as DiscreteDelta. I can individually reproduce the expected behaviour if I write code Assuming[q != qp, Simplify[KroneckerDelta[q, qp]]] code or code Assuming[q != qp, Simplify[DiscreteDelta[q, qp]]] code. But it doesn't work if I do something like code Assuming[q != q1, Simplify[b]]code. –  Paco Jan 22 '13 at 13:57
@Paco, use backticks for code! –  Andrew Jaffe Jan 22 '13 at 16:25

So your goal is to have a function similar to the delta function in that it can deal with purely symbolic variables instead of having attribute NumericFunction as does KroneckerDelta.

Here is a solution that doesn't require holding the expressions during the evaluation. What I do instead is to define MyDiscreteDelta only for situations where it occurs in an explicit sum, i.e., an expression whose Head has been converted from Sum to Plus. At that stage, you'll know that the dummy variable of Sum has been replaced by the actual values, which in your case includes symbolic names drawn from the set qpV, not just numbers.

I changed the definition slightly so that MyDiscreteDelta takes only one argument instead of two, more like the DiracDelta function. That way you can also insert more complicated expression as a condition in MyDiscreteDelta, such as MyDiscreteDelta[q^2 - qp^2] etc. This is also what your initial $\LaTeX$ example looked like.

To make the definitions active in the desired situations, I use TagSetDelayed. This can only be done using another auxiliary function that represents the multiplications in which MyDiscreteDelta might occur - which is a level below the summation:

ClearAll[MyDiscreteDelta, DiscreteDeltaTimes]

MyDiscreteDelta /: Times[any_, MyDiscreteDelta[expr_]] :=
DiscreteDeltaTimes[any, expr]

MyDiscreteDelta /: Plus[any_, MyDiscreteDelta[expr_]] :=
If[expr === 0, any + 1, any]

DiscreteDeltaTimes /: Plus[DiscreteDeltaTimes[any_, expr_], other_] :=
If[expr === 0, any, 0] + other

DiscreteDeltaTimes /:
MakeBoxes[DiscreteDeltaTimes[x_, y_], StandardForm] :=
RowBox[{ToBoxes[x], "\[ThinSpace]", "MyDiscreteDelta", "[",
ToBoxes[y], "]"}]

Clear[qV, qpV, f, g]

a = Sum[f[q] Sum[g[qp] MyDiscreteDelta[q - qp], {qp, qpV}], {q, qV}]


$\displaystyle\sum _q^{\text{qV}}$ f[q] $\displaystyle\sum _{\text{qp}}^{\text{qpV}}$ g[qp] MyDiscreteDelta[q-qp]

a/.{qpV->{q,q1,q2,q3,q4}}


$\displaystyle\sum _q^{\text{qV}}$ f[q] g[q]

In the first output, the MakeBoxes definition above gets applied so that the auxiliary function DiscreteDeltaTimes is hidden from view.

-