Plug-and-go Clebsch-Gordan computation in python?

Question

I started a little project in python, under the assumption that it would be easy to find a routine for numerically computing Clebsch-Gordan coefficients in some library such as scipy. When it came time to do it, I actually found that it was difficult to locate an appropriate library. Scipy doesn't seem to have it. I'm looking for a solution that is open-source, doesn't add a lot of lines of code to my own code, doesn't require me to write my own glue code for a C library, and is available in a debian library -- in other words, a really easy plug-and-go solution. There do seem to be bindings for GSL, but none of them seem to be packaged for debian.

The best solution I was able to come up with was this:

import numpy,sympy
from sympy.physics.quantum.cg import CG  

class Memoize:
  # https://stackoverflow.com/a/1988826/1142217
  def __init__(self, f):
    self.f = f
    self.memo = {}
  def __call__(self, *args):
    if not args in self.memo:
      self.memo[args] = self.f(*args)
    return self.memo[args]

def clebsch(j1,m1,j2,m2,j3,m3):
  """
  Computes <j1 m1 j2 m2 | j3 m3>, where all spins are given as double their values.
  """
  # https://docs.sympy.org/latest/modules/physics/quantum/cg.html
  # https://mattpap.github.io/scipy-2011-tutorial/html/numerics.html
  return CG(sympy.S(j1)/2,sympy.S(m1)/2,sympy.S(j2)/2,sympy.S(m2)/2,sympy.S(j3)/2,sympy.S(m3)/2).doit().evalf()

clebsch = Memoize(clebsch)

This uses the sympy symbolic math library to calculate the CG coefficient symbolically, then evaluates the resulting expression numerically. The performance is not actually too horrible due to memoization, but there is an annoying ~1 second hit on startup time because I have to load sympy. It also just seems goofy to use symbolic math for a numerical application.

I found some people's native python implementations online, but nothing that seemed like a high-quality implementation, and using one of them would have required cutting and pasting their code into mine, bloating the number of lines of code in what I had hoped would be a small, elegant program.

Any suggestions?

[EDIT] I opened a feature request in SciPy, and discussion has ensued about whether there is an implementation that can be copied into or adopted for their use, with an appropriate license.

Answer 1

As you mentioned, there doesn't seem to be a "reference implementation" in SciPy or a SciKit. They have been implemented in the Gnu Scientific Library, so you could call the GSL library from Python, either through PyGSL or by using ctypes. There also is a Fortran95 implementation with Python bindings in the SHTools (Spherical Harmonics) Tools package.

But if you believe this is lacking in SciPy, why not implement this yourself, either looking at GSL or SHTools and the references there in or, for example, based on this paper by Johansson and Forssen? You could share it afterwards with the SciPy community.