Get the eigenstate from PairStateInteractions.diagonalise

[yshen295]

Hi,
I have a problem about extracting the eigenstate after turning on sortEigenvectors=True, there will be an error ValueError: operands could not be broadcast together with shapes (432,250) (0,) . This is a toy model I build to find the eigenstate, and here is the code. Is there anything wrong with the library, or I was just not on the right track? Thanks!

def arc_eigenstate(a1,n1,a2,n2):
    if a1 == "Cs":
        atom1 = Caesium()
    elif a1 == "Rb":
        atom1 = Rubidium87()
    # n1 = 61
    l1 = 1
    j1 = 1.5
    m1 = 1.5
    if a2 == "Cs":
        atom2 = Caesium()
    elif a2 == "Rb":
        atom2 = Rubidium87()
    # n2 = 61
    l2 = 1
    j2 = 1.5
    m2 = 1.5
    

    Bz = 0.0001# #Tesla
    
    rmin = 1.5
    rmax = 6
    
    theta = 0 # Polar Angle [0-pi]
    phi = 0 # Azimuthal Angle [0-2pi]
    dn = 4 # Range of n to consider (n0-dn:n0+dn)
    dl = 4 # Range of l values
    deltaMax = 25e9  # Max pair-state energy difference [Hz]
    
    nEig = 250  # Number of eigenstates to extract
    
    
    calc = PairStateInteractions(atom1, n1, l1, j1, n2, l2, j2, m1, m2, s=0.5, 
                                 s2=0.5, atom2=atom2)
    
    r = np.linspace(rmin, rmax, 500)

    calc.defineBasis(theta, phi, dn, dl, deltaMax, Bz = Bz, 
                     progressOutput=True, debugOutput=False)
    
    calc.diagonalise(r, nEig, sortEigenvectors=True, progressOutput=True)
    
    i = np.argmin(np.abs(calc.r - 2.5))  # nearest to R = 5 μm
    

    #egvector = np.array(calc.composition)
    #print(egvector.shape)
    #for j in range(19):
       #print(egvector[i][j])
    R = 2.5  # Target distance            

    return

[nikolasibalic]

There was an issue in the library (incompatibility with new upstream likely). In version 3.10.2 the following works

Now this works

from arc import *
import numpy as np
import matplotlib.pyplot as plt
import matplotlib
matplotlib.use('qtagg')


def arc_eigenstate(a1,n1,a2,n2):
    if a1 == "Cs":
        atom1 = Caesium()
    elif a1 == "Rb":
        atom1 = Rubidium87()
    # n1 = 61
    l1 = 1
    j1 = 1.5
    m1 = 1.5
    if a2 == "Cs":
        atom2 = Caesium()
    elif a2 == "Rb":
        atom2 = Rubidium87()
    # n2 = 61
    l2 = 1
    j2 = 1.5
    m2 = 1.5
    

    Bz = 0.0001# #Tesla
    
    rmin = 1.5
    rmax = 6
    
    theta = 0 # Polar Angle [0-pi]
    phi = 0 # Azimuthal Angle [0-2pi]
    dn = 4 # Range of n to consider (n0-dn:n0+dn)
    dl = 4 # Range of l values
    deltaMax = 25e9  # Max pair-state energy difference [Hz]
    
    nEig = 6  # Number of eigenstates to extract
    
    
    calc = PairStateInteractions(atom1, n1, l1, j1, n2, l2, j2, m1, m2, s=0.5, 
                                 s2=0.5, atom2=atom2)
    
    r = np.linspace(rmin, rmax, 500)

    calc.defineBasis(theta, phi, dn, dl, deltaMax, Bz = Bz, 
                     progressOutput=True, debugOutput=False)
    
    calc.diagonalise(r, nEig, sortEigenvectors=True, progressOutput=True)

    #print(calc.x)

    #for j in range(len(calc.y[0])):
    j = 0
    x = []
    y = []
    for i in range(len(calc.r)):
        x.append(calc.r[i])
        y.append(calc.y[i][j])
        plt.plot(x,y,"b-")

    plt.xlabel(r"Inter-atomic distance ($\mu$m)")
    plt.ylabel(r"Pair-state relative energy, $\Delta/h$ (GHz)")
    plt.savefig("plot.png")  

    plt.close()

    calc.plotLevelDiagram()
    calc.showPlot()
    return

arc_eigenstate("Cs", 60, "Rb", 60)

Producing