Product involving Kronecker Deltas

Question

I am trying to do some products involving some objects made out of Kronecker deltas. For example, taking an object like $x_{abcd}=\delta_{ab}\delta_{cd}$, where all the indices run from $1$ to $N$, I would like to compute a product like
$$x_{abmn}x_{mncd},$$ where repeated indices are being summed over. Is there a way in Mathematica to fully compute such a product? I mean, it is quite easy to just get the a symbolic product like $$x_{abmn}x_{mncd}=\delta_{ab}\delta_{mn}\delta_{mn}\delta_{cd}.$$ But I would like to further carry out the sums over $m$ and $n$ to obtain $$x_{abmn}x_{mncd}=\delta_{ab}\delta_{cd}N.$$ Of course this is a very simple and quick example, but eventually I would calculate products involving hundreds of terms so any tips would be appreciated!

Edit: Here is an example code of just the symbolic product (What I have manage to do). There is not much to it, it is just a product, and it produces the right result.

x1=KroneckerDelta[a,c]*KroneckerDelta[m,n];
x2=KroneckerDelta[m,n]*KroneckerDelta[c,d];
x1*x2

The rest I do not know how to do. What I want the code to do is the following $$x_{abmn}x_{mncd}=\delta_{ab}\delta_{mn}\delta_{mn}\delta_{cd}=\sum_{n=1}^N\sum_{m=1}^N\delta_{ab}\delta_{mn}\delta_{mn}\delta_{cd}=\sum_{m=1}^N\delta_{ab}\delta_{mm}\delta_{cd}=\delta_{ab}\delta_{cd}N$$

Answer 1

Here's one way, using pattern matching. (I used nsize instead of N, since N is a reserved symbol). First define x, as you did:

x[a_, b_, c_, d_] := KroneckerDelta[a, b] KroneckerDelta[c, d];

Next define "contract using Einstein Summation" - i.e. sum on repeated indices. There are three rules, one for two different deltas with at least one repeated index, a second for powers, since we get delta^2 terms, and one to terminate the recursion:

contractES[
   Times[x___, KroneckerDelta[a_, b_], y___, KroneckerDelta[c_, b_], 
    z___]] := nsize contractES[Times[KroneckerDelta[a, c], x, y, z]];

