I have tried to adapt this answer to my problem of calculating some bosonic commutation relations, but there are still some issues.

The way I'm implementing the commutator is straightforward:

commutator[x_, y_] :=x**y-y**x

Example: if I want to compute $[b^\dagger b,ab^\dagger]$ I write

commutator[BosonC["b"]**BosonA["b"], BosonA["a"]**BosonC["b"]]

and the output is $ab^\dagger$ as it should be. However this fails when I compute $[a^\dagger a,ab^\dagger]$ (which should be $-ab^\dagger)$:

commutator[BosonC["a"]**BosonA["a"], BosonA["a"]**BosonC["b"]]
 Out: a†** a^2 ** b†- a ** b†** a†** a

How can I modify the code in this answer to have it work properly?

EDIT Building on the answers of @QuantumDot and @evanb, I came up with this solution. First I implement the commutator, with Distribute.

NCM[x___] = NonCommutativeMultiply[x];
SetAttributes[commutator, HoldAll]

NCM[] := 1
commutator[NCM[x___], NCM[y___]] := Distribute[NCM[x, y] - NCM[y, x]]
commutator[x_, y_] := Distribute[NCM[x, y] - NCM[y, x]]

Then I implement two tools, one for teaching Mathematica how to swap creation and annihilation operators and one is for operator ordering:

dag[x_] := ToExpression[ToString[x] ~~ "†"]

mode[x_] := With[{x†= dag[x]},
NCM[left___, x, x†, right___] := NCM[left, x†, x, right] + NCM[left, right]]

leftright[L_, R_] := With[{R† = dag[R], L† = dag[L]}, 
NCM[left___, pr : R | R†, pl : L | L†, right___] := NCM[left, pl, pr, right]]

Now I can use it like this: after evaluating the definitions I input (for instance)


And finally I can evaluate commutators, for instance

commutator[NCM[a†,a] + NCM[b†,b], NCM[a,b†]]
(* 0 *)
  • $\begingroup$ Can you write out the reasoning behind $[a^{\dagger} a, a b^{\dagger}] = -a b^{\dagger}$? It's not clear what rules you apply to get that, from this question or the linked one. $\endgroup$ – Patrick Stevens Sep 27 '15 at 14:46
  • 1
    $\begingroup$ The dagger is \[Dagger] or <escape> dg <escape> $\endgroup$ – Bob Hanlon Sep 27 '15 at 16:38
  • 1
    $\begingroup$ Can you post the full code that actually produces a ** b\[Dagger] from commutator[BosonC["b"]**BosonA["b"], BosonA["a"]**BosonC["b"]]? $\endgroup$ – Taiki Sep 28 '15 at 11:32
  • $\begingroup$ @PatrickStevens. See the relations here: bosons (not fermions). $\endgroup$ – march Sep 28 '15 at 17:39
  • $\begingroup$ perhaps you would like to have a look at SNEG package :) nrgljubljana.ijs.si/sneg $\endgroup$ – Sumit Sep 28 '15 at 17:55

The function NonCommutativeMultiply has too long of a name, so I make a short-hand version of it (NCP stands for non-commutative product):

NCP[x___] := NonCommutativeMultiply[x];

Now, here's the code

NCP[] := 1
NCP[left___, a, a†, right___] := NCP[left, a†, a, right] + NCP[left, right]
NCP[left___, b, b†, right___] := NCP[left, b†, b, right] + NCP[left, right]
NCP[left___, pl : a | a†, pr : b | b†, right___] := NCP[left, pr, pl, right]

Now your function:

SetAttributes[commutator, HoldAll]    
commutator[NCP[x___], NCP[y___]] := NCP[x, y] - NCP[y, x]

Let's give it a try:

commutator[NCP[b†, b], NCP[a, b†]]

b† ** a

  • $\begingroup$ Great answer too, I'm testing them out. Here I would need linearity to make it perfect. Maybe I can fiddle around with Distribute $\endgroup$ – Ziofil Sep 30 '15 at 15:32

You get get far by just implementing the linearity rules for the commutator (and that 0 and 1 "escape" the matrix product. The first two aren't strictly needed, but I include them for completeness:


NonCommutativeMultiply[H___, 1, T___] := NonCommutativeMultiply[H, T]

Now evaluating commutator[adag**a,a**bdag] gives


You can get further by implementing the fact that different oscillators commute. For example:

commutator[adag, bdag] = 0
commutator[a, bdag] = 0

Reevaluating commutator[adag**a,a**bdag] gives

commutator[adag, a] ** bdag ** a

Finally, you can implement the canonical relation:

commutator[adag, a] = 1;

To get what you want.

You can of course decorate / generalize what's above so that you don't have to write every relation out by hand.

  • $\begingroup$ This is great! How can I tell MMA that the commutator is antisymmetric (so I can avoid writing the rules for, say, both [a,a†] and [a†,a], etc..?) $\endgroup$ – Ziofil Sep 30 '15 at 15:26
  • $\begingroup$ Ah, that's a little bit tricky. You have to make sure you don't create a pattern that always binds, giving an infinite loop. One simple solution is commutator[A_, B_] := -commutator[B, A] /; Not[OrderedQ[{A, B}]]. You can come up with more complicated predicates if you want particular orderings, of course. The above works, for example, if you specify [a,adag]=-1, but it won't deduce that relation from [adag, a]=1. $\endgroup$ – evanb Oct 5 '15 at 3:49
  • $\begingroup$ Another approach is to use Simplify with a custom TransformationFunction. Assuming the definitions in my answer: tryOtherOrder[e_] := (e /. {commutator[A_, B_] :> -commutator[B, A]}) commutator[a, adag] (* doesn't change *) Simplify[commutator[a, adag], TransformationFunctions -> tryOtherOrder] (* gives -1*) $\endgroup$ – evanb Oct 5 '15 at 4:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.