For[i=0,i<50,i++, uu[i]= u /. FindRoot[u-Sin[u]/Sqrt[2] == .04 Pi i,{u, .04 Pi i}] ; pp[i]=N[2 ArcTan[ Cos[uu[i]/2], Sqrt[(Sqrt[2]+1)/(Sqrt[2]-1)] Sin[uu[i]/2]]] ] xy=Table[{.5 Cos[pp[i]]/(1+Cos[pp[i]]/Sqrt[2.]),.5 Sin[pp[i]]/(1+Cos[pp[i]]/Sqrt[2.])},{i,0,49}] ListPlot[%,AspectRatio->Automatic,PlotStyle->{RGBColor[1,0,0]}] vxy=Table[{-Sin[pp[i]],Cos[pp[i]]+1/Sqrt[2.]},{i,0,49}] ListPlot[%,AspectRatio->Automatic,PlotStyle->{RGBColor[1,0,0]}] --- e=N[1/Sqrt[2]] phi[u_]=2 ArcTan[ Cos[u/2], Sqrt[(1+e)/(1-e)] Sin[u/2]] r[u_]=(1-e^2)/(1+e Cos[phi[u]]) v[u_]=Sqrt[1+e^2+2 e Cos[phi[u]]] vr[u_]=Sqrt[v[u]^2-(1-e^2)^2/r[u]^2]Sign[Sin[u]] wt[u_]=25(u-e Sin[u])/Pi ParametricPlot[{wt[u],r[u]},{u,0,4 Pi},PlotRange->{0,2},PlotPoints->60] ParametricPlot[{wt[u],v[u]},{u,0,4 Pi},PlotRange->{0,2},PlotPoints->60] ParametricPlot[{wt[u],vr[u]},{u,0,4 Pi},PlotPoints->60] ParametricPlot[{wt[u],r[u]},{u,-uu[5],uu[5]},PlotRange->{.2,.85},PlotPoints->60] ParametricPlot[{wt[u],v[u]},{u,-uu[5],uu[5]},PlotRange->{.8,1.8},PlotPoints->60] ParametricPlot[{wt[u],vr[u]},{u,-uu[5],uu[5]},PlotRange->{-.8,.8},PlotPoints->60] --- Plot[-1/x +.25/x^2,{x,.1,2.5},PlotPoints->100,PlotStyle->{RGBColor[1,0,0]}] Plot[-.5,{x,0,2.5},PlotStyle->{RGBColor[1,0,1]}] Show[%,%%,PlotRange->{-1.25,.25}] Show[%%%,%%,PlotRange->{-1.25,.25}] ---- Plot[{-.5,-.5/4,-.5/9,-.5/16,-.5/25},{x,-.15,.15}] Plot[{-.5/4,-.5/9,-.5/16,-.5/25},{x,.85,1.15}] Plot[{-.5/9,-.5/16,-.5/25},{x,1.85,2.15}] Plot[{-.5/16,-.5/25},{x,2.85,3.15}] Plot[{-.5/25},{x,3.85,4.15}] Show[%%%%%,%%%%,%%%,%%,%,PlotRange->{-.6,.1},Ticks->{None,Automatic}] --- pr[n_,l_,r_]=2/(n^2 Sqrt[Pochhammer[n-l,2l+1]]) (2r/n)^l LaguerreL[n-l-1,2l+1,2r/n] Exp[-r/n] v[l_,r_]=l(l+1)/(2 r^2)-1/r Plot[{pr[1,0,x]-.5,pr[2,0,x]-.125,pr[3,0,x]-.5/9,pr[4,0,x]-.5/16,-.5,-.125,-.5/9,-.5/16, v[0,x]},{x,0,20},PlotRange->{-1,.5}, \ PlotStyle->{{RGBColor[0,0,0]},{RGBColor[0,0,0]},{RGBColor[0,0,0]},{RGBColor[0,0,0]}, {RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,0]}}] Integrate[pr[1,0,r] pr[1,0,r] r^2, {r,0,Infinity}] Integrate[pr[2,0,r] pr[2,0,r] r^2, {r,0,Infinity}] Integrate[pr[3,0,r] pr[3,0,r] r^2, {r,0,Infinity}] Integrate[pr[2,1,r] pr[2,1,r] r^2, {r,0,Infinity}] Integrate[pr[3,1,r] pr[3,1,r] r^2, {r,0,Infinity}] Integrate[pr[3,2,r] pr[3,2,r] r^2, {r,0,Infinity}] Integrate[pr[3,1,r] pr[2,1,r] r^2, {r,0,Infinity}] Plot[{-.5,-.125,-.5/9,-.5/16, v[0,x]},{x,0,20},PlotRange->{-1,.1}, \ PlotStyle->{{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,0]}}] Plot[{-.125,-.5/9,-.5/16, v[1,x]},{x,0,20},PlotRange->{-1,.1}, \ PlotStyle->{{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,1]},{RGBColor[1,0,0]}}] Plot[{-.5/16, v[3,x]},{x,0,20},PlotRange->{-1,.1}, \ PlotStyle->{{RGBColor[1,0,1]},{RGBColor[1,0,0]}}] Plot[pr[1,0,x],{x,0.01,5},PlotRange->{-.1,2}] Plot[pr[4,0,x],{x,0.01,50},PlotRange->{-.04,.25}] Plot[pr[3,0,x],{x,0.01,30},PlotRange->{-.06,.4}] Plot[pr[2,0,x],{x,0.01,15},PlotRange->{-.1,.7}] Plot[pr[4,1,x],{x,0.01,50}] Plot[pr[3,1,x],{x,0.01,30}] Plot[pr[2,1,x],{x,0.01,15}] Plot[pr[4,3,x],{x,0.01,50}] Plot[x pr[30,10,x],{x,0.01,2000}] Plot[pr[30,10,x],{x,0.01,2000},PlotRange->{-.00015,.0003}] Plot[{v[10,x],-.5/30^2},{x,0.01,2000},PlotRange->{-.005,.001}, PlotStyle->{{RGBColor[1,0,0]},{RGBColor[1,0,1]}}] ---- p[n_,l_,m_,r_,theta_,phi_]=2/(n^2 Sqrt[Pochhammer[n-l,2l+1]])(2r/n)^l LaguerreL[n-l-1,2l+1,2r/n] Exp[-r/n] SphericalHarmonicY[l,m,theta,phi] ContourPlot[Abs[p[5,2,2,Sqrt[x^2+z^2],ArcTan[z,x],0]]^2,{x,-60,60},{z,-60,60},PlotPoints->100, ContourLines->False,PlotRange->{0,.00001}] Plot[Abs[p[5,2,2,r,Pi/2,0]]^2,{r,0,60},PlotPoints->60,PlotRange->{0,.00008}] ContourPlot[Abs[p[5,2,2,Sqrt[x^2+y^2],Pi/2,0]]^2,{x,-60,60},{y,-60,60},PlotPoints->100, ContourLines->False,PlotRange->{0,.00001}] ContourPlot[Abs[p[5,2,0,Sqrt[x^2+z^2],ArcTan[z,x],0]]^2,{x,-60,60},{z,-60,60},PlotPoints->100, ContourLines->False,PlotRange->{0,.00004}] Plot[Abs[p[5,2,0,r,0,0]]^2,{r,0,60},PlotPoints->60,PlotRange->{0,.0002}] ContourPlot[Abs[p[5,2,1,Sqrt[x^2+z^2],ArcTan[z,x],0]]^2,{x,-60,60},{z,-60,60},PlotPoints->100, ContourLines->False,PlotRange->{0,.00002}] ContourPlot[Abs[p[4,1,1,Sqrt[x^2+z^2],ArcTan[z,x],0]]^2,{x,-35,35},{z,-35,35},PlotPoints->100, ContourLines->False,PlotRange->{0,.00005}] ContourPlot[Abs[p[4,1,1,Sqrt[x^2+y^2],Pi/2,0]]^2,{x,-35,35},{y,-35,35},PlotPoints->100, ContourLines->False,PlotRange->{0,.00005}] Plot[Abs[p[4,1,1,r,Pi/2,0]]^2,{r,0,35},PlotPoints->60,PlotRange->{0,.0004}] ContourPlot[Abs[p[4,1,0,Sqrt[x^2+z^2],ArcTan[z,x],0]]^2,{x,-35,35},{z,-35,35},PlotPoints->100, ContourLines->False,PlotRange->{0,.0001}] Plot[Abs[p[4,1,0,r,0,0]]^2,{r,0,35},PlotPoints->60,PlotRange->{0,.0008}] ---- pr[n_,l_,r_]=2/(n^2 Sqrt[Pochhammer[n-l,2l+1]]) (2r/n)^l LaguerreL[n-l-1,2l+1,2r/n] Exp[-r/n] Integrate[pr[1,0,r2]^2 r2 Integrate[pr[2,1,r1]^2 r1^2, {r1,0,r2}],{r2,0,Infinity}] =13/729 Integrate[pr[1,0,r2]^2 r2^2 Integrate[pr[2,1,r1]^2 r1, {r1,r2,Infinity}],{r2,0,Infinity}] =164/729 Total= 59/243 = 0.242798 Other way around: Integrate[pr[2,1,r2]^2 r2 Integrate[pr[1,0,r1]^2 r1^2, {r1,0,r2}],{r2,0,Infinity}] Integrate[pr[2,1,r2]^2 r2^2 Integrate[pr[1,0,r1]^2 r1, {r1,r2,Infinity}],{r2,0,Infinity}] -- Integrate[pr[1,0,r2]^2 r2 Integrate[pr[2,0,r1]^2 r1^2, {r1,0,r2}],{r2,0,Infinity}] =7/243 Integrate[pr[1,0,r2]^2 r2^2 Integrate[pr[2,0,r1]^2 r1, {r1,r2,Infinity}],{r2,0,Infinity}] =44/243 Total=17/81=0.209877 Other way around: Integrate[pr[2,0,r2]^2 r2 Integrate[pr[1,0,r1]^2 r1^2, {r1,0,r2}],{r2,0,Infinity}] Integrate[pr[2,0,r2]^2 r2^2 Integrate[pr[1,0,r1]^2 r1, {r1,r2,Infinity}],{r2,0,Infinity}] -- Integrate[pr[1,0,r2]^2 r2 Integrate[pr[1,0,r1]^2 r1^2, {r1,0,r2}],{r2,0,Infinity}] =5/16 Integrate[pr[1,0,r2]^2 r2^2 Integrate[pr[1,0,r1]^2 r1, {r1,r2,Infinity}],{r2,0,Infinity}] =5/16 Total=5/8=0.625 checked --Exchange: Integrate[pr[1,0,r2]pr[2,0,r2] r2 Integrate[pr[2,0,r1]pr[1,0,r1] r1^2, {r1,0,r2}],{r2,0,Infinity}] =8/729 Integrate[pr[1,0,r2]pr[2,0,r2] r2^2 Integrate[pr[2,0,r1]pr[1,0,r1] r1, {r1,r2,Infinity}],{r2,0,Infinity}] = must be same as above but get error: Integrate[pr[2,0,r1]pr[1,0,r1] r1, {r1,r2,Infinity}]= 4 Sqrt[2]/27 - Sqrt[2] r2^2 Exp [-3 r2/2]/3 - 4 Sqrt[2] Gamma[2, 0, 3 r2/2]/27 All except last term: Integrate[pr[1,0,r2]pr[2,0,r2] r2^2 (4 Sqrt[2]/27 - Sqrt[2] r2^2 Exp [-3 r2/2]/3),{r2,0,Infinity}] yields: -8/729 Gamma[2, 0, 3 r2/2] = 1-(1+3 r2/2)Exp[-3 r2/2] Integrate[pr[1,0,r2]pr[2,0,r2] r2^2 4 Sqrt[2] (1-(1+3 r2/2)Exp[-3 r2/2])/27,{r2,0,Infinity}] -16/729 (so with "-", is 16/729 SubTotal=-8/729+16/729=8/729 Total= 16/729= 0.0219479 Integrate[pr[1,0,r2]pr[2,1,r2] Integrate[pr[2,1,r1]pr[1,0,r1] r1^3, {r1,0,r2}],{r2,0,Infinity}] =56/2187 Integrate[pr[1,0,r2]pr[2,1,r2] r2^3 Integrate[pr[2,1,r1]pr[1,0,r1] , {r1,r2,Infinity}],{r2,0,Infinity}] =56/2187 Total= 112/2187 = 0.0512117 -(1+1/4)/2+(0.209877-0.0219479)/2 -(1+1/4)/2+(0.242798+0.0170706)/2 ----- Integrate[pr[2,1,r]pr[2,0,r]r^3,{r,0,Infinity}] ---- p[n_,l_,m_,r_,theta_,phi_]=2/(n^2 Sqrt[Pochhammer[n-l,2l+1]])(2r/n)^l LaguerreL[n-l-1,2l+1,2r/n] Exp[-r/n] SphericalHarmonicY[l,m,theta,phi] sp[x_,z_]=(p[2,0,0,Sqrt[x^2+z^2],ArcTan[z,x],0]+p[2,1,0,Sqrt[x^2+z^2],ArcTan[z,x],0])/Sqrt[2] sp2[x_,z_]=(2 - Sqrt[x^2 + z^2 ] + z)/(8 Sqrt[Pi] Exp[Sqrt[x^2 + z^2 ]/2]) sp1[x_,z_]=(2 - Sqrt[x^2 + z^2 ] - z)/(8 Sqrt[Pi] Exp[Sqrt[x^2 + z^2 ]/2]) ContourPlot[Abs[sp2[x,z]]^2,{x,-10,10},{z,-10,10},PlotPoints->100, ContourLines->False,PlotRange->{0,.01}] ContourPlot[Abs[sp1[x,z]]^2,{x,-10,10},{z,-10,10},PlotPoints->100, ContourLines->False,PlotRange->{0,.01}] Plot[Abs[sp2[0,r]]^2,{r,-10,10},PlotPoints->60,PlotRange->{0,.02}] Plot[Abs[sp1[0,r]]^2,{r,-10,10},PlotPoints->60,PlotRange->{0,.02}] --- Plot[{-30,-15,0,15,30},{x,-.15,.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{-45/2,-15/2,15/2,45/2},{x,.85,1.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{-15,0,15},{x,1.85,2.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{-15/2,15/2},{x,2.85,3.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{0},{x,3.85,4.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{-45/2,-15/2,15/2,45/2},{x,-.85,-1.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{-15,0,15},{x,-1.85,-2.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{-15/2,15/2},{x,-2.85,-3.15},PlotStyle->{RGBColor[1,0,0]}] Plot[{0},{x,-3.85,-4.15},PlotStyle->{RGBColor[1,0,0]}] Show[%%%%%%%%%,%%%%%%%%,%%%%%%%,%%%%%%,%%%%%,%%%%,%%%,%%,%] --- ephi[x_,y_,z_]=1/Sqrt[(x-1)^2+y^2+z^2]+1/Sqrt[(x+1)^2+y^2+z^2]+ 1/Sqrt[(y-1)^2+x^2+z^2]+1/Sqrt[(y+1)^2+x^2+z^2]+ 1/Sqrt[(z-1)^2+y^2+x^2]+1/Sqrt[(z+1)^2+y^2+x^2] Series[ephi[x,y,z],{x,0,6},{y,0,6},{z,0,6}] o4xyz[x_,y_,z_]=7(z^4+y^4+x^4)/2 -(21/2)(z^2 y^2 +y^2 x^2+z^2x^2) o6xyz[x_,y_,z_]=3(z^6+y^6+x^6)/4-45(z^4 y^2+ z^2y^4+z^4x^2+y^4x^2+z^2x^4+y^2x^4)/8 + 135(z^2y^2x^2)/2 o42[theta_,phi_]=o4xyz[Sin[theta] Cos[phi],Sin[theta] Sin[phi],Cos[theta]] o62[theta_,phi_]=o6xyz[Sin[theta] Cos[phi],Sin[theta] Sin[phi],Cos[theta]] Series[ephi[r Sin[theta] Cos[phi],r Sin[theta] Sin[phi], r Cos[theta]],{r,0,6}] Simplify[%] o4[theta_,phi_]=(7 (18 + 30 Cos[4 phi] - 20 Cos[4 phi - 2 theta] + 5 Cos[4 (phi - theta)] + 40 Cos[2 theta] + 70 Cos[4 theta] + 5 Cos[4 (phi + theta)] - 20 Cos[4 phi + 2 theta]))/256 Integrate[Sin[theta]Integrate[o42[theta,phi] SphericalHarmonicY[4,4,theta,phi] ,{phi,0,2 Pi}],{theta,0,Pi}] 4,4=35 Sqrt[Pi/70]/3 4,3=0 4,2=0 4,1=0 4,0=7 Sqrt[Pi]/3 4,-1=0 4,-2=0 4,-3=0 4,-4=35 Sqrt[Pi/70]/3 Simplify[(35 Sqrt[Pi/70]/3)(SphericalHarmonicY[4,4,theta,phi]+SphericalHarmonicY[4,-4,theta,phi])+ (7 Sqrt[Pi]/3) SphericalHarmonicY[4,0,theta,phi]-o4[theta,phi]] o6[theta_,phi_]=(3 (100 - 420 Cos[4 phi] - 231 Cos[4 phi - 6 theta] - 105 Cos[4 phi - 2 theta] + 546 Cos[4 (phi - theta)] + 210 Cos[2 theta] + 252 Cos[4 theta] + 462 Cos[6 theta] + 546 Cos[4 (phi + theta)] - 105 Cos[4 phi + 2 theta] - 231 Cos[4 phi + 6 theta]))/4096 Integrate[Sin[theta]Integrate[o6[theta,phi] SphericalHarmonicY[6,6,theta,phi] ,{phi,0,2 Pi}],{theta,0,Pi}] 6,6=0 6,5=0 6,4=-21 Sqrt[31 Pi/14]/26 6,3=0 6,2=0 6,1=0 6,0=3 Sqrt[13 Pi]/26 Simplify[(-21 Sqrt[13 Pi/14]/26)(SphericalHarmonicY[6,4,theta,phi]+SphericalHarmonicY[6,-4,theta,phi])+ (3 Sqrt[13 Pi]/26) SphericalHarmonicY[6,0,theta,phi]-o6[theta,phi]] Series[ephi[Sin[theta] Cos[phi]/r,Sin[theta] Sin[phi]/r, Cos[theta]/r],{r,0,8}] Simplify[%] phismall[x_,y_,z_]=6+7(z^4+y^4+x^4)/2 -(21/2)(z^2 y^2 +y^2 x^2+z^2x^2)+ 3(z^6+y^6+x^6)/4-45(z^4 y^2+ z^2y^4+z^4x^2+y^4x^2+z^2x^4+y^2x^4)/8 + 135(z^2y^2x^2)/2 phibig[x_,y_,z_]=(6+(7(z^4+y^4+x^4)/2 -(21/2)(z^2 y^2 +y^2 x^2+z^2x^2))/(x^2+y^2+z^2)^(4)+ (3(z^6+y^6+x^6)/4-45(z^4 y^2+ z^2y^4+z^4x^2+y^4x^2+z^2x^4+y^2x^4)/8 + 135(z^2y^2x^2)/2)/(x^2+y^2+z^2)^(6))/Sqrt[x^2+y^2+z^2] phiapprox[x_,y_,z_]:=If[x^2+y^2+z^2<1,phismall[x,y,z],phibig[x,y,z]] Plot[{phiapprox[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]],ephi[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]]},{x,0,2},PlotRange->{0,6}] Plot[{phiapprox[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]],ephi[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]]},{x,0,2},PlotRange->{0,6}, PlotStyle->{{RGBColor[1,0,0]},{RGBColor[0,0,0]}}] Plot[{phiapprox[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]],ephi[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]]},{x,.8,1.2}, PlotStyle->{{RGBColor[1,0,0]},{RGBColor[0,0,0]}}] Plot[ephi[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]],{x,0,2},PlotRange->{0,6}] Plot[{phiapprox[x/Sqrt[2],x/Sqrt[2],0],ephi[x/Sqrt[2],x/Sqrt[2],0]},{x,0,2},PlotRange->{0,6}, PlotStyle->{{RGBColor[1,0,0]},{RGBColor[0,0,0]}}] Plot[{phiapprox[x,0,0],ephi[x,0,0]},{x,0,2},PlotRange->{0,12}, PlotStyle->{{RGBColor[1,0,0]},{RGBColor[0,0,0]}}] ContourPlot[ephi[x,y,0],{x,-2,2},{y,-2,2},PlotPoints->100, Contours->25,PlotRange->{0,9},ContourShading->False] ContourPlot[ephi[x/Sqrt[2],x/Sqrt[2],z],{x,-2,2},{z,-2,2},PlotPoints->100, Contours->25,PlotRange->{0,9},ContourShading->False] ContourPlot[ephi[x/Sqrt[2],x/Sqrt[2],z],{x,-2,2},{z,-2,2},PlotPoints->100, Contours->{9,8.5,8,7.5,7,6.5,6.05,6,5.95,5.5,5,4.5,4,3.5,3,2.5,2},PlotRange->{0,9},ContourShading->False] ContourPlot[ephi[x,y,0],{x,-2,2},{y,-2,2},PlotPoints->100, Contours->{9,8.5,8,7.5,7,6.5,6.05,6,5.95,5.5,5,4.5,4,3.5,3,2.5,2},PlotRange->{0,9},ContourShading->False] ContourPlot[phiapprox[x,y,0],{x,-2,2},{y,-2,2},PlotPoints->100, Contours->{9,8.5,8,7.5,7,6.5,6.05,6,5.95,5.5,5,4.5,4,3.5,3,2.5,2},PlotRange->{0,9},ContourShading->False] ContourPlot[phiapprox[x/Sqrt[2],x/Sqrt[2],z],{x,-2,2},{z,-2,2},PlotPoints->100, Contours->{9,8.5,8,7.5,7,6.5,6.05,6,5.95,5.5,5,4.5,4,3.5,3,2.5,2},PlotRange->{0,9},ContourShading->False] ---- Integrate[Sin[theta]Integrate[SphericalHarmonicY[2,2,theta,phi]SphericalHarmonicY[2,2,theta,phi] \ (SphericalHarmonicY[4,-4,theta,phi]+SphericalHarmonicY[4,4,theta,phi]),{phi,0,2 Pi}],{theta,0,Pi}] int2[l1_,m1_,l2_,m2_,l3_,m3_]:=(-1)^(-m1)Sqrt[(2l1+1)(2l2+1)(2l3+1)/(4 Pi)]\ ThreeJSymbol[{l1,0},{l2,0},{l3,0}]ThreeJSymbol[{l1,-m1},{l2,m2},{l3,m3}] int2[2,-2,4,4,2,2]+int2[2,-2,4,-4,2,2] Sqrt[5/(14 Pi)] int2[2,2,4,0,2,2]=1/(14 Sqrt[Pi]) int2[2,1,4,0,2,1]=-2/(7 Sqrt[Pi]) int2[2,0,4,0,2,0] = 3/(7 Sqrt[Pi]) int2[2,-1,4,0,2,-1] = -2/(7 Sqrt[Pi]) int2[2,-2,4,0,2,-2] = 1/(14 Sqrt[Pi]) m40 = {{1/(14 Sqrt[Pi]),0,0,0,0},\ {0,-2/(7 Sqrt[Pi]),0,0,0},\ {0,0,3/(7 Sqrt[Pi]),0,0},\ {0,0,0,-2/(7 Sqrt[Pi]),0},\ {0,0,0,0,1/(14 Sqrt[Pi])}} m44={{0,0,0,0,Sqrt[5/(14 Pi)]},{0,0,0,0,0},{0,0,0,0,0},{0,0,0,0,0},{Sqrt[5/(14 Pi)],0,0,0,0}} m=(7Sqrt[Pi]/3)m40+Sqrt[35 Pi/18]m44 Eigensystem[m] ---- p[n_,l_,m_,r_,theta_,phi_]=2/(n^2 Sqrt[Pochhammer[n-l,2l+1]])(2r/n)^l LaguerreL[n-l-1,2l+1,2r/n] Exp[-r/n] SphericalHarmonicY[l,m,theta,phi] dxy[x_,y_,z_]=(p[3,2,2,Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]]- p[3,2,-2,Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]])/(Sqrt[2] I) dx2y2[x_,y_,z_]=(p[3,2,2,Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]]+ p[3,2,-2,Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]])/(Sqrt[2]) dxz[x_,y_,z_]=(p[3,2,1,Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]]+ p[3,2,-1,Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]])/(Sqrt[2]) N[dxy[4,4,0]] Plot[Abs[dxy[x/Sqrt[2],x/Sqrt[2],0]]^2,{x,0,15}] ContourPlot[Abs[dxy[x/Sqrt[2],x/Sqrt[2],z]]^2,{x,-16,16},{z,-16,16},PlotPoints->100, ContourLines->False,PlotRange->{0,.0006}] ContourPlot[Abs[dxy[x,y,0]]^2,{x,-16,16},{y,-16,16},PlotPoints->100, ContourLines->False,PlotRange->{0,.0006}] ContourPlot[Abs[dx2y2[x,y,0]]^2,{x,-16,16},{y,-16,16},PlotPoints->100, ContourLines->False,PlotRange->{0,.0006}] ContourPlot[Abs[dxz[x,0,z]]^2,{x,-16,16},{z,-16,16},PlotPoints->100, ContourLines->False,PlotRange->{0,.0006}] ----- p[n_,l_,m_,r_,theta_,phi_]=2/(n^2 Sqrt[Pochhammer[n-l,2l+1]])(2r/n)^l LaguerreL[n-l-1,2l+1,2r/n] Exp[-r/n] SphericalHarmonicY[l,m,theta,phi] px[r_,theta_,phi_]=(-p[2,1,1,r,theta,phi]+p[2,1,-1,r,theta,phi])/Sqrt[2] py[r_,theta_,phi_]=(p[2,1,1,r,theta,phi]+p[2,1,-1,r,theta,phi]) I/Sqrt[2] pz[r_,theta_,phi_]=p[2,1,0,r,theta,phi] s[r_,theta_,phi_]=p[2,0,0,r,theta,phi] teta[r_,theta_,phi_]=(-s[r,theta,phi]+px[r,theta,phi]+py[r,theta,phi]+pz[r,theta,phi])/2 teta2[x_,y_,z_]=Simplify[teta[Sqrt[x^2+y^2+z^2],ArcTan[z,Sqrt[x^2+y^2]],ArcTan[x,y]]] Plot[teta[r,.9553,Pi/4],{r,0,5}] Plot[teta[r,Pi-.9553,Pi/4+Pi],{r,0,5}] ContourPlot[Abs[teta2[x/Sqrt[2],x/Sqrt[2],z]]^2,{x,-8,8},{z,-8,8},PlotPoints->100, ContourLines->False,PlotRange->{0,.006}] Plot[teta2[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]],{x,-10,10}] Plot[Abs[teta2[x/Sqrt[3],x/Sqrt[3],x/Sqrt[3]]]^2,{x,-10,10}] ---- LaguerreL[5,9,x] -- Integrate[Exp[-2 a r]r^2,{r,0,Infinity}]=1/(4a^3) Integrate[Exp[-2 a r](-a^2+(-2+2a)/r)r^2,{r,0,Infinity}]=(-2+a)/(4a^2) Integrate[Exp[-2 a r2]r2 Integrate[Exp[-2 a r]r^2,{r,0,r2}],{r2,0,Infinity}]+ Integrate[Exp[-2 a r2]r2^2 Integrate[Exp[-2 a r]r,{r,r2,Infinity}],{r2,0,Infinity}]=5/(128a^5) a(a-2)+5a/(z 8) /. a->(1-5/(16 z)) ---- ephi[x_,y_,z_]=1/Sqrt[(x-1)^2+(y-1)^2+(z-1)^2]+1/Sqrt[(x+1)^2+(y+1)^2+(z-1)^2]+ 1/Sqrt[(y-1)^2+(x+1)^2+(z+1)^2]+1/Sqrt[(y+1)^2+(x-1)^2+(z+1)^2] Series[ephi[r Sin[theta] Cos[phi],r Sin[theta] Sin[phi], r Cos[theta]],{r,0,6}] Simplify[%] SeriesData[r, 0, {4/3^(1/2), 0, 0, (10*Cos[theta]*Sin[2*phi]*Sin[theta]^2)/(9*3^(1/2)), (7*(-18 - 30*Cos[4*phi] + 20*Cos[4*phi - 2*theta] - 5*Cos[4*(phi - theta)] - 40*Cos[2*theta] - 70*Cos[4*theta] - 5*Cos[4*(phi + theta)] + 20*Cos[4*phi + 2*theta]))/(5184*3^(1/2)), 0, (100 - 420*Cos[4*phi] - 231*Cos[4*phi - 6*theta] - 105*Cos[4*phi - 2*theta] + 546*Cos[4*(phi - theta)] + 210*Cos[2*theta] + 252*Cos[4*theta] + 462*Cos[6*theta] + 546*Cos[4*(phi + theta)] - 105*Cos[4*phi + 2*theta] - 231*Cos[4*phi + 6*theta])/(31104*3^(1/2))}, 0, 7, 1] Note: O3 term has no effect o4[theta_,phi_]=(7 (18 + 30 Cos[4 phi] - 20 Cos[4 phi - 2 theta] + 5 Cos[4 (phi - theta)] + 40 Cos[2 theta] + 70 Cos[4 theta] + 5 Cos[4 (phi + theta)] - 20 Cos[4 phi + 2 theta]))/(5184*3^(1/2)) Integrate[Sin[theta]Integrate[o4[theta,phi] SphericalHarmonicY[4,4,theta,phi] ,{phi,0,2 Pi}],{theta,0,Pi}] 4,4=(4*((35*Pi)/6)^(1/2))/243 4,3=0 4,2=0 4,1=0 4,0=(28*(Pi/3)^(1/2))/243 4,-1=0 4,-2=0 4,-3=0 4,-4=(4*((35*Pi)/6)^(1/2))/243