Linear Algebra

juobs.uweigvalsFunction
uweigvals(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 uwreal

  • b::Matrix{uwreal} : a matrix of uwreal, optional

It returns:

  • res = Vector{uwreal}: a vector where each elements is an eigenvalue
source
juobs.uweigvecsFunction
uweigvecs(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 uwreal

  • b::Matrix{uwreal} : a matrix of uwreal, optional

It returns:

  • res = Matrix{uwreal}: a matrix where each column is an eigenvector
source
juobs.uweigenFunction
uweigen(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 uwreal

  • b::Matrix{uwreal} : a matrix of uwreal, optional

It returns:

  • evals = Vector{uwreal}: a vector where each elements is an eigenvalue
  • evecs = Matrix{uwreal}: a matrix where the i-th column is the eigenvector of the i-th eigenvalue
source
juobs.get_matrixFunction
get_matrix(corr_diag::Vector{Array}, corr_upper::Vector{Array} )

This method returns an array of dim T where each element is a symmetrix matrix of dimension n of uwreal correlators at fixed time i=1..T. It takes as input:

corr_diag: vector of dimension n of correlators liying on the diagonal

corr_upper: vector of correlators liying on the upper diagonal.

Each correlator is an vector of uwreal variables of dimension T.

Example:

for i in 1:n
a[i,i]  = Vector{uwreal} # vector of uwreal variables of dimension T. They will constitute the diagonal elements of the matrices

for i in 1:n-1 
    for j in i+1:n
        a[i,j]  = Vector{uwreal} # vector of uwreal variables of dimension T. They will constitute the upper diagonal elements of the matrices. A matrix
        of dimension n*n has n(n-1)/2 upper diagonal elements.

Assume n=4

diagonal = Vector{Array}()
push!(diagonal, a[1,1],a[2,2],a[3,3],a[4,4])
upsize = Vector{Array}()
push!(upsize, a[1,2], a[1,3], a[1,4], a[2,3], a[2,4], a[3,4])

array_of_matrices = get_matrix(diagonal, upsize)

Julia> T-element Array{Array,1}

size(array_of_matrices)

Julia> (T,)

array_of_matrices[t] # t in 1:T

Julia> 4*4 Array{uwreal,2}
source
juobs.energiesFunction
energies(evals::Vector{Array})

Given a vector where each entry evals[t] is a uwreal array of eigenvalues, this method computes the effective energies of the first N states, where N=dim(evals[t]). The index t here runs from 1:T=lenght(evals), while the index i stands for the number of energy levels computed: i = length(evals[t]) It returns a vector array eff_en where each entry eff_en[t] contains the first N states energies as uwreal objects

source
juobs.getall_eigvalsFunction
getall_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 matrices

  • t0::Int64 : idex value at which the fixed matrix is taken

  • iter=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:

mat_array = get_matrix(diag, upper_diag)
evals = getall_eigvals(mat_array, 5)
source
juobs.getall_eigvecsFunction
getall_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 matrices

  • delta_t::Int64 : the fixed time shift t-t_0

  • iter=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)
source