Linear Algebra
juobs.uweigvals
— Functionuweigvals(a::Matrix{uwreal}; iter = 30)
uweigvals(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)
This function computes the eigenvalues of a matrix of uwreal objects. If a second matrix b is given as input, it returns the generalised eigenvalues instead. It takes as input:
a::Matrix{uwreal}
: a matrix of uwrealb::Matrix{uwreal}
: a matrix of uwreal, optionaliter=30
: optional flag to set the iterations of the qr algorithm used to solve the eigenvalue problem
It returns:
res = Vector{uwreal}
: a vector where each elements is an eigenvalue
a = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries
b = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries
res = uweigvals(a) ##eigenvalues
res1 = uweigvals(a,b) ## generalised eigenvalues
juobs.uweigvecs
— Functionuweigvecs(a::Matrix{uwreal}; iter = 30)
uweigvecs(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)
This function computes the eigenvectors of a matrix of uwreal objects. If a second matrix b is given as input, it returns the generalised eigenvectors instead. It takes as input:
a::Matrix{uwreal}
: a matrix of uwrealb::Matrix{uwreal}
: a matrix of uwreal, optionaliter=30
: the number of iterations of the qr algorithm used to extract the eigenvalues
It returns:
res = Matrix{uwreal}
: a matrix where each column is an eigenvector
Examples:
a = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries
b = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries
res = uweigvecs(a) ##eigenvectors in column
res1 = uweigvecs(a,b) ## generalised eigenvectors in column
juobs.uweigen
— Functionuweigen(a::Matrix{uwreal}; iter = 30)
uweigen(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)
This function computes the eigenvalues and the eigenvectors of a matrix of uwreal objects. If a second matrix b is given as input, it returns the generalised eigenvalues and eigenvectors instead. It takes as input:
a::Matrix{uwreal}
: a matrix of uwrealb::Matrix{uwreal}
: a matrix of uwreal, optionaliter=30
: the number of iterations of the qr algorithm used to extract the eigenvalues
It returns:
evals = Vector{uwreal}
: a vector where each elements is an eigenvalueevecs = Matrix{uwreal}
: a matrix where the i-th column is the eigenvector of the i-th eigenvalue
Examples:
a = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries
b = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries
eval, evec = uweigen(a)
eval1, evec1 = uweigvecs(a,b)
juobs.get_matrix
— Functionget_matrix(diag::Vector{Array}, upper::Vector{Array} )
This function allows the user to build an array of matrices, where each matrix is a symmetric matrix of correlators at a given timeslice.
It takes as input:
diag::Vector{Array}
: vector of correlators. Each correlator will constitute a diagonal element of the matrices A[i] where i runs over the timeslices, i.e.A[i][1,1] = diag[1], .... A[i][n,n] = diag[n]
Given n=length(diag)
, the matrices will have dimension n*nupper::Vector{Array}
: vector of correlators liying on the upper diagonal position.A[i][1,2] = upper[1], .. , A[i][1,n] = upper[n-1],
A[i][2,3] = upper[n], .. , A[i][n-1,n] = upper[n(n-1)/2]
Given n,length(upper)=n(n-1)/2
The method returns an array of symmetric matrices of dimension n for each timeslice
Examples:
## load data
pp_data = read_mesons(path, "G5", "G5")
pa_data = read_mesons(path, "G5", "G0G5")
aa_data = read_mesons(path, "G0G5", "G0G5")
## create Corr struct
corr_pp = corr_obs.(pp_data)
corr_pa = corr_obs.(pa_data)
corr_aa = corr_obs.(aa_data) # array of correlators for different \mu_q combinations
## set up matrices
corr_diag = [corr_pp[1], corr_aa[1]]
corr_upper = [corr_pa[1]]
matrices = [corr_diag, corr_upper]
Julia> matrices[i]
pp[i] ap[i]
pa[i] aa[i]
## where i runs over the timeslices
juobs.energies
— Functionenergies(evals::Vector{Array}; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
This method computes the energy level from the eigenvalues according to:
$E_i(t) = \log(λ(t) / λ(t+1))$
where i=1,..,n
with n=length(evals[1])
and t=1,..,T
total time slices. It returns a vector array en where each entry en[i][t] contains the i-th states energy at time t
Examples:
## load data
pp_data = read_mesons(path, "G5", "G5")
pa_data = read_mesons(path, "G5", "G0G5")
aa_data = read_mesons(path, "G0G5", "G0G5")
## create Corr struct
corr_pp = corr_obs.(pp_data)
corr_pa = corr_obs.(pa_data)
corr_aa = corr_obs.(aa_data) # array of correlators for different \mu_q combinations
## set up matrices
corr_diag = [corr_pp[1], corr_aa[1]]
corr_upper = [corr_pa[1]]
matrices = [corr_diag, corr_upper]
## solve the GEVP
evals = getall_eigvals(matrices, 5) #where t_0=5
en = energies(evals)
Julia> en[i] # i-th state energy at each timeslice
juobs.getall_eigvals
— Functiongetall_eigvals(a::Vector{Matrix}, t0; iter=30 )
This function solves a GEVP problem, returning the eigenvalues, for a list of matrices, taking as generalised matrix the one at index t0, i.e:
$C(t_i)v_i = λ_i C(t_0) v_i$, with i=1:lenght(a)
It takes as input:
a::Vector{Matrix}
: a vector of matricest0::Int64
: idex value at which the fixed matrix is takeniter=30
: the number of iterations of the qr algorithm used to extract the eigenvalues
It returns:
res
= Vector{Vector{uwreal}}
where res[i]
are the generalised eigenvalues of the i-th matrix of the input array.
Examples:
## load data
pp_data = read_mesons(path, "G5", "G5")
pa_data = read_mesons(path, "G5", "G0G5")
aa_data = read_mesons(path, "G0G5", "G0G5")
## create Corr struct
corr_pp = corr_obs.(pp_data)
corr_pa = corr_obs.(pa_data)
corr_aa = corr_obs.(aa_data) # array of correlators for different \mu_q combinations
## set up matrices
corr_diag = [corr_pp[1], corr_aa[1]]
corr_upper = [corr_pa[1]]
matrices = [corr_diag, corr_upper]
## solve the GEVP
evals = getall_eigvals(matrices, 5) #where t_0=5
Julia>
juobs.getall_eigvecs
— Functiongetall_eigvecs(a::Vector{Matrix}, delta_t; iter=30 )
This function solves a GEVP problem, returning the eigenvectors, for a list of matrices.
$C(t_i)v_i = λ_i C(t_i-\delta_t) v_i$, with i=1:lenght(a)
Here delta_t
is the time shift within the two matrices of the problem, and is kept fixed. It takes as input:
a::Vector{Matrix}
: a vector of matricesdelta_t::Int64
: the fixed time shift t-t_0iter=30
: the number of iterations of the qr algorithm used to extract the eigenvalues
It returns:
res = Vector{Matrix{uwreal}}
where each res[i]
is a matrix with the eigenvectors as columns Examples:
mat_array = get_matrix(diag, upper_diag)
evecs = getall_eigvecs(mat_array, 5)