lattice/pyerrors
1
0
Fork 0
mirror of https://github.com/fjosw/pyerrors.git synced 2025-05-14 19:43:41 +02:00
Code Issues Projects Releases Packages Wiki Activity

!feat: GEVP now returns all eigenvectors instead of just the ones for

Browse source 
the specified state.
This commit is contained in:
Fabian Joswig 2022-05-16 11:44:02 +01:00
parent ba054fa11c
commit 92b19cba9f
2 changed files with 10 additions and 15 deletions
Download patch file Download diff file Expand all files Collapse all files

18
pyerrors/correlators.py
View file

@ -241,8 +241,8 @@ class Corr:
if self.N == 1: if self.N == 1:
raise Exception("Trying to symmetrize a correlator matrix, that already has N=1.") raise Exception("Trying to symmetrize a correlator matrix, that already has N=1.")
def GEVP(self, t0, ts=None, state=0, sorted_list="Eigenvalue"): def GEVP(self, t0, ts=None, sorted_list="Eigenvalue"):
"""Solve the general eigenvalue problem on the current correlator """Solve the generalized eigenvalue problem on the current correlator and returns the corresponding eigenvectors.
Parameters Parameters
---------- ----------
@ -251,8 +251,6 @@ class Corr:
ts : int ts : int
fixed time G(t_s)v= lambda G(t_0)v if return_list=False fixed time G(t_s)v= lambda G(t_0)v if return_list=False
If return_list=True and sorting=Eigenvector it gives a reference point for the sorting method. If return_list=True and sorting=Eigenvector it gives a reference point for the sorting method.
state : int
The state one is interested in, ordered by energy. The lowest state is zero.
sorted_list : string sorted_list : string
if this argument is set, a list of vectors (len=self.T) is returned. If it is left as None, only one vector is returned. if this argument is set, a list of vectors (len=self.T) is returned. If it is left as None, only one vector is returned.
"Eigenvalue" - The eigenvector is chosen according to which eigenvalue it belongs individually on every timeslice. "Eigenvalue" - The eigenvector is chosen according to which eigenvalue it belongs individually on every timeslice.
@ -278,7 +276,7 @@ class Corr:
Gt[i, j] = symmetric_corr[ts][i, j].value Gt[i, j] = symmetric_corr[ts][i, j].value
sp_vecs = _GEVP_solver(Gt, G0) sp_vecs = _GEVP_solver(Gt, G0)
sp_vec = sp_vecs[state] sp_vec = [sp_vecs[s] for s in range(self.N)]
return sp_vec return sp_vec
elif sorted_list in ["Eigenvalue", "Eigenvector"]: elif sorted_list in ["Eigenvalue", "Eigenvector"]:
if sorted_list == "Eigenvalue" and ts is not None: if sorted_list == "Eigenvalue" and ts is not None:
@ -294,7 +292,7 @@ class Corr:
sp_vecs = _GEVP_solver(Gt, G0) sp_vecs = _GEVP_solver(Gt, G0)
if sorted_list == "Eigenvalue": if sorted_list == "Eigenvalue":
sp_vec = sp_vecs[state] sp_vec = [sp_vecs[s] for s in range(self.N)]
all_vecs.append(sp_vec) all_vecs.append(sp_vec)
else: else:
all_vecs.append(sp_vecs) all_vecs.append(sp_vecs)
@ -304,7 +302,7 @@ class Corr:
if (ts is None): if (ts is None):
raise Exception("ts is required for the Eigenvector sorting method.") raise Exception("ts is required for the Eigenvector sorting method.")
all_vecs = _sort_vectors(all_vecs, ts) all_vecs = _sort_vectors(all_vecs, ts)
all_vecs = [a[state] if a is not None else None for a in all_vecs] all_vecs = [[a[s] if a is not None else None for a in all_vecs] for s in range(self.N)]
else: else:
raise Exception("Unkown value for 'sorted_list'.") raise Exception("Unkown value for 'sorted_list'.")
@ -328,7 +326,7 @@ class Corr:
"Eigenvector" - Use the method described in arXiv:2004.10472 [hep-lat] to find the set of v(t) belonging to the state. "Eigenvector" - Use the method described in arXiv:2004.10472 [hep-lat] to find the set of v(t) belonging to the state.
The reference state is identified by its eigenvalue at t=ts The reference state is identified by its eigenvalue at t=ts
""" """
vec = self.GEVP(t0, ts=ts, state=state, sorted_list=sorted_list) vec = self.GEVP(t0, ts=ts, sorted_list=sorted_list)[state]
return self.projected(vec) return self.projected(vec)
def Hankel(self, N, periodic=False): def Hankel(self, N, periodic=False):
@ -1176,9 +1174,7 @@ class Corr:
if basematrix.N != self.N: if basematrix.N != self.N:
raise Exception('basematrix and targetmatrix have to be of the same size.') raise Exception('basematrix and targetmatrix have to be of the same size.')
evecs = [] evecs = basematrix.GEVP(t0proj, tproj, sorted_list=None)[:Ntrunc]
for i in range(Ntrunc):
evecs.append(basematrix.GEVP(t0proj, tproj, state=i, sorted_list=None))
tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object) tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
rmat = [] rmat = []

7
tests/correlators_test.py
View file

@ -231,11 +231,10 @@ def test_matrix_corr():
corr_mat = pe.Corr(np.array([[corr_aa, corr_ab], [corr_ab, corr_aa]])) corr_mat = pe.Corr(np.array([[corr_aa, corr_ab], [corr_ab, corr_aa]]))
corr_mat.item(0, 0) corr_mat.item(0, 0)
vec_0 = corr_mat.GEVP(0, 1, sorted_list=None) vecs = corr_mat.GEVP(0, 1, sorted_list=None)
vec_1 = corr_mat.GEVP(0, 1, state=1, sorted_list=None)
corr_0 = corr_mat.projected(vec_0) corr_0 = corr_mat.projected(vecs[0])
corr_1 = corr_mat.projected(vec_1) corr_1 = corr_mat.projected(vecs[0])
assert np.all([o == 0 for o in corr_0 - corr_aa]) assert np.all([o == 0 for o in corr_0 - corr_aa])
assert np.all([o == 0 for o in corr_1 - corr_aa]) assert np.all([o == 0 for o in corr_1 - corr_aa])