# Define operator algebra

My objective is to work out commutators like this $$[f(x,y)\partial_x^2+g(x,y)\partial_x+h(x,y),a(x,y)\partial_x^2+b(x,y)\partial_x+c(x,y)]$$ or $$[f(x,y)\partial_x^2+g(x,y)\partial_x+h(x,y),a(x,y)\partial_y^2+b(x,y)\partial_y+c(x,y)]$$ by starting from simpler commutator identities, for example, using $$[\partial_x,f(x,z)]=f_x(x,z)$$ it is easy to derive $$[\partial^2_x,f(x,z)]=2f_x(x,z)\partial_x+f_{xx}(x,z)$$ and so on, until I get to the first two commutators that I want. But how can I do that with Mathematica?

My difficulty is that I don't know how to teach Mathematica to use identities like the third equation to derive the fourth and up to the first two.

• – Jens Sep 15 '13 at 20:07

Define the 2 operators whose commutator you wish to derive. I have ignored the dependence on $y$ to simplify the expressions without losing any essential content.

op1 = (f[x] Dt[#, {x, 2}] + g[x] Dt[#, x] + h[x] #) &;
op2 = (a[x] Dt[#, {x, 2}] + b[x] Dt[#, x] + c[x] #) &;


Apply the commutator to a test function $u(x)$, and then separately collect the $u(x), u^\prime (x), u^{\prime\prime}(x), u^{\prime\prime\prime}(x)$ terms.

op1[op2[u[x]]] - op2[op1[u[x]]] // Collect[#, {u[x], u'[x], u''[x], u'''[x]}] &

(* the result is a bit too long to quote here *)


You can then read off the operator-valued version of the commutator by simply removing the test function $u(x)$ from this result.

Alternatively, a fully operator-valued derivation of the commutator could follow the lines I outline below, which was a quick bit of interactive operator manipulation that I concocted.

Define the 2 operators, using op[...] as a container to hold an "operator product", which allows us to control the order in which the various operators (and c-numbers) occur.

o1 = op[f, dx, dx] + op[g, dx] + op[h];
o2 = op[a, dx, dx] + op[b, dx] + op[c];


Form the commutator.

op[o1, o2] - op[o2, o1]


Distribute the op over Plus.

% //. op[u_ + v_, w_] :> op[u, w] + op[v, w]
% //. op[w_, u_ + v_] :> op[w, u] + op[w, v]


Flatten the nested op.

% /. op[op[u__], op[v__]] :> op[u, v]


Move the dx operators to the right.

% //. op[u___, dx, v_?(# =!= dx &), w___] :> op[u, v, dx, w] + op[u, d[v, x], w]


Move the c-numbers outside the op.

% //. op[u_?(# =!= dx &), v___] :> u op[v]


Tidy up the notation.

% /. d[d[u_, v_], v_] :> d[u, {v, 2}]
% /. op[] -> 1


Collect terms by derivative.

% // Collect[#, op[__]] &


This gives the same result as the quick derivation earlier.

• I believe the Distribute and Flatten steps can be simplified by using Flat and Distribute -- may I edit your answer? If I get it wrong you can of course revert the edit. Also _?(# =!= dx &) could be replaced with Except[dx] I believe. – Mr.Wizard Sep 15 '13 at 20:49
• I see what you mean about your suggestions. Though which approach is to be preferred? - (1) explicit but inefficient patterns/replacements, or (2) core functions which are more efficient. I usually use (1) when I knock code together quickly, so I can be sure about what is going on, and then later I refine this to (2) — so my long-term “library code” looks more like (2) than (1). It would be instructive to mention both approaches, so I welcome your edits/additions. – Stephen Luttrell Sep 16 '13 at 0:28
• I like your argument and I don't think I should replace your code. I also don't want to post a second answer merely changing a couple of things. I'm not sure how I should add things; does a section at the bottom illustrating these alternatives seem appropriate to you? – Mr.Wizard Sep 16 '13 at 1:03
• Because you intend to show how a subset of patterns/replacements can be finessed by using the core functions appropriately, a separate answer would be appropriate, and it would give you more space for your material. Actually, to properly discuss the (1) versus (2) issue (see my previous comment) would warrant a whole separate question/answer. – Stephen Luttrell Sep 16 '13 at 10:56