contractES[
  Times[x___, KroneckerDelta[a_, b_]^n_?(# > 1 &), y___]] := 
 nsize^2 contractES[Times[x, y, KroneckerDelta[a, b]^(n - 2)]];

contractES[x___] := x;

KroneckerDelta has the Orderless Attribute, so these patterns suffice even if the indices are in a different order in the input.

Now runs tests, first for x:

x[a, b, m, n] x[m, n, c, d]

Out[7]= KroneckerDelta[a, b] KroneckerDelta[c, d] KroneckerDelta[m, 
  n]^2

And now for your example:

contractES[x[a, b, m, n] x[m, n, c, d]]

Out[8]= nsize^2 KroneckerDelta[a, b] KroneckerDelta[c, d]

I get nsize^2 because we get one factor of nsize for each index, and we contract over both m and n. So my nsize^2 is basically what you're calling "N".

Finally do a couple more cases, to show that the order of the indices doesn't matter:

contractES /@ {x[a, b, a, b] x[c, d, a, b], 
  x[m, n, a, b] x[c, d, m, n]}

Out[9]= {nsize^2 KroneckerDelta[a, b] KroneckerDelta[c, d], 
 nsize^2 KroneckerDelta[a, b] KroneckerDelta[c, d]}

And one last one with more x's:

contractES[x[a, b, m, n] x[m, n, c, d] x[c, d, g, h]]

Out[10]= nsize^4 KroneckerDelta[a, b] KroneckerDelta[g, h]

EDIT: Request for additional functionality.

As requested, here's code to (crudely) cover the case where you want differently dimensioned deltas. The basic idea is to build objects with constructors for both a group of products, all containing objects of the same size, and then call another constructor for each symbolic index, using TagSet to store the dimension common to this group. I also added code to allow for mixed integer/symbolic indices, and validity checking, of sorts.

You don't have to create all the "groups" at once. You do have to make all the objects in a given group the same size, since I pass in a single nMax parameter for each call. Assign the output from the constructor to either a symbol or an indexed symbol, as in the examples.

One important limitation: make sure your index names are unique. I didn't add any protection in the index constructor code, so if you repeat a symbol name, you'll overwrite the size[index] value stored, which will affect objects you create earlier. I suppose it would be possible to check to see if an index has already undergone processing, throw an error message and abort, but I'll leave that to you if you want it.

Another less critical limitation: I can't do validity checking on the indices when you either supply a symbolic size or a symbolic index. I..e if you define nMax = 4 and use 'a' for your index, or choose nMax = m and use '8' as an index, I have no way of knowing if it's in range or not, so that's for you to ensure.

Predicates to test for valid indices:

kdIndex::invindx = "Index `1` is invalid for size `2`";

validIndexQ[nMax_Integer][i_Integer] := IntervalMemberQ[Interval[{1,nMax}],i];

validIndexQ[_][_Symbol] = True;

validIndexQ[n_][i_] := (Message[kdIndex::invindx,i,n];False);

Index and group constructors:

makeKDIndex[nMax_][a_Symbol] /; validIndexQ[nMax][a] := size[a]^=nMax;

makeKroneckerIndexGroup[nMax_][texp:Times[expr__]] :=
    Module[{indices={}},
        indices = Union@Flatten@Cases[expr, KroneckerDelta[a_,b_]:>{a,b},Infinity];
        makeKDIndex[nMax][#]&/@indices;
        texp
    ];

Transformation rules for contraction. I added two new rules to cover the case that we collapse down to a single delta raised to a power with no other terms, and hence no multiplication. The other "Times" based rules don't match for that case.

contractES[Times[x___,KroneckerDelta[a_,b_Symbol],y___,KroneckerDelta[c_,b_Symbol],z___]] := 
    size[b] contractES[Times[KroneckerDelta[a,c],x,y,z]];

contractES[Times[x___,KroneckerDelta[a_Integer,b_Symbol]^n_?(#>1&),y___]] := 
    size[b] contractES[Times[x,y,KroneckerDelta[a,b]^(n-2)]];

contractES[Times[x___,KroneckerDelta[a_Symbol,b_Symbol]^n_?(#>1&),y___]] := 
    size[a]size[b] contractES[Times[x,y,KroneckerDelta[a,b]^(n-2)]];

contractES[KroneckerDelta[a_Integer,b_Symbol]^n_?(#>1&)] :=
    size[b] contractES[KroneckerDelta[a,b]^(n-2)];

contractES[KroneckerDelta[a_Symbol,b_Symbol]^n_?(#>1&)] := 
    size[a] size[b] contractES[KroneckerDelta[a,b]^(n-2)];

contractES[x___] := x;

And some tests:

Construct some groups. Assign to some symbol or indexed symbol.

x[1] = makeKroneckerIndexGroup[4][
  KroneckerDelta[a, b] KroneckerDelta[c, b]]
Out[13]= KroneckerDelta[a, b] KroneckerDelta[b, c]

x[2] = makeKroneckerIndexGroup[5][KroneckerDelta[3,d]KroneckerDelta[e,3]KroneckerDelta[3,e]]
Out[14]= KroneckerDelta[3, d] KroneckerDelta[3, e]^2

x[3] = makeKroneckerIndexGroup[6][
  KroneckerDelta[q, r] KroneckerDelta[r, q]]
Out[16]= KroneckerDelta[q, r]^2

x[4] = makeKroneckerIndexGroup[mMax][
  KroneckerDelta[s, t] KroneckerDelta[u, s]]
Out[17]= KroneckerDelta[s, t] KroneckerDelta[s, u]

Check the sizes:

size /@ {a, b, c, d, e, q, r, s, t, u}
Out[23]= {4, 4, 4, 5, 5, 6, 6, mMax, mMax, mMax}

Some contractions:

contractES /@ {x[1] x[2], x[2] x[3], x[3], x[4], x[4] x[1] x[4]}
Out[24]= {20 KroneckerDelta[3, d] KroneckerDelta[a, c], 
 180 KroneckerDelta[3, d], 36, mMax KroneckerDelta[t, u], 
 4 mMax^4 KroneckerDelta[a, c]}

You can mix numbers (and other symbols) in with the deltas, as in the earlier answer.

contractES[x[1] 4.2 x[2] 5.03 x[3]]
Out[20]= 15210.7 KroneckerDelta[3, d] KroneckerDelta[a, c]

contractES[x[3] x[3] x[1] x[2] x[1] x[2]]
Out[21]= 41472000

Answer 2

Define a function and $N$ (recall that N is reserved/already defined)

x[a_,b_,c_,d_]:=KroneckerDelta[a,b]KroneckerDelta[c,d]
nmax = 10;

Then we can write

\!\(\*UnderoverscriptBox[\(\[Sum]\), \(n = 1\), \(nmax\)]
  \(\*UnderoverscriptBox[\(\[Sum]\), \(m = 1\), \(nmax\)]
  \((x[a,b,m,n] x[m,n,c,d])\)\)\)

I copy-pasted that after typing escsumesc, but an equivalent input is

Sum[x[a,b,m,n] x[m,n,c,d], {m, nmax}, {n, nmax}]

Answer 3

Its standard job of replacement regardless of any algebra

  { Sum[x___  * KroneckerDelta[n_,p_]  KroneckerDelta[n_,q_] ,
d1___, {n_ , ___},  d2___ ]  :>  Sum[x  KroneckerDelta[p,q] ,d1,d2] , ...} 

plus same rule for all permutations of n,p,q in the arguments of the deltas