bootstrap.hankel_summed: GEVP method based on Hankel matrices.

View source: R/hankel.R

bootstrap.hankel_summedR Documentation

GEVP method based on Hankel matrices.

Description

Alternative method to determine energy levels from correlation matrices. A so-called Hankel matrix is generated from an input cf object and a generalised eigenvalue problem is solved then. This is the function to call. It will perform a bootstrap analysis.

Usage

bootstrap.hankel_summed(cf, t0values = c(1:(N - 2 * n - deltat)),
  deltat = 1, n = 2, N = cf$Time/2 + 1)

Arguments

cf

object of type cf

t0values

Integer vector. The t0 values to sum over. Default is c(1:max). All elements must be larger or equal to zero and smaller or equal than max=N-2*n-deltat

deltat

Integer. value of deltat used in the hankel GEVP. Default is 1.

n

Integer. Size of the Hankel matrices to generate, default is 2.

N

Integer. Maximal time index in correlation function to be used in Hankel matrix

Details

See vignette(name="hankel", package="hadron")

Value

List object of class "hankel.summed". The eigenvalues are stored in a numeric vector t0, the corresonding samples in t. The reference input times t0values is stored as t0values in the returned list. In addition, deltat is stored in the returned list.

See Also

Other hankel: bootstrap.hankel(), gevp.hankel_summed(), gevp.hankel(), hankel2cf(), hankel2effectivemass(), plot_hankel_spectrum()

Examples


data(correlatormatrix)
correlatormatrix <- bootstrap.cf(correlatormatrix, boot.R=99, boot.l=1, seed=132435)
t0 <- 4
correlatormatrix.gevp <- bootstrap.gevp(cf=correlatormatrix, t0=t0, element.order=c(1,2,3,4))
pc1 <- gevp2cf(gevp=correlatormatrix.gevp, id=1)
pc1.hankel <- bootstrap.hankel_summed(cf=pc1, t0=c(1:15), n=2)

hadron documentation built on Sept. 9, 2022, 5:06 p.m.