Maple User Group Answers

[Anfang] [Hauptseite] [Suchen] [LDFM]

Kronecker delta   (4.11.95)

Holger Friedrich

Is it possible to define a Kronecker delta in a way that would allow Maple to compute infinite sums of expressions involving the delta?

Of course, I can define it as a procedure,

> delta := proc(j, k)
>               if j = k then RETURN(1) else RETURN(0) fi
> end:

but in this way I cannot even compute the simple and obvious sum sum(delta(i, 0), i=0..infinity) = 1 since, it would seem, Maple cannot find a closed-form solution for the sum of an expression that I can't provide to it in closed form in the first place!

> sum('delta(i, 0)', 'i'=0..98);
                                     1

> sum('delta(i, 0)', 'i'=0..99);
Error, (in sum/hgpossib2) division by zero

> sum('delta(i, 0)', 'i'=0..infinity);
Error, (in sum/hgpossib2) division by zero

Such kinds of sums could come in handy when I try to feed some elementary quantum-mechanical problems to Maple. Finite sums suffice for the sort of things I am trying out right now, but maybe infinite sums would be useful for other problems...

Ross Taylor

Yes, I have dealt with this problem in another context (differentiation of arbitrary sums with respect to indexed variables).

Here is an illustration:

> read `c:/maple/sum/tdtest`:

> A:= sum(x[i]*x[j]*delta[i,j],i=1..infinity);

             infinity
              -----
               \
        A :=    )     x[i] x[j] delta[i, j]
               /
              -----
              i = 1

> `simplify/delta`(");

                           2
                       x[j]

The file tdtest can be obtained via anonymous ftp from omnigate.clarkson.edu in the /pub/maple directory.

Let me know if it works (or does not work for you).

Robert Israel

Holger Friedrich (hf@glare.in-chemnitz.de) wrote:

> Of course, things might be a little easier if Maple had a function
> that could look for occurrences of a *function* and return the
> *arguments* of the function call.  I could not get match() to do the
> trick. 

Oh, but it does have one: "indets" with the "specfunc" type.

> S:= sum('delta(n, 0)^2', 'n'=0..infinity):
> indets(S, specfunc(anything,delta));

                                 {delta(n, 0)}

You might also try defining a type:

> `type/deltaproduct`:= proc(x)
       if type(x,`*`) then type([op(x)], list({deltaproduct,deltacoeff }))
    elif type(x,`^`) then type(op(1,x), deltaproduct) and type (op(2,x), posint)
       else  type(x, specfunc(anything,delta))
       fi end;
> `type/deltacoeff`:= x -> not(hastype(x,specfunc(anything,delta)));

Now "type(expr, deltaproduct)" should return true if expr is the sort of thing you want to be able to sum.

I think your "delta" procedure should also return 1 if the two arguments are equal, even if they aren't numeric.

> delta := proc(j, k)
      if j=k then 1     
      elif  type([j, k], [numeric, numeric]) then  0          
      else 'delta'(args)
      fi
      end;

Thanks to all who responded to my inquiry about the Kronecker delta. The reply by Michael Monagan was particularly helpful. If anyone else is interested, here is a raw solution that already works in a few simple cases. Thus, the basic ideas should make sense.

> restart;
> delta := proc(j, k)
>     option `Workaround avoiding sum/hgpossib2 suggested by Michael Monagan`;
>     if type([j, k], [numeric, numeric])
>         then if j = k then 1 else 0 fi
>         else 'delta'(args)
>     fi
> end:
> `sum/delta` := proc()
>     local x, s;
>     if isdelta(args[1]) then
>         x := solve(op(1, args[1]) = op(2, args[1]), args[2]);
>         RETURN(Heaviside(args[2] - x - 1))
>     fi;
>     if isproductwithdelta(args[1], 's') then
>         x := solve(op(1, op(s, args[1])) = op(2, op(s, args[1])), args[2]);
>         RETURN(Heaviside(args[2] - x - 1) * subs(args[2]=x, args[1]))
>     fi;
>     ERROR(`sorry, input expression too complicated :-(`)
> end:
> isdelta := proc(expr)
>     local i;
>     if not has(expr, delta) then RETURN(false) fi;
>     for i to nops(expr) do if has(op(i, expr), delta) then RETURN(false) fi od;
>     true;
> end:
> isproductwithdelta := proc(expr, j)
>     local i;
>     j := `nowhere`;
>     if not type(expr, `*`) then RETURN(false) fi;
>     for i to nops(expr) do if isdelta(op(i, expr)) then j := i; RETURN(true) fi od;
>     false
> end:
> sum('2 * (delta(n, 0) + delta(n, 1))', 'n'=0..infinity);
                                     4

> sum('delta(n, 5) * n^2', 'n'=0..infinity);
                                    25

> sum(expand('(delta(n, 0) + delta(n, 1)) * (delta(n, 2) + delta(n, 3))'), 'n'=0..infinity);
                                     0

> sum('delta(n, 0) * delta(m, 0) * delta(n, m)', 'n'=0..infinity);
                          delta(m, 0) delta(0, m)

> sum('delta(n, 0)^2', 'n'=0..infinity);
Error, (in sum/delta) sorry, input expression too complicated :-(

Holger Friedrich

You were right, of course. This just did not matter when I was investigating situations in which the arguments to the deltas always evaluated to specific integers in the end.

If I want to get the final answers for something like delta(n, n) or delta(n, n+1) without specifying a specific value for n, I have to do what you suggested... and even take it one step further:

delta := proc(j, k)
        if not type(j - k, numeric) then RETURN('procname(args)') fi; 
        if j = k then 1 else 0 fi;        
end:

This version returns zero for delta(n, n+1). After all, there is already enough information to figure out that it is zero.

Thanks again for your suggestion. It was helpful.

[Anfang] [Hauptseite] [Suchen] [LDFM]