Replacing values of a function

Question

I have recently starting using Mathematica and have recently come to what seems as an impasse in my understanding of the language. If this is too "tell me how to do it" I would certainly understand.

I have written a program that applies commutator logic to a sequence of creation and annihilation operators to produce a string of Kronecker deltas. I then need to apply the Kronecker delta restrictions to my operators. This is purely symbolic algebra, no real numbers need be produced at the end of the day.

For example,

(* KD[a,b] is a Kronecker delta with indices a and b. *)
(* KD Attributes: Orderless,NumericFunction *)
(* A[{x},{y}] is an operator with creation indices x and annihilation indices y*)

In:
A[{x},{y}] (KD[x,a1]*KD[y,a2]-KD[x,b1]*KD[y,b2])

Out:
A[{a1},{a2}]-A[{b1},{b2}]

What I have so far:

In:
rules=torules[(KD[x,a1]*KD[y,a2]]
A[{x},{y}]/.rules

Out:
{x->a1,y->a2}
A[{a1},{a2}]

This, however, is not quite what I need and converting from KD's to the rules list seems overly complicated. The main issues with this is that I cannot generate operator expressions that have plus within them. The structure of KD is not easily changed, but the input of the operator expressions is quite flexible.

Another take on this, I am not sure if it will help, but its my current direction.

In:
SetAttributes[A, NumericFunciton]
A[{x},{y}] (KD[x,a1]*KD[y,a2]-KD[x,b1]*KD[y,b2])

Out:
A[{x},{y}]*KD[x,a1]*KD[y,a2]- A[{x},{y}]*KD[x,b1]*KD[y,b2]

This should be a fairly simple replacement pattern, but I cannot get any further.

Answer 1

Here are several rules that can help to implement the Einstein convention:

deltaSimplifyRules =
  {
   (a_ + b_)*x_ :> (a*x + b*x),
   KD[i_Symbol , j_]*
      x_ /; (MemberQ[Attributes[i], Temporary] && 
       MemberQ[x, i, Infinity]) :> (x /. i :> j),
   KD[i_ , j_Symbol]*
      x_ /; (MemberQ[Attributes[j], Temporary] && 
       MemberQ[x, j, Infinity]) :> (x /. j :> i),
   KD[i_Symbol , i_] :>  3 /; MemberQ[Attributes[i], Temporary],
   KD[i_Symbol , j_]^2 /; MemberQ[Attributes[i], Temporary] :> KD[j, j],
   KD[x__] :> Sort[KD[x]]
   };

simplifyKD[x_] := x //. deltaSimplifyRules

now you can use simplifyKD to simplify your expression:

simplifyKD[
 Module[{x, y}, 
  A[{x}, {y}] (KD[x, a1]*KD[y, a2] - KD[x, b1]*KD[y, b2])  ]]

A[{a1}, {a2}] - A[{b1}, {b2}]

Module is used to create a temporary variable. Any temporary variable is summed over by these rules.

Answer 2

Assuming you want to normal-order products of creation and annihilation operators, followed by tidying up the Kronecker deltas and the ordering of indices, you could do something like this (using an obvious notation):

Clear[opprod, del];

Define the basic commutation rule for moving annihilation operators to the right.

opprod[u___, an[i_], cr[j_], v___] := opprod[u, cr[j], an[i], v] + del[i, j] opprod[u, v];

Tidy up a corner case.

opprod[] = id;

Define rules for simplifying and tidying Kronecker deltas - there are other possibilities that could be used here.

del[i_, i_] = 1;
del[i_, j_] /; ! OrderedQ[{i, j}] := del[j, i];
del /: del[i_, j_]^_ := del[i, j];
del /: del[i_, j_] del[j_, k_] := del[i, j] del[i, k];
del /: del[i_, k_] del[j_, k_] := del[i, j] del[i, k];

Define a rule for standardising the indices between Kronecker deltas and operator products.

opprod /: del[i_, j_] opprod[u___, (x : (an | cr))[j_], v___] :=
  del[i, j] opprod[u, x[i], v];

Define a rule for tidying up indices in operator products - this puts the indices of commuting operators in a canonical order.

opprod[u___, (x : (an | cr))[i_], (x : (an | cr))[j_], v___] /; !OrderedQ[{i, j}] :=
  opprod[u, x[j], x[i], v];

Now try out these definitions on a non-trivial operator product.

opprod[an[k], an[j], cr[i], an[i], cr[j], cr[k]] // Expand

(* id del[i, j] + id del[i, k] + 2 id del[i, j] del[i, k] + 
opprod[cr[i], an[i]] + 2 del[i, j] opprod[cr[i], an[i]] + 
2 del[i, k] opprod[cr[i], an[i]] + 
6 del[i, j] del[i, k] opprod[cr[i], an[i]] + 
del[j, k] opprod[cr[i], an[i]] + del[i, k] opprod[cr[j], an[j]] + 
del[i, j] opprod[cr[k], an[k]] + 
opprod[cr[i], cr[j], an[i], an[j]] + 
2 del[i, k] opprod[cr[i], cr[j], an[i], an[j]] + 
2 del[j, k] opprod[cr[i], cr[j], an[i], an[j]] + 
opprod[cr[i], cr[k], an[i], an[k]] + 
2 del[i, j] opprod[cr[i], cr[k], an[i], an[k]] + 
opprod[cr[i], cr[j], cr[k], an[i], an[j], an[k]] *)

I have quoted this result, not because it is necessarily interesting, but so you can verify that the above definitions are working as I had intended. I have cross-checked that this result is the same as what I get using another home-grown method that I have "up my sleeve". If it turns out to be incorrect then please tell me!

It is also very useful to use the Notation package to format the output in a human-readable way. For instance, you could try this (it looks much less messy when you enter it using the Notation package's palette):

<< Notation`

Notation[ParsedBoxWrapper[
  SubscriptBox["C", 
  "i_"]]\[DoubleLongLeftRightArrow]ParsedBoxWrapper[
  RowBox[{"cr", "[", "i_", "]"}]]]
Notation[ParsedBoxWrapper[
  SubscriptBox["A", 
  "i_"]]\[DoubleLongLeftRightArrow]ParsedBoxWrapper[
  RowBox[{"an", "[", "i_", "]"}]]]
InfixNotation[ParsedBoxWrapper["\[SmallCircle]"], opprod]
Notation[ParsedBoxWrapper[
  SubscriptBox["\[Delta]", 
  RowBox[{"i_", ",", 
  "j_"}]]]\[DoubleLongLeftRightArrow]ParsedBoxWrapper[
  RowBox[{"del", "[", RowBox[{"i_", ",", "j_"}], "]"}]]]