Message Boards Message Boards

0
|
11759 Views
|
27 Replies
|
0 Total Likes
View groups...
Share
Share this post:
GROUPS:

Recursive law to calculate derivatives of Spherical Harmonics

Posted 9 years ago

Hello

Can you please help with the following?

Assume we know the following identity, which yields the derivative of the product of Spherical Harmonic with some function of radial coordinate, f[r].

(\[PartialD]/\[PartialD]z)(f[r] SphericalHarmonicY[l, 
     m, \[Theta], \[CurlyPhi]]) = 
 Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 1))]
    SphericalHarmonicY[l + 1, m, \[Theta], \[CurlyPhi]] (\!\(
\*SubscriptBox[\(\[PartialD]\), \(r\)]\(f[r]\)\) - l f[r]/r) + 
  Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 1))] (\!\(
\*SubscriptBox[\(\[PartialD]\), \(r\)]\(f[r]\)\) + (l + 1) f[r]/r)

The relation tells that the derivative is a linear combination of two SphericalHarmonics with lowered and increased first argument,multiplied by some differential operator which acts on the function f[r].

I would like to take advantage of the recursive definition and understand how to define the derivative operator in such a way that applying it twice (or more) will give a second (or higher) derivative.

Can you please help to write such an operator in a proper way? I have attached short nb file showing my attempt..

Thank you, Shimon

Attachments:
POSTED BY: Shimon Rubin
27 Replies

In the meanwhile you might have figured it out, but nevertheless

(* re-create the identity (A.37): The point is keeping the m invariant, i.e. using identities (A.35) and (A.36)
for D[LegendreP[l,m,Cos[\[CurlyTheta]]],\[CurlyTheta]] together with 
(A.20) and (A.21) for Sin[\[CurlyTheta]] D[LegendreP[l,m,Cos[\[CurlyTheta]]] to do that. 
Did not work out, so use SphericalHarmonicY as an index function as before. *)
Off[SphericalHarmonicY::argrx]
Clear[kp, km, \[Xi], singleQ, \[Zeta]]
km[l_, m_] := Sqrt[(l + m) (l - m)/((2 l + 1) (2 l - 1))]
kp[l_, m_] := km[l + 1, m]
singleQ[G_, f_] := Length[Select[Level[f, {1}], Head[#] == G &]] == 1 /; Head[f] === Times
singleQ[G_, f_] :=  And @@ (singleQ[G, #] & /@ (List @@ f)) /; Head[f] === Plus
singleQ[G_, f_] := False /; Head[f] =!= Times && Head[f] =!= Plus
\[Zeta][G_, y_, f_] := ((D[f,  y] /. {G[l_, m_, x__] ->  kp[l, m] G[l + 1, m, x] + km[l, m] G[l - 1, m, x]}) + 
    (f /. {G[l_, m_, x__] -> -(l/r) kp[l, m] G[l + 1, m, x] + (l + 1)/r km[l, m] G[l - 1, m, x]})) /; ! FreeQ[f, G] 
    && FreeQ[Denominator[Together[f]], G] && singleQ[G, Expand[f]] && FreeQ[D[f, y], Derivative[__][G]]


\[Zeta][SphericalHarmonicY, r, 
  f[r] SphericalHarmonicY[\[Lambda], \[Mu], \[Alpha], \[Beta]]] // TraditionalForm



Simplify[\[Zeta][SphericalHarmonicY, 
   r, \[Zeta][SphericalHarmonicY, r, 
    f[r] SphericalHarmonicY[\[Lambda], \[Mu], \[Alpha], \[Beta]]]]] // TraditionalForm



Simplify[\[Zeta][SphericalHarmonicY, 
   r, \[Zeta][SphericalHarmonicY, 
    r, \[Zeta][SphericalHarmonicY, r, 
     f[r] SphericalHarmonicY[\[Lambda], \[Mu], \[Alpha], \[Beta]]]]]] // TraditionalForm

this last one gives

3 calls

POSTED BY: Udo Krause

Even if the formula I wrote has a typo

It has a typo, indeed. Bethe/Salpeter Q.M. of One- and Two-Electron Atoms, Plenum/Rosetta, 1977,p. 349, formula (A.37) reads

bethe<em>salpetercorrect

do you see the difference? And now the verification works in principle

In[13]:= FullSimplify[
 Together[Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 
            1))] SphericalHarmonicY[l + 1, 
       m, \[Theta], \[CurlyPhi]] (f'[r] - l f[r]/r) + 
       Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 
            1))] SphericalHarmonicY[l - 1, 
       m, \[Theta], \[CurlyPhi]] (f'[
         r] + (l + 1) f[r]/r) - (Cos[\[Theta]] D[#, r] - 
       Sin[\[Theta]]/r D[#, \[Theta]]) &[
   f[r] SphericalHarmonicY[l, m, \[Theta], \[CurlyPhi]]]], 
 Element[l | m, Integers]]

Out[13]= $Aborted

In[15]:= FullSimplify[
 Coefficient[
  Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 
            1))] SphericalHarmonicY[l + 1, 
       m, \[Theta], \[CurlyPhi]] (f'[r] - l f[r]/r) + 
       Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 
            1))] SphericalHarmonicY[l - 1, 
       m, \[Theta], \[CurlyPhi]] (f'[
         r] + (l + 1) f[r]/r) - (Cos[\[Theta]] D[#, r] - 
       Sin[\[Theta]]/r D[#, \[Theta]]) &[
   f[r] SphericalHarmonicY[l, m, \[Theta], \[CurlyPhi]]], f[r]], 
 Element[l | m, Integers]]

