docs: def calculate_gamma_matrix(magnetic_states, Omega=1):
      r"""Calculate the matrix of decay between states.
  
      This function calculates the matrix $\gamma_{ij}$ of decay rates between
      states |i> and |j> (in the units specified by the Omega argument).
  
      >>> g=State("Rb",87,5,0,1/Integer(2))
      >>> e=State("Rb",87,5,1,3/Integer(2))
      >>> magnetic_states=make_list_of_states([g,e],"magnetic")
  
      To return the rates in rad/s:
      >>> print calculate_gamma_matrix(magnetic_states)
      [[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -12702506.296014734, -15878132.870018415, -15878132.870018415, -0.0, -19053759.4440221, -9526879.72201105, -3175626.5740036834, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -12702506.296014734, -15878132.870018415, -0.0, -15878132.870018415, -0.0, -9526879.72201105, -12702506.296014734, -9526879.72201105, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -12702506.296014734, -0.0, -15878132.870018415, -15878132.870018415, -0.0, -0.0, -3175626.5740036834, -9526879.72201105, -19053759.4440221, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -3810751.8888044204, -0.0, -0.0, -12702506.296014734, -6351253.148007367, -0.0, -0.0, -0.0, -38107518.8880442, -12702506.296014734, -2540501.2592029474, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -1905375.9444022102, -1905375.9444022102, -0.0, -6351253.148007367, -3175626.5740036834, -9526879.72201105, -0.0, -0.0, -0.0, -25405012.592029467, -20324010.07362358, -7621503.77760884, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -635125.3148007367, -2540501.259202947, -635125.3148007367, -0.0, -9526879.72201105, -0.0, -9526879.72201105, -0.0, -0.0, -0.0, -15243007.55521768, -22864511.33282652, -15243007.55521768, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -0.0, -1905375.9444022102, -1905375.9444022102, -0.0, -0.0, -9526879.72201105, -3175626.5740036834, -6351253.148007367, -0.0, -0.0, -0.0, -7621503.77760884, -20324010.07362358, -25405012.592029467, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -0.0, -0.0, -3810751.8888044204, -0.0, -0.0, -0.0, -6351253.148007367, -12702506.296014734, -0.0, -0.0, -0.0, -0.0, -2540501.2592029474, -12702506.296014734, -38107518.8880442], [12702506.296014734, 12702506.296014734, 12702506.296014734, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [15878132.870018415, 15878132.870018415, 0.0, 3810751.8888044204, 1905375.9444022102, 635125.3148007367, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [15878132.870018415, 0.0, 15878132.870018415, 0.0, 1905375.9444022102, 2540501.259202947, 1905375.9444022102, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 15878132.870018415, 15878132.870018415, 0.0, 0.0, 635125.3148007367, 1905375.9444022102, 3810751.8888044204, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [19053759.4440221, 0.0, 0.0, 12702506.296014734, 6351253.148007367, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [9526879.72201105, 9526879.72201105, 0.0, 6351253.148007367, 3175626.5740036834, 9526879.72201105, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [3175626.5740036834, 12702506.296014734, 3175626.5740036834, 0.0, 9526879.72201105, 0.0, 9526879.72201105, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 9526879.72201105, 9526879.72201105, 0.0, 0.0, 9526879.72201105, 3175626.5740036834, 6351253.148007367, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 19053759.4440221, 0.0, 0.0, 0.0, 6351253.148007367, 12702506.296014734, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 38107518.8880442, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 12702506.296014734, 25405012.592029467, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 2540501.2592029474, 20324010.07362358, 15243007.55521768, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 7621503.77760884, 22864511.33282652, 7621503.77760884, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 15243007.55521768, 20324010.07362358, 2540501.2592029474, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 25405012.592029467, 12702506.296014734, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 38107518.8880442, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]]
  
      To return the rates in 10^6 rad /s:
      >>> gamma = calculate_gamma_matrix(magnetic_states,Omega=1e6)
      >>> gamma
      [[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -12.702506296014734, -15.878132870018417, -15.878132870018417, -0.0, -19.053759444022102, -9.526879722011051, -3.1756265740036835, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -12.702506296014734, -15.878132870018417, -0.0, -15.878132870018417, -0.0, -9.526879722011051, -12.702506296014734, -9.526879722011051, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -12.702506296014734, -0.0, -15.878132870018417, -15.878132870018417, -0.0, -0.0, -3.1756265740036835, -9.526879722011051, -19.053759444022102, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -3.8107518888044205, -0.0, -0.0, -12.702506296014734, -6.351253148007367, -0.0, -0.0, -0.0, -38.107518888044204, -12.702506296014734, -2.5405012592029474, -0.0, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -1.9053759444022103, -1.9053759444022103, -0.0, -6.351253148007367, -3.1756265740036835, -9.526879722011051, -0.0, -0.0, -0.0, -25.40501259202947, -20.32401007362358, -7.62150377760884, -0.0, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -0.6351253148007368, -2.540501259202947, -0.6351253148007368, -0.0, -9.526879722011051, -0.0, -9.526879722011051, -0.0, -0.0, -0.0, -15.24300755521768, -22.86451133282652, -15.24300755521768, -0.0, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -0.0, -1.9053759444022103, -1.9053759444022103, -0.0, -0.0, -9.526879722011051, -3.1756265740036835, -6.351253148007367, -0.0, -0.0, -0.0, -7.62150377760884, -20.32401007362358, -25.40501259202947, -0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.0, -0.0, -0.0, -3.8107518888044205, -0.0, -0.0, -0.0, -6.351253148007367, -12.702506296014734, -0.0, -0.0, -0.0, -0.0, -2.5405012592029474, -12.702506296014734, -38.107518888044204], [12.702506296014734, 12.702506296014734, 12.702506296014734, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [15.878132870018417, 15.878132870018417, 0.0, 3.8107518888044205, 1.9053759444022103, 0.6351253148007368, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [15.878132870018417, 0.0, 15.878132870018417, 0.0, 1.9053759444022103, 2.540501259202947, 1.9053759444022103, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 15.878132870018417, 15.878132870018417, 0.0, 0.0, 0.6351253148007368, 1.9053759444022103, 3.8107518888044205, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [19.053759444022102, 0.0, 0.0, 12.702506296014734, 6.351253148007367, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [9.526879722011051, 9.526879722011051, 0.0, 6.351253148007367, 3.1756265740036835, 9.526879722011051, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [3.1756265740036835, 12.702506296014734, 3.1756265740036835, 0.0, 9.526879722011051, 0.0, 9.526879722011051, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 9.526879722011051, 9.526879722011051, 0.0, 0.0, 9.526879722011051, 3.1756265740036835, 6.351253148007367, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 19.053759444022102, 0.0, 0.0, 0.0, 6.351253148007367, 12.702506296014734, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 38.107518888044204, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 12.702506296014734, 25.40501259202947, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 2.5405012592029474, 20.32401007362358, 15.24300755521768, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 7.62150377760884, 22.86451133282652, 7.62150377760884, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 15.24300755521768, 20.32401007362358, 2.5405012592029474, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 25.40501259202947, 12.702506296014734, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 38.107518888044204, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]]
  
      Let us test if all D2 lines decay at the expected rate (6.065 MHz):
      >>> Gamma =[ sum([ gamma[i][j] for j in range(i)])/2/pi for i in range(len(magnetic_states))][8:]
      >>> for Gammai in Gamma: print Gammai
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
      6.065
  
      """
      Ne = len(magnetic_states)
  
      II = magnetic_states[0].i
  
      gamma = [[0.0 for j in range(Ne)] for i in range(Ne)]
      for i in range(Ne):
          for j in range(i):
              ei = magnetic_states[i]
              ej = magnetic_states[j]
              einsteinAij = Transition(ei, ej).einsteinA
  
              if einsteinAij != 0:
                  ji = ei.j; jj = ej.j
                  fi = ei.f; fj = ej.f
                  mi = ei.m; mj = ej.m
  
                  gammaij = (2.0*ji+1)
                  gammaij *= (2.0*fi+1)
                  gammaij *= (2.0*fj+1)
  
                  gammaij *= float(wigner_6j(ji, fi, II, fj, jj, 1)**2)
                  gammaij *= sum([float(wigner_3j(fj, 1, fi, -mj, q, mi)**2)
                                  for q in [-1, 0, 1]])
  
                  gammaij *= einsteinAij/Omega
                  gammaij = float(gammaij)
  
                  gamma[i][j] = gammaij
                  gamma[j][i] = -gammaij
  
      return gamma