CODEWITHSUNDEEP
python-3.xnumpyfor-loopmatplotlib
M00N KNIGHT
M00N KNIGHT

Reputation: 137

Plot for Function not looping

I'm trying to plot h vs. eig_gs but it only returns one value. First of all Hamiltonian is defined as follows:

def Hamiltonian(alpha,h):

    Sx = np.array([[0,1],[1,0]])
    Sy = np.array([[0,-1j],[1j,0]])
    Sz = np.array([[1,0],[0,-1]])
    I  = np.array([[1,0],[0,1]])

    H = ( (alpha*np.kron(np.kron(Sx,Sx),I))
       + (alpha*np.kron(np.kron(Sy,Sy),I))
       + (alpha*np.kron(np.kron(I,Sx),Sx))
       + (alpha*np.kron(np.kron(I,Sy),Sy))
       + (h*np.kron(np.kron(I,Sz),I)) )

    return H

Which returns an 8x8 matrix I can use for my work. My question is what is wrong in my code to stop it from looping over all oh h rather than just the last value? I have tried sticking H inside the for loop but this does not change the value and i have tried writing is as 

H = Hamiltonian(1,h.size)

(inside the for loop) But this does not solve the problem

# Computation of eigenvalues from density matrix

h = np.array([0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
eig_gs = np.zeros(h.size, dtype = 'complex')

for i in range(h.size):
    H = Hamiltonian(1,h.size)
    eigvals, eigvecs = LA.eigh(H)

    # Density of the GS
    g_state = eigvecs[:,0]
    rho_gs = np.outer(g_state, g_state.conjugate())

    # Expectation value of the energy GS density matrix
    eig_gs[i] = np.trace(rho_gs.dot(H))

plt.plot(h, eig_gs.real)
plt.show()

result

If I put the plt.plot in the for loop I get the following (which doesn't make sense):

enter image description here

Upvotes: 0

Views: 33

Answers (1)

Sheldore
Sheldore

Reputation: 39072

You need to put Hamiltonian inside the for loop and pass the corresponding h value while computing Hamiltonian. After appending all the eigenvalues, you can plot them outside the for loop. The result below shows the discrete energy spectrum obtained by solving your eigenvalue problem. The discrete energy levels can be visualized best using horizontal lines plotted using axhline

for i in range(h.size):
    H = Hamiltonian(1,h[i]) # Update the Hamiltonian
    eigvals, eigvecs = LA.eigh(H)
    # Density of the GS
    g_state = eigvecs[:,0]
    rho_gs = np.outer(g_state, g_state.conjugate())
    # Expectation value of the energy GS density matrix
    eig_gs[i] = np.trace(rho_gs.dot(H))

for i in range(h.size):
    plt.axhline(eig_gs[i].real)

enter image description here

Upvotes: 1

Related Questions