from __future__ import print_function
import numpy as np
def make_local_ops():
ops = {}
# c^+
cdag = np.zeros((2,2))
cdag[1, 0] = 1
ops['c^+'] = cdag
# c
c = cdag.transpose()
ops['c'] = c
# I
ops['I'] = np.identity(2)
# F (for fermionic anticommutation)
ops['F'] = np.diag((1.0, -1.0))
return ops
def make_matrix(*ops):
r = 1.0
for op in ops[::-1]:
r = np.kron(op, r)
return r
def hermite_conjugate(mat):
return mat.conj().T
def solve_hubbard_1site(U, mu=0):
local_ops = make_local_ops()
cdag = local_ops['c^+']
I = local_ops['I']
F = local_ops['F']
# creation operators
Cdag = {}
Cdag['1u'] = make_matrix(I, cdag) # F represents fermionic anticommutation
Cdag['1d'] = make_matrix(cdag, F)
C = {} # annihilation operators
N = {} # number operators
for key, cdag in Cdag.items():
C[key] = hermite_conjugate(cdag)
N[key] = cdag @ C[key]
hamil = U * N['1u'] @ N['1d'] - mu * (N['1u'] + N['1d'])
print("H =\n", hamil)
eigval, eigvec = np.linalg.eigh(hamil)
return eigval, eigvec
def main():
U = 10.0
mu = 2
E, vec = solve_hubbard_1site(U, mu)
print("\nE =\n", E)
print("\nEigenvectors =")
for i in range(vec.shape[1]):
print(vec[:, i])
if __name__ == '__main__':
main()