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$ Commented Sep 27, 2015 at 14:46
  • 1
    $\begingroup$ The dagger is \[Dagger] or <escape> dg <escape> $\endgroup$
    – Bob Hanlon
    Commented Sep 27, 2015 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
    Commented Sep 28, 2015 at 11:32
  • $\begingroup$ @PatrickStevens. See the relations here: bosons (not fermions). $\endgroup$
    – march
    Commented Sep 28, 2015 at 17:39
  • $\begingroup$ perhaps you would like to have a look at SNEG package :) nrgljubljana.ijs.si/sneg $\endgroup$
    – Sumit
    Commented Sep 28, 2015 at 17:55

2 Answers 2


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
    Commented Sep 30, 2015 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
    Commented Sep 30, 2015 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
    Commented Oct 5, 2015 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
    Commented Oct 5, 2015 at 4:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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