How to write a differential operator in Mathematica

Question

I have a very basic question. I am trying to code the LHS of the differential equation: $$\Big[\frac{1}{\sin\theta}\frac{d}{d\theta}\big(\sin\theta\frac{d}{d\theta}\big)-\frac{m^2}{\sin^2\theta}+n(n+1)\Big]f(\theta)=0$$ into Mathematica. As you can see, in the first term, the differential operator is nested and is to properly act on the function $f(\theta)$.

I tried to use #1, but (lol) I don't know how to use it!

This is what I have:

Sin[θ]^(-1) D[Sin[θ] D[#1, θ], θ] - m^2*#1/Sin[θ]^2 + n (n + 1) #1

What next?!

Also: I don't want to just stick in f[θ] into where I have #1, because eventually I would like to be able to make change of variables!

Answer 1

Why not just define a Function ?

ClearAll[m, n, x, foo]
oper = Function[{f, \[Theta]}, (1/
      Sin[\[Theta]] D[Sin[\[Theta]] D[f[\[Theta]], \[Theta]]] - 
     m^2/Sin[\[Theta]]^2 f[\[Theta]] + n (n + 1) f[\[Theta]])];
oper[foo, x]
foo[x_] := Sin[x]*Cos[x]*x^4
oper[foo, x]

Answer 2

Well, I hope you're not pulling my leg here (at 1000 rep one would think you have met pure (anonymous) functions before, or have at least learned where to get help: type & and press F1). But here it goes:

(*Assign the operator in the form of an anonymous function  to oper*)
oper = Sin[\[Theta]]^(-1) D[Sin[\[Theta]] D[#1, \[Theta]], \[Theta]] - 
   m^2*#1/Sin[\[Theta]]^2 + n (n + 1) #1 &;

(*use it*)
oper [f[\[Theta]]]

Answer 3

The differential operator in this question is itself indexed by two variables m and n. This arises after expressing the Laplace operator in spherical coordinates (see the answer by b.gatessucks, but that answer was incomplete because it dropped the $\phi$ dependence too early) and then assuming that the function on which we operate is proportional to $\exp(im\phi)$.

In the other answers that define a function for this operator, the parameters m and n appear as global variables inside the function. That's inconvenient in practice because these indices correspond to physical parameters (the angular momentum components) that are adjustable. So m and n need to appear in the list of variables on which the operator depends.

The next thing on which the operator depends is the independent variable $\theta$ with respect to which we want to take derivatives. It would be reasonable to make this work along the lines of D, which means in particular that it should handle a case like this:

foo[x_] := Sin[θ] Sin[x]
D[foo[θ], θ]

2 Sin[θ] Cos[θ]

This will not work if you define an operator function like this:

myOp = Function[{f, θ}, D[f[θ], θ]];
myOp[foo[θ], θ]

$\sin ^2(\theta )'(\theta )+(2 \sin (\theta )\cos (\theta ))(\theta )$

The problem is that myOp expects to be passed a function instead of an expression. This is not how D works. The purpose of the second argument of D[f,x] is to identify the variable appearing in f that is to be varied. The above method in myOp is basically what Nasser used in his answer. It performs a differentiation with respect to the first slot in the first argument assuming it is a function, and then return the result with the variable $\theta$ inserted in that slot. That's not wrong, it's just an interpretation of differentiation that differs from what D does.

An approach more consistent with D is this definition:

newOp = Function[{f, θ}, D[f, θ]];
newOp[foo[θ], θ]

2 Cos[θ] Sin[θ]

But in this question, the variable θ plays a role that is better grouped together with the indices m and n, and separated from the function f on which the operator acts. Now that I identified how I want to think about the variables on which the differential operator depends, the definition I would use is as follows:

Clear[operNM, n, m, var]
operNM[n_, m_, var_] = Function[{f}, D[Sin[var] D[f, var], var]/Sin[var] - 
    (m^2 f)/Sin[var]^2 + n (n + 1) f
  ];

This defines operNM[n, m, var] as an operator to which you then supply an expression that it operates on:

operNM[n,m,θ][foo[θ]]

$\csc (\theta ) \left(4 \sin (\theta ) \cos^2(\theta )-2 \sin ^3(\theta )\right)-m^2+n(n+1) \sin ^2(\theta )$

Having the operator act on expressions instead of expecting its argument to be a function is more practical especially when you deal with functions of several variables. In quantum mechanics, the name of the coordinate in a many-particle wave function is tied to the index of the particle it represents, and function slots are then not a convenient way to keep track of that information. In operNM you can choose the desired differentiation variable independently of where in the expression it appears.

Answer 4

If you are on version 9 you could use the built-in functions (previously available in the VectorAnalysis package). The following might get you close to what you need :

Laplacian[f[theta], {r, theta, phi}, "Spherical"] // Simplify
(* (Cot[theta] Derivative[1][f][theta] + (f^\[Prime]\[Prime])[theta])/r^2 *)