sphericalharmonics {cooltools}

Spherical Harmonics

Description

Evaluates spherical harmonics Y, either in the real-valued or complex-valued basis.

Usage

sphericalharmonics(l, m, x, basis = "real")

Arguments

l

degree of the spherical harmonic, accurate to about l=500; (0=monopole, 1=dipole, 2=quadrupole, 3=octupole, 4=hexadecapole,...)

m

order of the spherical harmonic (-l,-l+1,...,+l)

x

either an n-by-2 matrix specifying the polar angle theta (0...pi) and azimuthal angle phi (0...2*pi); or an n-by-3 matrix specifying the 3D coordinates of n vectors (whose normalization is irrelevant).

basis

a string specifying the type of spherical harmonics; this has to be either "complex" for the standard complex-valued harmonics with Condon-Shortley phase convention, or "real" (default) for the standard real-valued harmonics.

Value

Returns an n-vector of the spherical harmonics; for points x=c(0,0,0), a value of 0 is returned

Author(s)

Danail Obreschkow

Examples

## Check orthonormalization of all spherical harmonics up to 3rd degree

# make indices l and m up to 3rd degree
l = c(0,rep(1,3),rep(2,5),rep(3,7))
m = c(0,seq(-1,1),seq(-2,2),seq(-3,3))

# check orthonormalization for all pairs
for (i in seq(16)) {
  for (j in seq(16)) {

    # compute scalar product
    f = function(theta,phi) {
      Yi = sphericalharmonics(l[i],m[i],cbind(theta,phi))
      Yj = sphericalharmonics(l[j],m[j],cbind(theta,phi))
      return(Re(Yi*Conj(Yj))*sin(theta))
    }
    g = Vectorize(function(phi) integrate(f,0,pi,phi)$value)
    scalar.product = integrate(g,0,2*pi)$value

    # compare scalar product to expected value
    ok = abs(scalar.product-(i==j))<1e-6
    cat(sprintf('(l=%1d,m=%+1d|l=%1d,m=%+1d)=%5.3f  %s\n',l[i],m[i],l[j],m[j],
                scalar.product+1e-10,ifelse(ok,'ok','wrong')))
  }
}