Out[15]= (E^(-I \[CurlyPhi]) (Sqrt[Gamma[1 + l - m]] Sqrt[
      Gamma[2 + l + m]]
       Sin[\[Theta]] SphericalHarmonicY[l, 
       1 + m, \[Theta], \[CurlyPhi]] + 
     E^(I \[CurlyPhi]) Sqrt[Gamma[l - m]] Sqrt[
      Gamma[1 + l + 
        m]] ((1 + l) Sqrt[((l - m) (l + m))/(-1 + 4 l^2)]
          SphericalHarmonicY[-1 + l, m, \[Theta], \[CurlyPhi]] + 
        m Cos[\[Theta]] SphericalHarmonicY[l, 
          m, \[Theta], \[CurlyPhi]] - 
        l Sqrt[((1 + l - m) (1 + l + m))/(3 + 4 l (2 + l))]
          SphericalHarmonicY[1 + l, 
          m, \[Theta], \[CurlyPhi]])))/(r Sqrt[Gamma[l - m]] Sqrt[
   Gamma[1 + l + m]])

In[16]:= fsf = %15;

In[24]:= fsf /. {l -> RandomInteger[73], 
  m -> RandomInteger[99], \[Theta] -> 
   RandomReal[{0, \[Pi]}], \[CurlyPhi] -> RandomReal[{0, 2 \[Pi]}]}

Out[24]= (1.18791*10^-14 - 3.56613*10^-15 I)/r

In[17]:= FullSimplify[
 Coefficient[
  Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 
            1))] SphericalHarmonicY[l + 1, 
       m, \[Theta], \[CurlyPhi]] (f'[r] - l f[r]/r) +
       Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 
            1))] SphericalHarmonicY[l - 1, 
       m, \[Theta], \[CurlyPhi]] (f'[
         r] + (l + 1) f[r]/r) - (Cos[\[Theta]] D[#, r] - 
       Sin[\[Theta]]/r D[#, \[Theta]]) &[
   f[r] SphericalHarmonicY[l, m, \[Theta], \[CurlyPhi]]], f'[r]], 
 Element[l | m, Integers]]

Out[17]= Sqrt[((l - m) (l + m))/(-1 + 4 l^2)]
   SphericalHarmonicY[-1 + l, m, \[Theta], \[CurlyPhi]] - 
 Cos[\[Theta]] SphericalHarmonicY[l, m, \[Theta], \[CurlyPhi]] + 
 Sqrt[((1 + l - m) (1 + l + m))/(3 + 4 l (2 + l))]
   SphericalHarmonicY[1 + l, m, \[Theta], \[CurlyPhi]]

In[18]:= fsfs = %17;

In[20]:= fsfs /. {l -> RandomInteger[73], 
  m -> RandomInteger[99], \[Theta] -> 
   RandomReal[{0, \[Pi]}], \[CurlyPhi] -> RandomReal[{0, 2 \[Pi]}]}

Out[20]= 0.

I say in principle because the fsfs term gave in all runs practically zero, but the fsf term might give also Indeterminate depending on the chosen quantum numbers and angles. That needs inspection. At least now there is a valuable candidate for the identity. Even more, there is also a possible candidate to solve p.o.'s original question because every summand on the r.h.s. now contains a SphericalHarmonicY guiding the multiple application of the rule. Go ahead!

POSTED BY: Udo Krause

I say in principle because the fsfs term gave in all runs practically zero, but the fsf term might give also Indeterminate depending on the chosen quantum numbers and angles. That needs inspection.

The inspection consists of course in remembering that $Y_l^m(\vartheta,\varphi)$ is defined for $l=0,1,2,\dots$ and $ -l\leq m \leq l$ and the numerical tests must take that into account

In[363]:= fsf /. Thread[Rule[{l, m, \[Theta], \[CurlyPhi]}, 
   Flatten[{Block[{o = RandomInteger[99]}, {o, RandomInteger[{-o+1, o-1}]}], RandomReal[{0, \[Pi]}], RandomReal[{0, 2 \[Pi]}]}]]]

Out[363]= (1.3764*10^-14 + 1.26875*10^-14 I)/r


In[364]:= fsfs /. Thread[Rule[{l, m, \[Theta], \[CurlyPhi]}, 
   Flatten[{Block[{o = RandomInteger[99]}, {o, RandomInteger[{-o+1, o-1}]}], RandomReal[{0, \[Pi]}], RandomReal[{0, 2 \[Pi]}]}]]]

Out[364]= -2.77556*10^-17 - 3.46945*10^-18 I

to be intelligible. The tests can be done without specifying the angles:

In[446]:= fsf /. Thread[Rule[{l, m}, Block[{o = RandomInteger[199]}, {o, RandomInteger[{-o + 1, o - 1}]}]]] // FullSimplify
Out[446]= 0

In[447]:= fsfs /. Thread[Rule[{l, m}, Block[{o = RandomInteger[199]}, {o, RandomInteger[{-o + 1, o - 1}]}]]] // FullSimplify
Out[447]= 0
POSTED BY: Udo Krause
Posted 9 years ago

Thank you.

Yes. Had a mistake when rewrote the identity from the book into the message. Of course it doesn't matter to the question I asked regarding how to define recursive operators and the answer you provided.

Regarding identities in the book. identity A.37 is correct (proved it). Identities A.38, and A.39 are wrong. To correct these, one has to multiply the right hand side by -1.

Overall these relations are very useful to obtain different relations involving Spherical harmonics (often encountered in Quantum Mechanics). For instance, using these relations, it can be shown that three dimensional Laplacian of any Spherical Harmonic times a function of just the radial coordinates yields Spherical harmonics multiplied by a differential operator which acts on the function of the radial coordinate, explicitly given by

((\[PartialD]^2/\[PartialD]x^2)+(\[PartialD]^2/\[PartialD]y^2)+(\
\[PartialD]^2/\[PartialD]z^2)) (f[r] SphericalHarmonicY[l, 
     m, \[Theta], \[CurlyPhi]]) = 
 SphericalHarmonicY[l, 
   m, \[Theta], \[CurlyPhi]] ((\[PartialD]^2/\[PartialD]r^2)+(2/
      r) (\[PartialD]/\[PartialD]r)-((l (l + 1) )/r^2)) f[r]

A first step in the proof is to recognize that the Laplacian can be rewritten as

(\[PartialD]^2/\[PartialD]x^2)+(\[PartialD]^2/\[PartialD]y^2)+(\
\[PartialD]^2/\[PartialD]z^2) = ((\[PartialD]/\[PartialD]x)+I \
\[PartialD]/\[PartialD]y) ((\[PartialD]/\[PartialD]x)-I \[PartialD]/\
\[PartialD]y) + \[PartialD]^2/\[PartialD]z^2

And the second step is to use the relations above.

POSTED BY: Shimon Rubin

identity taken from the book by Bethe and Salpeter, "Quantum Mechanics of One and Two Electron Systems".

Okay, let's check this with the book; I ordered the 1977 edition from the library.

Even if the formula I wrote has a typo

A typo could be of great benefit to accomplish the task at hand.

POSTED BY: Udo Krause

I'm sorry but I'm still having problems to verify the original relation

po_relation

It doesn't work out in general, so it was tried to work it out per coefficient of f[r],D[f[r],r] where it didn't work out again. Finally the coefficients have been checked numerically, again unsuccessful:

 In[47]:= FullSimplify[Coefficient[
      Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 
                1))] SphericalHarmonicY[l + 1, 
           m, \[Theta], \[Phi]] (f'[r] - l f[r]/r) + 
           Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 1))] (f'[
             r] + (l + 1) f[r]/r) - (Cos[\[Theta]] D[#, r] - 
           Sin[\[Theta]]/r D[#, \[Theta]]) &[f[r] SphericalHarmonicY[l, m, \[Theta], \[Phi]]], f[r]], 
     Element[l | m, Integers]]

    Out[47]= (1/r)((1 + l) Sqrt[((l - m) (l + m))/(-1 + 4 l^2)] + 
      m Cos[\[Theta]] SphericalHarmonicY[l, m, \[Theta], \[Phi]] + (
      E^(-I \[Phi]) Sqrt[Gamma[1 + l - m]] Sqrt[Gamma[2 + l + m]]
        Sin[\[Theta]] SphericalHarmonicY[l, 1 + m, \[Theta], \[Phi]])/(
      Sqrt[Gamma[l - m]] Sqrt[Gamma[1 + l + m]]) - 
      l Sqrt[((1 + l - m) (1 + l + m))/(3 + 4 l (2 + l))]
        SphericalHarmonicY[1 + l, m, \[Theta], \[Phi]])

    In[51]:= fsf = %47;

    In[94]:= fsf /.  {l -> RandomInteger[73], m -> RandomInteger[99], 
\[Theta] -> RandomReal[{0, \[Pi]}], \[Phi] -> RandomReal[{0, 2 \[Pi]}]}

    Out[94]= (22.7394 + 0. I)/r

    In[48]:= FullSimplify[Coefficient[
      Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 
                1))] SphericalHarmonicY[l + 1, 
           m, \[Theta], \[Phi]] (f'[r] - l f[r]/r) + 
           Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 1))] (f'[
             r] + (l + 1) f[r]/r) - (Cos[\[Theta]] D[#, r] - 
           Sin[\[Theta]]/r D[#, \[Theta]]) &[f[r] SphericalHarmonicY[l, m, \[Theta], \[Phi]]], f'[r]], 
     Element[l | m, Integers]]

    Out[48]= Sqrt[((l - m) (l + m))/(-1 + 4 l^2)] - 
     Cos[\[Theta]] SphericalHarmonicY[l, m, \[Theta], \[Phi]] + 
     Sqrt[((1 + l - m) (1 + l + m))/(3 + 4 l (2 + l))]
       SphericalHarmonicY[1 + l, m, \[Theta], \[Phi]]

    In[54]:= fsfs = %48;

    In[85]:= fsfs /. {l -> RandomInteger[73], m -> RandomInteger[99], \[Theta] -> 
       RandomReal[{0, \[Pi]}], \[Phi] -> RandomReal[{0, 2 \[Pi]}]}

    Out[85]= 0.366093 - 0.00046913 I

An identity should turn the coefficients of the linear independent variables f[r], D[f[r],r] into zero.

The intention was to use the knowledge gained in verification of the relation to define an operator applicable in any order.

POSTED BY: Udo Krause
Posted 9 years ago

I still didn't look in the numerical check for the identity, but I am pretty sure the identity is correct. In any case, I will check this identity taken from the book by Bethe and Salpeter, "Quantum Mechanics of One and Two Electron Systems". Even if the formula I wrote has a typo the question how to define an operator that "knows" how to read the indices of the function (e.g. Spherical Harmonic) he acts on is relevant to me.

POSTED BY: Shimon Rubin

This general question is possibly easier to answer. Let's call the operator gop and let be it linear. Then

gop[H[y] G[l, x]] = D[H[y] (l G[l + 1, x] + (l + 1) G[l - 1, x]), y]
= l D[H[y] G[l + 1, x], y] + (l + 1) D[H[y] G[l - 1, x], y]
= l D[H[y], y] G[l + 1, x] + (l + 1) D[H[y], y] G[l - 1, x]

gop[gop[H[y] G[l, x]]] = l gop[D[H[y], y] G[l + 1, x]] + (l + 1) gop[D[H[y], y] G[l - 1, x]]
= l ((l + 1) D[H[y], {y, 2}] G[l + 2, x] + (l + 2) D[H[y], {y, 2}] G[
       l, x]) + (l + 1) ((l - 1) D[H[y], {y, 2}] G[l, x] + 
     l D[H[y], {y, 2}] G[l - 2, x])

In fact gop operates on it's index function G and then does a differentiation. Because the differentiation is also linear, one should make the index function as well as the differentiation variable explicit:

Clear[gop, singleQ]
singleQ[G_, f_] := Length[Select[Level[f, {1}], Head[#] == G &]] == 1 /; Head[f] === Times
singleQ[G_, f_] := And @@ (singleQ[G, #] & /@ (List @@ f)) /; Head[f] === Plus
singleQ[G_, f_] := False /; Head[f] =!= Times && Head[f] =!= Plus
gop[G_, y_, f_] := D[f /. G[l_, x__] -> l G[l + 1, x] + (l + 1) G[l - 1, x], y] /; ! FreeQ[f, G] && 
    FreeQ[Denominator[Together[f]], G] && singleQ[G, Expand[f]] && FreeQ[D[f, y], Derivative[__][G]]

Things one isn't able to pin down are excluded. So one tries to implement the In-doubt-never-Rule. Check against the above double application of gop:

In[302]:= gop[G, y, gop[G, y, H[y] G[l, x]]]

Out[302]= ((1 + l) (l G[-2 + l, x] + (-1 + l) G[l, x]) + 
   l ((2 + l) G[l, x] + (1 + l) G[2 + l, x])) (H^\[Prime]\[Prime])[y]

In[309]:= %302 - (l ((l + 1) D[H[y], {y, 2}] G[l + 2, x] + (l + 2) D[
         H[y], {y, 2}] G[l, x]) + (l + 
       1) ((l - 1) D[H[y], {y, 2}] G[l, x] + 
       l D[H[y], {y, 2}] G[l - 2, x])) // Simplify

Out[309]= 0

Well, some more test cases:

In[289]:= gop[G, z, gop[G, y, H[y] Q[y] G[l, x]]]

Out[289]= 0

In[304]:= gop[G, x, gop[G, y, H[y] Q[y] G[l, x]]]

Out[304]= gop[G, x, ((1 + l) G[-1 + l, x] + l G[1 + l, x]) Q[
    y] Derivative[1][H][y] + ((1 + l) G[-1 + l, x] + l G[1 + l, x]) H[
    y] Derivative[1][Q][y]]

In[305]:= gop[G, x, gop[G, x, gop[G, y, H[y] Q[y] G[l, x]]]]

Out[305]= gop[G, x, 
 gop[G, x, ((1 + l) G[-1 + l, x] + l G[1 + l, x]) Q[y] Derivative[1][
     H][y] + ((1 + l) G[-1 + l, x] + l G[1 + l, x]) H[y] Derivative[
     1][Q][y]]]

In[307]:= gop[G, y, 
 gop[G, y, gop[G, y, H[y] G[l, x] + Q[y] P[y] G[l, x]]]]

Out[307]= 
3 ((1 + l) (l ((-1 + l) G[-3 + l, x] + (-2 + l) G[-1 + l, x]) + (-1 + 
          l) ((1 + l) G[-1 + l, x] + l G[1 + l, x])) + 
    l ((2 + l) ((1 + l) G[-1 + l, x] + l G[1 + l, x]) + (1 + 
          l) ((3 + l) G[1 + l, x] + (2 + l) G[3 + l, x]))) Derivative[
   1][Q][y] (P^\[Prime]\[Prime])[y] + 
 3 ((1 + l) (l ((-1 + l) G[-3 + l, x] + (-2 + l) G[-1 + l, x]) + (-1 +
           l) ((1 + l) G[-1 + l, x] + l G[1 + l, x])) + 
    l ((2 + l) ((1 + l) G[-1 + l, x] + l G[1 + l, x]) + (1 + 
          l) ((3 + l) G[1 + l, x] + (2 + l) G[3 + l, x]))) Derivative[
   1][P][y] (Q^\[Prime]\[Prime])[
   y] + ((1 + 
       l) (l ((-1 + l) G[-3 + l, x] + (-2 + l) G[-1 + l, x]) + (-1 + 
          l) ((1 + l) G[-1 + l, x] + l G[1 + l, x])) + 
    l ((2 + l) ((1 + l) G[-1 + l, x] + l G[1 + l, x]) + (1 + 
          l) ((3 + l) G[1 + l, x] + (2 + l) G[3 + l, x]))) 
\!\(\*SuperscriptBox[\(H\), 
TagBox[
RowBox[{"(", "3", ")"}],
Derivative],
MultilineFunction->None]\)[
   y] + ((1 + 
       l) (l ((-1 + l) G[-3 + l, x] + (-2 + l) G[-1 + l, x]) + (-1 + 
          l) ((1 + l) G[-1 + l, x] + l G[1 + l, x])) + 
    l ((2 + l) ((1 + l) G[-1 + l, x] + l G[1 + l, x]) + (1 + 
          l) ((3 + l) G[1 + l, x] + (2 + l) G[3 + l, x]))) Q[y] 
\!\(\*SuperscriptBox[\(P\), 
TagBox[
RowBox[{"(", "3", ")"}],
Derivative],
MultilineFunction->None]\)[
   y] + ((1 + 
       l) (l ((-1 + l) G[-3 + l, x] + (-2 + l) G[-1 + l, x]) + (-1 + 
          l) ((1 + l) G[-1 + l, x] + l G[1 + l, x])) + 
    l ((2 + l) ((1 + l) G[-1 + l, x] + l G[1 + l, x]) + (1 + 
          l) ((3 + l) G[1 + l, x] + (2 + l) G[3 + l, x]))) P[y] 
\!\(\*SuperscriptBox[\(Q\), 
TagBox[
RowBox[{"(", "3", ")"}],
Derivative],
MultilineFunction->None]\)[y]

Despite the restriction to linear expressions in the index function G and despite the inability to detect the index function as well as the differentiation variable automatically the definitions might be useful. Please check for errors before usage.

Mix the index function a bit:

In[310]:= gop[G, y, gop[G, y, H[y] G[l, x] + Q[y] P[y] G[l, a, b, c]]]

Out[310]= 
2 ((1 + l) (l G[-2 + l, a, b, c] + (-1 + l) G[l, a, b, c]) + 
    l ((2 + l) G[l, a, b, c] + (1 + l) G[2 + l, a, b, c])) Derivative[
   1][P][y] Derivative[1][Q][
   y] + ((1 + l) (l G[-2 + l, x] + (-1 + l) G[l, x]) + 
    l ((2 + l) G[l, x] + (1 + l) G[2 + l, x])) (H^\[Prime]\[Prime])[
   y] + ((1 + l) (l G[-2 + l, a, b, c] + (-1 + l) G[l, a, b, c]) + 
    l ((2 + l) G[l, a, b, c] + (1 + l) G[2 + l, a, b, c])) Q[y] (
   P^\[Prime]\[Prime])[
   y] + ((1 + l) (l G[-2 + l, a, b, c] + (-1 + l) G[l, a, b, c]) + 
    l ((2 + l) G[l, a, b, c] + (1 + l) G[2 + l, a, b, c])) P[y] (
   Q^\[Prime]\[Prime])[y]
POSTED BY: Udo Krause
Posted 9 years ago

Thank you very much for the detailed answered. I am trying to digest all details and focus on the following piece which gives the correct answer when the operator 'gop' is applied twice

Clear[gop, singleQ]
singleQ[G_, f_] := 
 Length[Select[Level[f, {1}], Head[#] == G &]] == 1 /; 
  Head[f] === Times
singleQ[G_, f_] := 
 And @@ (singleQ[G, #] & /@ (List @@ f)) /; Head[f] === Plus
singleQ[G_, f_] := False /; Head[f] =!= Times && Head[f] =!= Plus
gop[G_, y_, f_] := 
 D[f /. G[l_, x__] -> l G[l + 1, x] + (l + 1) G[l - 1, x], 
   y] /; ! FreeQ[f, G] && FreeQ[Denominator[Together[f]], G] && 
   singleQ[G, Expand[f]] && FreeQ[D[f, y], Derivative[__][G]]
gop[G, y, gop[G, y, H[y] G[l, x]]]

Can you please tell what is the role of each one of the three lines above, which start with singleQ[G, f] := ...?

POSTED BY: Shimon Rubin
  • !FreeQ[f, G] do not evaluate gop if the index function G does not appear in the expression f
  • FreeQ[Denominator[Together[f]], G] do not evaluate gop if the index function G does apear in the denominator of the expression f
  • singleQ[G, Expand[f]] do not evaluate gop if the index function G does not appear linear in the expression f
  • FreeQ[D[f, y], Derivative[__][G]] do not evaluate gop if the index function G depends on the differentiation variable y

Usually, the evaluation of something is only meaningful if some conditions apply. Despite or better because of the very far reaching rule evaluation machinery of Mathematica it is often preferable to state what was in one's mind as one wrote an expression.

POSTED BY: Udo Krause
Posted 9 years ago

Thank you very much for your answer. Trying to bring all details together.

POSTED BY: Shimon Rubin
Posted 9 years ago

Thank you very much for the code. Indeed, the function shRubinD you have defined doesn't give the correct result when applied twice. It doesn't take into account that the indices 'l' and 'm' have changed.

My question isn't limited to SphericalHarmonics and can be refined as following. Assume we have the following recursive definition of a symbolic Operator which operates on a function H[y] Subscript[G, l][x]

Operator[H[y] Subscript[G, l][x]] = 
 l D[H[y], y] Subscript[G, l + 1][x] + 
(l + 1) D[H[y], y] Subscript[G, l - 1][x]

Here, 'x' and 'y' are two different arguments while 'l' is an index. How can one "teach" Mathematica to increase and decrease the index of the function G, as well as to take proper combination of the H function and yield a symbolic result due to applying the Operator arbitrary number of times? Particularly, this will enable to obtain the result of Operator^2.

POSTED BY: Shimon Rubin

It doesn't take into account that the indices 'l' and 'm' have changed.

Right, right: A possible SphericalHarmonicY member in f must be detected by shRubinD and treated afterwards.

POSTED BY: Udo Krause

Doesn't this tend to do the job?

In[47]:= Clear[shRubinD]
shRubinD[f_, x_, l_, m_, y_, z_] := 
  Sqrt[((l + m + 1) (l - m + 1))/((2 l + 3) (2 l + 
           1))] SphericalHarmonicY[l + 1, m, y, 
      z] (D[f, x] - l f/x) + 
    Sqrt[((l + m) (l - m))/((2 l + 1) (2 l - 1))] (D[f, 
        x] + (l + 1) f/x) /; 
   AtomQ[l] && AtomQ[m] && AtomQ[x] && AtomQ[y] && 
    AtomQ[z] && ! FreeQ[f, x];
shRubinD[f_, x_, n_Integer, l_, m_, y_, z_] := 
 Nest[shRubinD[#, x, l, m, y, z] &, f, n] /; n > 1

In[50]:= shRubinD[Sin[x], x, p, y, \[CurlyPhi], \[CurlyTheta]]

Out[50]= Sqrt[((p - y) (p + y))/((-1 + 2 p) (1 + 2 p))] (Cos[
     x] + ((1 + p) Sin[x])/x) + 
 Sqrt[((1 + p - y) (1 + p + y))/((1 + 2 p) (3 + 2 p))] (Cos[x] - (
    p Sin[x])/x) SphericalHarmonicY[1 + p, 
   y, \[CurlyPhi], \[CurlyTheta]]

In[52]:= shRubinD[Sin[x], x, 2, p, 
  y, \[CurlyPhi], \[CurlyTheta]] // FullSimplify

Out[52]= (1/(x^2))(((p - y) (p + 
     y) (2 (1 + p) x Cos[x] + (p + p^2 - x^2) Sin[x]))/(-1 + 4 p^2) - 
  Sqrt[((p - y) (p + y))/(-1 + 4 p^2)]
    Sqrt[((1 + p - y) (1 + p + y))/(
   3 + 4 p (2 + p))] (-2 x Cos[x] + (1 + 2 p (1 + p) + 2 x^2) Sin[
       x]) SphericalHarmonicY[1 + p, 
    y, \[CurlyPhi], \[CurlyTheta]] + ((1 + p - y) (1 + p + 
     y) (-2 p x Cos[x] + (p + p^2 - x^2) Sin[x]) SphericalHarmonicY[
    1 + p, y, \[CurlyPhi], \[CurlyTheta]]^2)/(3 + 4 p (2 + p)))

Well, it is extraordinarily ugly because it has SphericalHarmonicY as inmate. But you seem to intend that.

POSTED BY: Udo Krause

Possibly it's not (yet) correct because the second shRubinD does not take on the changed l and m indices of the SphericalHarmonicY from the previous derivative. Please check that.

POSTED BY: Udo Krause
Posted 9 years ago

Thanks. Yet, the problem I like to solve seems to demand other solution.

As far as I see, the main difficulty is that the function needs to recognize the value of the first argument of the Spherical Harmonic function and then modify it accordingly. I don't know how to access the argument of the Spherical Harmonic function (or any other) in case it is multiplied by another function (which is relevant here according to the identity).

POSTED BY: Shimon Rubin

You should be able to do that with pattern matching.

In[3]:= h[f_[x_]*g[x_]] := f[x]*g[x + 1]^2

In[4]:= h[f[x]*g[x]]

Out[4]= f[x] g[1 + x]^2
POSTED BY: Frank Kampas
Posted 9 years ago

Thanks. Can you please share more details? It is hard for me to follow.

Which pattern matching function ti use? How do I access the argument of Spherical Harmonic function which is multiplied by another function?

POSTED BY: Shimon Rubin

Putting an underscore after an expression makes it a named pattern. That's how I accessed the argument of g in the above example.

POSTED BY: Frank Kampas
Posted 9 years ago

Thanks for good intention, nevertheless I couldn't understand any of your answers.

POSTED BY: Shimon Rubin

Since Mathematica already knows the derivatives of the spherical harmonics with respect to its arguments, why don't you use the chain rule to get the derivative with respect to z?

In[1]:= D[SphericalHarmonicY[l, m, \[Theta], phi], \[Theta]]

Out[1]= m Cot[\[Theta]] SphericalHarmonicY[l, m, \[Theta], phi] + (
 E^(-I phi) Sqrt[Gamma[1 + l - m]] Sqrt[Gamma[2 + l + m]]
   SphericalHarmonicY[l, 1 + m, \[Theta], phi])/(
 Sqrt[Gamma[l - m]] Sqrt[Gamma[1 + l + m]])
POSTED BY: Frank Kampas
Posted 9 years ago

Thanks. That was the first thing I tried. However, taking the second derivative of the Spherical Harmonic with respect to the variable 'z' directly (using the chain rule), yields a very long expression which isn't reduced upon using FullSimplify

Simplify[-Sin[\[Theta]]/
  r D[-Sin[\[Theta]]/
    r D[SphericalHarmonicY[l, m, \[Theta], phi], \[Theta]], \[Theta]]]

(E^(-2 I phi) (E^(I phi) Sqrt[
      Gamma[-1 + l - 
        m]] (E^(I phi) m (-1 + (1 + m) Cos[\[Theta]]^2) Sqrt[
         Gamma[l - m]] Sqrt[Gamma[1 + l + m]]
          SphericalHarmonicY[l, m, \[Theta], phi] + (1 + m) Sqrt[
         Gamma[1 + l - m]] Sqrt[Gamma[2 + l + m]]
          Sin[2 \[Theta]] SphericalHarmonicY[l, 1 + m, \[Theta], 
          phi]) + Sqrt[Gamma[l - m]] Sqrt[Gamma[1 + l - m]] Sqrt[
      Gamma[3 + l + m]]
       Sin[\[Theta]]^2 SphericalHarmonicY[l, 2 + m, \[Theta], 
       phi]))/(r^2 Sqrt[Gamma[-1 + l - m]] Sqrt[Gamma[l - m]] Sqrt[
   Gamma[1 + l + m]])

Unfortunately using FullSimplify doesn't use the internal relations between Spherical Harmonics and the compact relation for the derivative I wrote above isn't reproduced. Since I am interested in the simplest symbolic expression and not just on numeric value, it doesn't resolve the issue.

POSTED BY: Shimon Rubin

I'm puzzled by your expression for the derivative with respect to z, calculated with the chain rule. Since z is a function of r and theta, shouldn't there be more than one term?

POSTED BY: Frank Kampas
Posted 9 years ago

Just used the chain rule, yielding that the derivative with respect to the Cartesian 'z' coordinate in Spherical coordinates is given by

 Cos[\[Theta]] (\[PartialD]/\[PartialD]r)-(Sin[\[Theta]]/
   r) \[PartialD]/\[PartialD]\[Theta]

Next if interested in operating on 'r' independent function such as Spherical Harmonic just drop the partial derivative with respect to 'r'.

POSTED BY: Shimon Rubin

You can define a function that takes a function as an argument.

In[27]:= Clear[f, g]

In[28]:= g[f[x_]] := f[x + 1]^2;

In[29]:= g[f[2]]

Out[29]= f[3]^2
POSTED BY: Frank Kampas

Probably better to define your function, rather than redefine the built-in derivative.

POSTED BY: Frank Kampas
Posted 9 years ago

Thanks for the reply. Can you please elaborate what do you mean by define your function?

POSTED BY: Shimon Rubin
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract