Commit 46040edf authored by Alessandro 's avatar Alessandro

Merge branch 'master' of gitlab.ift.uam-csic.es:jugarrio/juobs

parents 0f06a290 4ca6bf15
......@@ -2,3 +2,5 @@
analysis/plat\.txt
*.pdf
analysis/*.txt
......@@ -2,17 +2,23 @@
[[ADerrors]]
deps = ["BDIO", "FFTW", "FastGaussQuadrature", "ForwardDiff", "LaTeXStrings", "LinearAlgebra", "Nettle", "PGFPlotsX", "Plots", "Printf", "QuadGK", "Roots", "Statistics", "Test", "UnicodePlots"]
git-tree-sha1 = "f1c9954c17c0e87f185d532163e9e630efd12fe3"
repo-rev = "9c41ec2e161e4826bfe0752be8bf2a0bc025759c"
git-tree-sha1 = "2d3b701de1b206f497c632413b5a783fbd65d714"
repo-rev = "master"
repo-url = "https://gitlab.ift.uam-csic.es/alberto/aderrors.jl"
uuid = "5e92007d-7bf1-471c-8ceb-4591b8b567a9"
version = "0.1.0"
[[AbstractFFTs]]
deps = ["LinearAlgebra"]
git-tree-sha1 = "8ed9de2f1b1a9b1dee48582ad477c6e67b83eb2c"
git-tree-sha1 = "485ee0867925449198280d4af84bdb46a2a404d0"
uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c"
version = "1.0.0"
version = "1.0.1"
[[Adapt]]
deps = ["LinearAlgebra"]
git-tree-sha1 = "ffcfa2d345aaee0ef3d8346a073d5dd03c983ebe"
uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
version = "3.2.0"
[[ArgCheck]]
git-tree-sha1 = "dedbbb2ddb876f899585c4ec4433265e3017215a"
......@@ -197,9 +203,9 @@ version = "4.3.1+4"
[[FFTW]]
deps = ["AbstractFFTs", "FFTW_jll", "IntelOpenMP_jll", "Libdl", "LinearAlgebra", "MKL_jll", "Reexport"]
git-tree-sha1 = "c31e446bf3b12aad2ec8fc500fe19528c148d811"
git-tree-sha1 = "1b48dbde42f307e48685fa9213d8b9f8c0d87594"
uuid = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
version = "1.3.1"
version = "1.3.2"
[[FFTW_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
......@@ -275,9 +281,9 @@ version = "6.1.2+6"
[[GR]]
deps = ["Base64", "DelimitedFiles", "GR_jll", "HTTP", "JSON", "LinearAlgebra", "Pkg", "Printf", "Random", "Serialization", "Sockets", "Test", "UUIDs"]
git-tree-sha1 = "b90b826782cb3ac5b7a7f41b3fd0113180257ed4"
git-tree-sha1 = "12d971c928b7ecf19b748a2c7df6a365690dbf2c"
uuid = "28b8d3ca-fb5f-59d9-8090-bfdbd6d07a71"
version = "0.53.0"
version = "0.55.0"
[[GR_jll]]
deps = ["Artifacts", "Bzip2_jll", "Cairo_jll", "FFMPEG_jll", "Fontconfig_jll", "GLFW_jll", "JLLWrappers", "JpegTurbo_jll", "Libdl", "Libtiff_jll", "Pixman_jll", "Pkg", "Qt_jll", "Zlib_jll", "libpng_jll"]
......@@ -287,9 +293,9 @@ version = "0.53.0+0"
[[GeometryBasics]]
deps = ["EarCut_jll", "IterTools", "LinearAlgebra", "StaticArrays", "StructArrays", "Tables"]
git-tree-sha1 = "a28d728c2d825285fe27f38ca322399d35d1a5b9"
git-tree-sha1 = "4d4f72691933d5b6ee1ff20e27a102c3ae99d123"
uuid = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
version = "0.3.6"
version = "0.3.9"
[[Gettext_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Libiconv_jll", "Pkg", "XML2_jll"]
......@@ -309,10 +315,10 @@ uuid = "42e2da0e-8278-4e71-bc24-59509adca0fe"
version = "1.0.0"
[[HTTP]]
deps = ["Base64", "Dates", "IniFile", "MbedTLS", "Sockets"]
git-tree-sha1 = "c7ec02c4c6a039a98a15f955462cd7aea5df4508"
deps = ["Base64", "Dates", "IniFile", "MbedTLS", "NetworkOptions", "Sockets", "URIs"]
git-tree-sha1 = "c9f380c76d8aaa1fa7ea9cf97bddbc0d5b15adc2"
uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"
version = "0.8.19"
version = "0.9.5"
[[IOCapture]]
deps = ["Logging"]
......@@ -386,6 +392,18 @@ git-tree-sha1 = "3a0084cec7bf157edcb45a67fac0647f88fe5eaf"
uuid = "23fbe1c1-3f47-55db-b15f-69d7ec21a316"
version = "0.14.7"
[[LazyArtifacts]]
deps = ["Pkg"]
git-tree-sha1 = "4bb5499a1fc437342ea9ab7e319ede5a457c0968"
uuid = "4af54fe1-eca0-43a8-85a7-787d91b784e3"
version = "1.3.0"
[[LeastSquaresOptim]]
deps = ["FiniteDiff", "ForwardDiff", "LinearAlgebra", "Optim", "Printf", "SparseArrays", "Statistics", "SuiteSparse"]
git-tree-sha1 = "b5a1931bf37616820c4ef71629b6308d67f0b393"
uuid = "0fc2ff8b-aaa3-5acd-a817-1944a5e08891"
version = "0.8.1"
[[LibGit2]]
deps = ["Printf"]
uuid = "76f85450-5226-5b5a-8eaa-529ad045b433"
......@@ -467,10 +485,10 @@ uuid = "2fda8390-95c7-5789-9bda-21331edee243"
version = "0.12.0"
[[MKL_jll]]
deps = ["IntelOpenMP_jll", "Libdl", "Pkg"]
git-tree-sha1 = "eb540ede3aabb8284cb482aa41d00d6ca850b1f8"
deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "Pkg"]
git-tree-sha1 = "c253236b0ed414624b083e6b72bfe891fbd2c7af"
uuid = "856f044c-d86e-5d09-b602-aeab76dc8ba7"
version = "2020.2.254+0"
version = "2021.1.1+1"
[[MacroTools]]
deps = ["Markdown", "Random"]
......@@ -531,6 +549,11 @@ git-tree-sha1 = "b3b169cec36b850cfacbefc0b0a5f68a6089eb70"
uuid = "4c82536e-c426-54e4-b420-14f461c4ed8b"
version = "3.4.1+2"
[[NetworkOptions]]
git-tree-sha1 = "ed3157f48a05543cce9b241e1f2815f7e843d96e"
uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908"
version = "1.2.0"
[[Ogg_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "a42c0f138b9ebe8b58eba2271c5053773bde52d0"
......@@ -614,9 +637,9 @@ uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
[[PlotThemes]]
deps = ["PlotUtils", "Requires", "Statistics"]
git-tree-sha1 = "c6f5ea535551b3b16835134697f0c65d06c94b91"
git-tree-sha1 = "a3a964ce9dc7898193536002a6dd892b1b5a6f1d"
uuid = "ccf2f8ad-2431-5c83-bf29-c5338b663b6a"
version = "2.0.0"
version = "2.0.1"
[[PlotUtils]]
deps = ["ColorSchemes", "Colors", "Dates", "Printf", "Random", "Reexport", "Statistics"]
......@@ -626,9 +649,9 @@ version = "1.0.10"
[[Plots]]
deps = ["Base64", "Contour", "Dates", "FFMPEG", "FixedPointNumbers", "GR", "GeometryBasics", "JSON", "Latexify", "LinearAlgebra", "Measures", "NaNMath", "PlotThemes", "PlotUtils", "Printf", "REPL", "Random", "RecipesBase", "RecipesPipeline", "Reexport", "Requires", "Scratch", "Showoff", "SparseArrays", "Statistics", "StatsBase", "UUIDs"]
git-tree-sha1 = "4797acb266b8d9ff316f4581924e71c6709f152d"
git-tree-sha1 = "142dd04f5060c04de91cc10ca76dffb291a02426"
uuid = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
version = "1.10.1"
version = "1.10.6"
[[PositiveFactorizations]]
deps = ["LinearAlgebra"]
......@@ -679,9 +702,9 @@ version = "1.1.1"
[[RecipesPipeline]]
deps = ["Dates", "NaNMath", "PlotUtils", "RecipesBase"]
git-tree-sha1 = "9ea2f5bf1b26918b16e9f885bb8e05206bfc2144"
git-tree-sha1 = "c4d54a78e287de7ec73bbc928ce5eb3c60f80b24"
uuid = "01d81517-befc-4cb6-b9ec-a95719d0359c"
version = "0.2.1"
version = "0.3.1"
[[Reexport]]
deps = ["Pkg"]
......@@ -709,9 +732,9 @@ version = "0.2.2+1"
[[Roots]]
deps = ["Printf"]
git-tree-sha1 = "8f743e4f4368d1d753f3806bf635899dad6b4847"
git-tree-sha1 = "369e25546984dff5df351bc056fccc30de615080"
uuid = "f2b01f46-fcfa-551c-844a-d8ac1e96c665"
version = "1.0.7"
version = "1.0.8"
[[SHA]]
uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"
......@@ -777,10 +800,10 @@ uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c"
version = "0.9.6"
[[StructArrays]]
deps = ["DataAPI", "Tables"]
git-tree-sha1 = "ad1f5fd155426dcc879ec6ede9f74eb3a2d582df"
deps = ["Adapt", "DataAPI", "Tables"]
git-tree-sha1 = "26ea43b4be7e919a2390c3c0f824e7eb4fc19a0a"
uuid = "09ab397b-f2b6-538f-b94a-2f83cf4a842a"
version = "0.4.2"
version = "0.5.0"
[[SuiteSparse]]
deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"]
......@@ -794,14 +817,19 @@ version = "1.0.0"
[[Tables]]
deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "TableTraits", "Test"]
git-tree-sha1 = "240d19b8762006ff04b967bdd833269ad642d550"
git-tree-sha1 = "a716dde43d57fa537a19058d044b495301ba6565"
uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
version = "1.2.2"
version = "1.3.2"
[[Test]]
deps = ["Distributed", "InteractiveUtils", "Logging", "Random"]
uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
[[URIs]]
git-tree-sha1 = "7855809b88d7b16e9b029afd17880930626f54a2"
uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4"
version = "1.2.0"
[[UUIDs]]
deps = ["Random", "SHA"]
uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"
......@@ -983,9 +1011,9 @@ version = "1.2.11+18"
[[Zstd_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "6f1abcb0c44f184690912aa4b0ba861dd64f11b9"
git-tree-sha1 = "2c1332c54931e83f8f94d310fa447fd743e8d600"
uuid = "3161d3a3-bdf6-5164-811a-617609db77b4"
version = "1.4.5+2"
version = "1.4.8+0"
[[libass_jll]]
deps = ["Artifacts", "Bzip2_jll", "FreeType2_jll", "FriBidi_jll", "JLLWrappers", "Libdl", "Pkg", "Zlib_jll"]
......
......@@ -7,6 +7,7 @@ version = "0.1.0"
ADerrors = "5e92007d-7bf1-471c-8ceb-4591b8b567a9"
BDIO = "375f315e-f2c4-11e9-2ef9-134f02f79e27"
LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
LeastSquaresOptim = "0fc2ff8b-aaa3-5acd-a817-1944a5e08891"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
LsqFit = "2fda8390-95c7-5789-9bda-21331edee243"
Optim = "429524aa-4258-5aef-a3af-852621145aeb"
......
using juobs, ADerrors, DelimitedFiles, PyPlot, LaTeXStrings
const path = "/home/javier/Lattice/charm/production_2"
const path_plat = "/home/javier/Lattice/juobs/analysis/plat.txt"
const ensembles = ["H400", "N202", "N200", "N203", "N300", "J303"]
const deg = [true, true, false, false, true, false]
const L = [32, 48, 48, 48, 48, 64]
const beta = [3.46 , 3.55, 3.55, 3.55, 3.70, 3.70]
const R = [["H400r001", "H400r002"], "N202r001", ["N200r000", "N200r001"], ["N203r000", "N203r001"], "N300r002", "J303r003"]
include("/home/javier/Lattice/juobs/analysis/functions.jl")
include("/home/javier/Lattice/juobs/constants/juobs_const.jl")
m_ll = Vector{uwreal}(undef, length(ensembles))
m_ss = Vector{uwreal}(undef, length(ensembles))
m_hh = Vector{Vector{uwreal}}(undef, length(ensembles))
mu_pp = Vector{Vector{Vector{Float64}}}(undef, length(ensembles))
for iens = 1:length(ensembles)
pp = read_dat(R[iens], "G5", "G5")
a0p = read_dat(R[iens], "G5", "G0G5")
rew = read_rew(R[iens])
pp_obs = corr_obs.(pp, L=L[iens], rw=rew)
a0p_obs = corr_obs.(a0p, L=L[iens], rw=rew)
mu_pp[iens] = getfield.(pp_obs, :mu)
m = comp_pcac(a0p_obs, pp_obs, deg[iens], ensembles[iens])
m_ll[iens] = get_ll(mu_pp[iens], m, deg[iens])
m_ss[iens] = deg[iens] ? m_ll[iens] : get_ss(mu_pp[iens], m, deg[iens])
m_hh[iens] = get_hh(mu_pp[iens], m, deg[iens])
end
mm = get_mu.(mu_pp, deg)
mul_pp = getindex.(mm, 1)
mus_pp = getindex.(mm, 2)
muh_pp = getindex.(mm, 3)
cot_ll = Vector{uwreal}(undef, length(ensembles))
cot_ss = Vector{uwreal}(undef, length(ensembles))
cot_hh = Vector{Vector{uwreal}}(undef, length(ensembles))
for iens = 1:length(ensembles)
cot_ll[iens] = za(beta[iens]) * m_ll[iens] / mul_pp[iens]
cot_ss[iens] = za(beta[iens]) * m_ss[iens] / mus_pp[iens]
cot_hh[iens] = fill(za(beta[iens]), 3) .* m_hh[iens] ./ muh_pp[iens]
end
\ No newline at end of file
......@@ -8,14 +8,19 @@ function read_dat(rep::Vector{String}, g1::String="G5", g2::String="G5")
p = joinpath.(path, rep, "info")
f = [filter(x-> contains(x, ".dat"), readdir(p[k]))[1] for k = 1:length(p)]
f = joinpath.(p, f)
aux = read_mesons.(f, g1, g2)
res = []
for j = 1:length(aux[1])
aux2 = [aux[i][j] for i = 1:length(aux)]
push!(res, aux2)
return read_mesons(f, g1, g2)
end
function read_rew(rep::String)
p = joinpath(path, rep, "rew")
f = filter(x-> contains(x, ".dat"), readdir(p))
f = joinpath(p, f[1])
try
return read_ms1(f)
catch
return read_ms1(f, v="1.4")
end
return res
end
read_rew(rep::Vector{String}) = read_rew.(rep)
function get_mu(mu_list::Vector{Vector{Float64}}, deg::Bool)
mu_sorted = unique(sort(minimum.(mu_list)))
......@@ -48,6 +53,16 @@ function get_ls(mu_list, obs::Vector{uwreal}, deg::Bool)
end
end
function get_ss(mu_list, obs::Vector{uwreal}, deg::Bool)
mul, mus, muh = get_mu(mu_list, deg)
for k = 1:length(mu_list)
mu = mu_list[k]
if mus in mu && mu[1] == mu[2] #s-s
return obs[k]
end
end
end
function get_lh(mu_list, obs::Vector{uwreal}, deg::Bool)
mul, mus, muh = get_mu(mu_list, deg)
obs_lh = Vector{uwreal}(undef, 0)
......@@ -72,6 +87,17 @@ function get_sh(mu_list, obs::Vector{uwreal}, deg::Bool)
return obs_sh
end
function get_hh(mu_list, obs::Vector{uwreal}, deg::Bool)
mul, mus, muh = get_mu(mu_list, deg)
obs_hh = Vector{uwreal}(undef, 0)
for k = 1:length(mu_list)
mu = mu_list[k]
if mu[1] == mu[2] && !(mul in mu) && !(mus in mu)
push!(obs_hh, obs[k])
end
end
return obs_hh
end
function select_plateau(ens::String, mu_list, deg::Bool)
mul, mus, muh = get_mu(mu_list, deg)
f = readdlm(path_plat)
......@@ -116,6 +142,11 @@ function comp_meff(pp_obs::Vector{juobs.Corr}, deg::Bool, ens::String; pl::Bool=
plat = select_plateau(ens, mu_list, deg)
return meff.(pp_obs, plat, pl=pl)
end
function comp_pcac(a0p_obs::Vector{juobs.Corr}, pp_obs::Vector{juobs.Corr}, deg::Bool, ens::String; pl::Bool=false)
mu_list = getfield.(pp_obs, :mu)
plat = select_plateau(ens, mu_list, deg)
return mpcac.(a0p_obs, pp_obs, plat, pl=pl)
end
function comp_f(pp_obs::Vector{juobs.Corr}, m::Vector{uwreal}, deg::Bool, ens::String; pl::Bool=false)
mu_list = getfield.(pp_obs, :mu)
plat = select_plateau(ens, mu_list, deg)
......
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Home · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li class="is-active"><a class="tocitem" href="index.html">Home</a><ul class="internal"><li><a class="tocitem" href="#Contents"><span>Contents</span></a></li></ul></li><li><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="index.html">Home</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="index.html">Home</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs" title="Edit on GitLab"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitLab</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="DOCUMENTATION"><a class="docs-heading-anchor" href="#DOCUMENTATION">DOCUMENTATION</a><a id="DOCUMENTATION-1"></a><a class="docs-heading-anchor-permalink" href="#DOCUMENTATION" title="Permalink"></a></h1><h2 id="Contents"><a class="docs-heading-anchor" href="#Contents">Contents</a><a id="Contents-1"></a><a class="docs-heading-anchor-permalink" href="#Contents" title="Permalink"></a></h2><ul><li><a href="reader.html#Reader">Reader</a></li><li><a href="tools.html#Tools">Tools</a></li><li><a href="obs.html#Observables">Observables</a></li><li><a href="linalg.html#Linear-Algebra">Linear Algebra</a></li></ul></article><nav class="docs-footer"><a class="docs-footer-nextpage" href="reader.html">Reader »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Wednesday 3 March 2021 11:34">Wednesday 3 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Home · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li class="is-active"><a class="tocitem" href="index.html">Home</a><ul class="internal"><li><a class="tocitem" href="#Contents"><span>Contents</span></a></li></ul></li><li><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="index.html">Home</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="index.html">Home</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs" title="Edit on GitLab"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitLab</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="DOCUMENTATION"><a class="docs-heading-anchor" href="#DOCUMENTATION">DOCUMENTATION</a><a id="DOCUMENTATION-1"></a><a class="docs-heading-anchor-permalink" href="#DOCUMENTATION" title="Permalink"></a></h1><h2 id="Contents"><a class="docs-heading-anchor" href="#Contents">Contents</a><a id="Contents-1"></a><a class="docs-heading-anchor-permalink" href="#Contents" title="Permalink"></a></h2><ul><li><a href="reader.html#Reader">Reader</a></li><li><a href="tools.html#Tools">Tools</a></li><li><a href="obs.html#Observables">Observables</a></li><li><a href="linalg.html#Linear-Algebra">Linear Algebra</a></li></ul></article><nav class="docs-footer"><a class="docs-footer-nextpage" href="reader.html">Reader »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Tuesday 23 March 2021 11:06">Tuesday 23 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
......@@ -78,4 +78,4 @@ evals = getall_eigvals(matrices, 5) #where t_0=5
Julia&gt;</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.getall_eigvecs" href="#juobs.getall_eigvecs"><code>juobs.getall_eigvecs</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">getall_eigvecs(a::Vector{Matrix}, delta_t; iter=30 )</code></pre><p>This function solves a GEVP problem, returning the eigenvectors, for a list of matrices.</p><p><span>$C(t_i)v_i = λ_i C(t_i-\delta_t) v_i$</span>, with i=1:lenght(a)</p><p>Here <code>delta_t</code> is the time shift within the two matrices of the problem, and is kept fixed. It takes as input:</p><ul><li><p><code>a::Vector{Matrix}</code> : a vector of matrices</p></li><li><p><code>delta_t::Int64</code> : the fixed time shift t-t_0</p></li><li><p><code>iter=30</code> : the number of iterations of the qr algorithm used to extract the eigenvalues </p></li></ul><p>It returns:</p><ul><li><code>res = Vector{Matrix{uwreal}}</code></li></ul><p>where each <code>res[i]</code> is a matrix with the eigenvectors as columns Examples:</p><pre><code class="language-">mat_array = get_matrix(diag, upper_diag)
evecs = getall_eigvecs(mat_array, 5)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="obs.html">« Observables</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Wednesday 3 March 2021 11:34">Wednesday 3 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
evecs = getall_eigvecs(mat_array, 5)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="obs.html">« Observables</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Tuesday 23 March 2021 11:06">Tuesday 23 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Observables · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="index.html">Home</a></li><li><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li class="is-active"><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="obs.html">Observables</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="obs.html">Observables</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs" title="Edit on GitLab"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitLab</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Observables"><a class="docs-heading-anchor" href="#Observables">Observables</a><a id="Observables-1"></a><a class="docs-heading-anchor-permalink" href="#Observables" title="Permalink"></a></h1><article class="docstring"><header><a class="docstring-binding" id="juobs.meff" href="#juobs.meff"><code>juobs.meff</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false )
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Observables · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="index.html">Home</a></li><li><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li class="is-active"><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="obs.html">Observables</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="obs.html">Observables</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs" title="Edit on GitLab"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitLab</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Observables"><a class="docs-heading-anchor" href="#Observables">Observables</a><a id="Observables-1"></a><a class="docs-heading-anchor-permalink" href="#Observables" title="Permalink"></a></h1><article class="docstring"><header><a class="docstring-binding" id="juobs.meff" href="#juobs.meff"><code>juobs.meff</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false)</code></pre><p>Computes effective mass for a given correlator corr at a given plateau <code>plat</code>. Correlator can be passed as an <code>Corr</code> struct or <code>Vector{uwreal}</code>.</p><p>The flags <code>pl</code> and <code>data</code> allow to show the plots and return data as an extra result.</p><pre><code class="language-">data = read_mesons(path, &quot;G5&quot;, &quot;G5&quot;)
meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)</code></pre><p>Computes effective mass for a given correlator corr at a given plateau <code>plat</code>. Correlator can be passed as an <code>Corr</code> struct or <code>Vector{uwreal}</code>.</p><p>The flags <code>pl</code> and <code>data</code> allow to show the plots and return data as an extra result.</p><pre><code class="language-">data = read_mesons(path, &quot;G5&quot;, &quot;G5&quot;)
corr_pp = corr_obs.(data)
m = meff(corr_pp[1], [50, 60], pl=false)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.dec_const_pcvc" href="#juobs.dec_const_pcvc"><code>juobs.dec_const_pcvc</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu::Vector{Float64}, y0::Int64 ; pl::Bool=true, data::Bool=false)meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false)
m = meff(corr_pp[1], [50, 60], pl=false)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.mpcac" href="#juobs.mpcac"><code>juobs.mpcac</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">mpcac(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false)</code></pre><p>Computes decay constant using the PCVC relation for twisted mass fermions. The decay constant is computed in the plateau <code>plat</code>. Correlator can be passed as an <code>Corr</code> struct or <code>Vector{uwreal}</code>. If it is passed as a uwreal vector, vector of twisted masses <code>mu</code> and source position <code>y0</code> must be specified.</p><p>The flags <code>pl</code> and <code>data</code> allow to show the plots and return data as an extra result.</p><pre><code class="language-">data = read_mesons(path, &quot;G5&quot;, &quot;G5&quot;)
corr_pp = corr_obs.(data)
mpcac(a0p::Corr, pp::Corr, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)</code></pre><p>Computes the bare PCAC mass for a given correlator <code>a0p</code> and <code>pp</code> at a given plateau <code>plat</code>. Correlator can be passed as an <code>Corr</code> struct or <code>Vector{uwreal}</code>.</p><p>The flags <code>pl</code> and <code>data</code> allow to show the plots and return data as an extra result. The <code>ca</code> variable allows to compute <code>mpcac</code> using the improved axial current.</p><pre><code class="language-">data_pp = read_mesons(path, &quot;G5&quot;, &quot;G5&quot;)
data_a0p = read_mesons(path, &quot;G5&quot;, &quot;G0G5&quot;)
corr_pp = corr_obs.(data_pp)
corr_a0p = corr_obs.(data_a0p)
m12 = mpcac(corr_a0p, corr_pp, [50, 60], pl=false)
p0 = 9.2056
p1 = -13.9847
g2 = 1.73410
ca = -0.006033 * g2 *( 1 + exp(p0 + p1/g2))
m12 = mpcac(corr_a0p, corr_pp, [50, 60], pl=false, ca=ca)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.dec_const" href="#juobs.dec_const"><code>juobs.dec_const</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">dec_const(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, y0::Int64; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const(a0p::Corr, pp::Corr, plat::Vector{Int64}, m::uwreal; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)</code></pre><p>Computes the bare decay constant using <span>$A_0P$</span> and <span>$PP$</span> correlators . The decay constant is computed in the plateau <code>plat</code>. Correlator can be passed as an <code>Corr</code> struct or <code>Vector{uwreal}</code>. If it is passed as a uwreal vector, effective mass <code>m</code> and source position <code>y0</code> must be specified.</p><p>The flags <code>pl</code> and <code>data</code> allow to show the plots and return data as an extra result. The <code>ca</code> variable allows to compute <code>dec_const</code> using the improved axial current.</p><p><strong>The method assumes that the source is close to the boundary.</strong> It takes the following ratio to cancel boundary effects. <span>$R = \frac{f_A(x_0, y_0)}{\sqrt{f_P(T-y_0, y_0)}} * e^{m (x_0 - T/2)}$</span></p><pre><code class="language-">data_pp = read_mesons(path, &quot;G5&quot;, &quot;G5&quot;, legacy=true)
data_a0p = read_mesons(path, &quot;G5&quot;, &quot;G0G5&quot;, legacy=true)
corr_pp = corr_obs.(data_pp, L=32)
corr_a0p = corr_obs.(data_a0p, L=32)
m = meff(corr_pp[1], [50, 60], pl=false)
beta = 3.46
p0 = 9.2056
p1 = -13.9847
g2 = 6 / beta
ca = -0.006033 * g2 *( 1 + exp(p0 + p1/g2))
f = dec_const(corr_a0p[1], corr_pp[1], [50, 60], m, pl=true, ca=ca)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.dec_const_pcvc" href="#juobs.dec_const_pcvc"><code>juobs.dec_const_pcvc</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu::Vector{Float64}, y0::Int64 ; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)</code></pre><p>Computes decay constant using the PCVC relation for twisted mass fermions. The decay constant is computed in the plateau <code>plat</code>. Correlator can be passed as an <code>Corr</code> struct or <code>Vector{uwreal}</code>. If it is passed as a uwreal vector, vector of twisted masses <code>mu</code> and source position <code>y0</code> must be specified.</p><p>The flags <code>pl</code> and <code>data</code> allow to show the plots and return data as an extra result.</p><p><strong>The method assumes that the source is in the bulk.</strong></p><pre><code class="language-">data = read_mesons(path, &quot;G5&quot;, &quot;G5&quot;)
corr_pp = corr_obs.(data, L=32)
m = meff(corr_pp[1], [50, 60], pl=false)
f = dec_const_pcvc(corr_pp[1], [50, 60], m, pl=false)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.comp_t0" href="#juobs.comp_t0"><code>juobs.comp_t0</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2)
f = dec_const_pcvc(corr_pp[1], [50, 60], m, pl=false)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.comp_t0" href="#juobs.comp_t0"><code>juobs.comp_t0</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Vector{Matrix{Float64}}, Nothing}=nothing, npol::Int64=2)</code></pre><p>Computes <code>t0</code> using the energy density of the action <code>Ysl</code>(Yang-Mills action). <code>t0</code> is computed in the plateau <code>plat</code>. A polynomial interpolation in <code>t</code> is performed to find <code>t0</code>, where <code>npol</code> is the degree of the polynomial (linear fit by default)</p><p>The flag <code>pl</code> allows to show the plot.</p><pre><code class="language-">#Single replica
comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Vector{Matrix{Float64}}, Nothing}=nothing, npol::Int64=2, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)</code></pre><p>Computes <code>t0</code> using the energy density of the action <code>Ysl</code>(Yang-Mills action). <code>t0</code> is computed in the plateau <code>plat</code>. A polynomial interpolation in <code>t</code> is performed to find <code>t0</code>, where <code>npol</code> is the degree of the polynomial (linear fit by default)</p><p>The flag <code>pl</code> allows to show the plot.</p><pre><code class="language-">#Single replica
Y = read_ms(path)
rw = read_ms(path_rw)
......@@ -25,4 +54,4 @@ rw2 = read_ms(path_rw2)
t0 = comp_t0([Y1, Y2], [38, 58], L=32, pl=true)
t0_r = comp_t0(Y, [38, 58], L=32, rw=[rw1, rw2], pl=true)
</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="tools.html">« Tools</a><a class="docs-footer-nextpage" href="linalg.html">Linear Algebra »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Wednesday 3 March 2021 11:34">Wednesday 3 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="tools.html">« Tools</a><a class="docs-footer-nextpage" href="linalg.html">Linear Algebra »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Tuesday 23 March 2021 11:06">Tuesday 23 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Reader · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="index.html">Home</a></li><li class="is-active"><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="reader.html">Reader</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="reader.html">Reader</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs" title="Edit on GitLab"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitLab</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Reader"><a class="docs-heading-anchor" href="#Reader">Reader</a><a id="Reader-1"></a><a class="docs-heading-anchor-permalink" href="#Reader" title="Permalink"></a></h1><article class="docstring"><header><a class="docstring-binding" id="juobs.read_mesons" href="#juobs.read_mesons"><code>juobs.read_mesons</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Reader · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="index.html">Home</a></li><li class="is-active"><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="reader.html">Reader</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="reader.html">Reader</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs" title="Edit on GitLab"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitLab</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Reader"><a class="docs-heading-anchor" href="#Reader">Reader</a><a id="Reader-1"></a><a class="docs-heading-anchor-permalink" href="#Reader" title="Permalink"></a></h1><article class="docstring"><header><a class="docstring-binding" id="juobs.read_mesons" href="#juobs.read_mesons"><code>juobs.read_mesons</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)
read_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)</code></pre><p>This function read a mesons dat file at a given path and returns a vector of <code>CData</code> structures for different masses and Dirac structures. Dirac structures <code>g1</code> and/or <code>g2</code> can be passed as string arguments in order to filter correaltors. ADerrors id can be specified as argument. If is not specified, the <code>id</code> is fixed according to the ensemble name (example: &quot;H400&quot;-&gt; id = 400)</p><p>*For the old version (without smearing, distance preconditioning and theta) set legacy=true Examples:</p><pre><code class="language-">read_mesons(path)
read_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)</code></pre><p>This function read a mesons dat file at a given path and returns a vector of <code>CData</code> structures for different masses and Dirac structures. Dirac structures <code>g1</code> and/or <code>g2</code> can be passed as string arguments in order to filter correaltors. ADerrors id can be specified as argument. If is not specified, the <code>id</code> is fixed according to the ensemble name (example: &quot;H400&quot;-&gt; id = &quot;H400&quot;)</p><p>*For the old version (without smearing, distance preconditioning and theta) set legacy=true.</p><p>Examples:</p><pre><code class="language-">read_mesons(path)
read_mesons(path, &quot;G5&quot;)
read_mesons(path, nothing, &quot;G5&quot;)
read_mesons(path, &quot;G5&quot;, &quot;G5&quot;)
read_mesons(path, &quot;G5&quot;, &quot;G5&quot;, id=1)
read_mesons(path, &quot;G5_d2&quot;, &quot;G5_d2&quot;, legacy=true)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.read_ms" href="#juobs.read_ms"><code>juobs.read_ms</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_ms(path::String; id::Union{Int64, Nothing}=nothing, dtr::Int64=1)</code></pre><p>Reads openQCD ms dat files at a given path. This method return YData: </p><ul><li><p><code>t(t)</code>: flow time values</p></li><li><p><code>Ysl(icfg, x0, t)</code>: the time-slice sums of the densities of the Yang-Mills action </p></li><li><p><code>vtr</code>: vector that contains trajectory number</p></li><li><p><code>id</code>: ensmble id</p></li></ul><p>Examples:</p><pre><code class="language-">Y = read_ms(path)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.read_ms1" href="#juobs.read_ms1"><code>juobs.read_ms1</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_ms1(path::String; v::String=&quot;1.2&quot;)</code></pre><p>Reads openQCD ms1 dat files at a given path. This method returns a matrix <code>W[irw, icfg]</code> that contains the reweighting factors, where <code>irw</code> is the <code>rwf</code> index and icfg the configuration number. The function is compatible with the output files of openQCD v=1.2, 1.4 and 1.6. Version can be specified as argument.</p><p>Examples:</p><pre><code class="language-">read_ms1(path)
read_mesons(path, &quot;G5&quot;, &quot;G5&quot;, id=&quot;H100&quot;)
read_mesons(path, &quot;G5_d2&quot;, &quot;G5_d2&quot;, legacy=true)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.read_ms" href="#juobs.read_ms"><code>juobs.read_ms</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_ms(path::String; id::Union{String, Nothing}=nothing, dtr::Int64=1, obs::String=&quot;Y&quot;)</code></pre><p>Reads openQCD ms dat files at a given path. This method return YData: </p><ul><li><p><code>t(t)</code>: flow time values</p></li><li><p><code>obs(icfg, x0, t)</code>: the time-slice sums of the densities of the observable (Wsl, Ysl or Qsl)</p></li><li><p><code>vtr</code>: vector that contains trajectory number</p></li><li><p><code>id</code>: ensmble id</p></li></ul><p><code>dtr</code> = <code>dtr_cnfg</code> / <code>dtr_ms</code>, where <code>dtr_cnfg</code> is the number of trajectories computed before saving the configuration. <code>dtr_ms</code> is the same but applied to the ms.dat file.</p><p>Examples:</p><pre><code class="language-">Y = read_ms(path)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.read_ms1" href="#juobs.read_ms1"><code>juobs.read_ms1</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_ms1(path::String; v::String=&quot;1.2&quot;)</code></pre><p>Reads openQCD ms1 dat files at a given path. This method returns a matrix <code>W[irw, icfg]</code> that contains the reweighting factors, where <code>irw</code> is the <code>rwf</code> index and icfg the configuration number. The function is compatible with the output files of openQCD v=1.2, 1.4 and 1.6. Version can be specified as argument.</p><p>Examples:</p><pre><code class="language-">read_ms1(path)
read_ms1(path, v=&quot;1.4&quot;)
read_ms1(path, v=&quot;1.6&quot;)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.read_md" href="#juobs.read_md"><code>juobs.read_md</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">read_md(path::String)</code></pre><p>Reads openQCD pbp.dat files at a given path. This method returns a matrix <code>md[irw, icfg]</code> that contains the derivatives <span>$dS/dm$</span>, where <span>$md[irw=1] = dS/dm_l$</span> and <span>$md[irw=2] = dS/dm_s$</span></p><p><span>$Seff = -tr(log(D+m))$</span></p><p><span>$dSeff/ dm = -tr((D+m)^-1)$</span></p><p>Examples:</p><pre><code class="language-">md = read_md(path)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article><article class="docstring"><header><a class="docstring-binding" id="juobs.truncate_data!" href="#juobs.truncate_data!"><code>juobs.truncate_data!</code></a><span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia">truncate_data!(data::YData, nc::Int64)
......@@ -24,4 +24,4 @@ truncate_data!(Y, nc)
dat = read_mesons([path1, path2], &quot;G5&quot;, &quot;G5&quot;)
Y = read_ms.([path1_ms, path2_ms])
truncate_data!(dat, [nc1, nc2])
truncate_data!(Y, [nc1, nc2])</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="index.html">« Home</a><a class="docs-footer-nextpage" href="tools.html">Tools »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Wednesday 3 March 2021 11:34">Wednesday 3 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
truncate_data!(Y, [nc1, nc2])</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="index.html">« Home</a><a class="docs-footer-nextpage" href="tools.html">Tools »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Tuesday 23 March 2021 11:06">Tuesday 23 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Search · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="index.html">Home</a></li><li><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="search.html">Search</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="search.html">Search</a></li></ul></nav><div class="docs-right"><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article><p id="documenter-search-info">Loading search...</p><ul id="documenter-search-results"></ul></article><nav class="docs-footer"><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Wednesday 3 March 2021 11:34">Wednesday 3 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body><script src="search_index.js"></script><script src="assets/search.js"></script></html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Search · juobs Documentation</title><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit">juobs Documentation</span></div><form class="docs-search" action="search.html"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="index.html">Home</a></li><li><a class="tocitem" href="reader.html">Reader</a></li><li><a class="tocitem" href="tools.html">Tools</a></li><li><a class="tocitem" href="obs.html">Observables</a></li><li><a class="tocitem" href="linalg.html">Linear Algebra</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="search.html">Search</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="search.html">Search</a></li></ul></nav><div class="docs-right"><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article><p id="documenter-search-info">Loading search...</p><ul id="documenter-search-results"></ul></article><nav class="docs-footer"><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Tuesday 23 March 2021 11:06">Tuesday 23 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body><script src="search_index.js"></script><script src="assets/search.js"></script></html>
var documenterSearchIndex = {"docs":
[{"location":"reader.html#Reader","page":"Reader","title":"Reader","text":"","category":"section"},{"location":"reader.html","page":"Reader","title":"Reader","text":"read_mesons\nread_ms\nread_ms1\nread_md\ntruncate_data!","category":"page"},{"location":"reader.html#juobs.read_mesons","page":"Reader","title":"juobs.read_mesons","text":"read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)\n\nread_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)\n\nThis function read a mesons dat file at a given path and returns a vector of CData structures for different masses and Dirac structures. Dirac structures g1 and/or g2 can be passed as string arguments in order to filter correaltors. ADerrors id can be specified as argument. If is not specified, the id is fixed according to the ensemble name (example: \"H400\"-> id = 400)\n\n*For the old version (without smearing, distance preconditioning and theta) set legacy=true Examples:\n\nread_mesons(path)\nread_mesons(path, \"G5\")\nread_mesons(path, nothing, \"G5\")\nread_mesons(path, \"G5\", \"G5\")\nread_mesons(path, \"G5\", \"G5\", id=1)\nread_mesons(path, \"G5_d2\", \"G5_d2\", legacy=true)\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.read_ms","page":"Reader","title":"juobs.read_ms","text":"read_ms(path::String; id::Union{Int64, Nothing}=nothing, dtr::Int64=1)\n\nReads openQCD ms dat files at a given path. This method return YData: \n\nt(t): flow time values\nYsl(icfg, x0, t): the time-slice sums of the densities of the Yang-Mills action \nvtr: vector that contains trajectory number\nid: ensmble id\n\nExamples:\n\nY = read_ms(path)\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.read_ms1","page":"Reader","title":"juobs.read_ms1","text":"read_ms1(path::String; v::String=\"1.2\")\n\nReads openQCD ms1 dat files at a given path. This method returns a matrix W[irw, icfg] that contains the reweighting factors, where irw is the rwf index and icfg the configuration number. The function is compatible with the output files of openQCD v=1.2, 1.4 and 1.6. Version can be specified as argument.\n\nExamples:\n\nread_ms1(path)\nread_ms1(path, v=\"1.4\")\nread_ms1(path, v=\"1.6\")\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.read_md","page":"Reader","title":"juobs.read_md","text":"read_md(path::String)\n\nReads openQCD pbp.dat files at a given path. This method returns a matrix md[irw, icfg] that contains the derivatives dSdm, where mdirw=1 = dSdm_l and mdirw=2 = dSdm_s\n\nSeff = -tr(log(D+m))\n\ndSeff dm = -tr((D+m)^-1)\n\nExamples:\n\nmd = read_md(path)\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.truncate_data!","page":"Reader","title":"juobs.truncate_data!","text":"truncate_data!(data::YData, nc::Int64)\n\ntruncate_data!(data::Vector{YData}, nc::Vector{Int64})\n\ntruncate_data!(data::Vector{CData}, nc::Int64)\n\ntruncate_data!(data::Vector{Vector{CData}}, nc::Vector{Int64})\n\nTruncates the output of read_mesons and read_ms taking the first nc configurations.\n\nExamples:\n\n#Single replica\ndat = read_mesons(path, \"G5\", \"G5\")\nY = read_ms(path)\ntruncate_data!(dat, nc)\ntruncate_data!(Y, nc)\n\n#Two replicas\ndat = read_mesons([path1, path2], \"G5\", \"G5\")\nY = read_ms.([path1_ms, path2_ms])\ntruncate_data!(dat, [nc1, nc2])\ntruncate_data!(Y, [nc1, nc2])\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#Linear-Algebra","page":"Linear Algebra","title":"Linear Algebra","text":"","category":"section"},{"location":"linalg.html","page":"Linear Algebra","title":"Linear Algebra","text":"uweigvals\nuweigvecs\nuweigen\nget_matrix\nenergies\ngetall_eigvals\ngetall_eigvecs","category":"page"},{"location":"linalg.html#juobs.uweigvals","page":"Linear Algebra","title":"juobs.uweigvals","text":"uweigvals(a::Matrix{uwreal}; iter = 30)\n\nuweigvals(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)\n\nThis 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:\n\na::Matrix{uwreal} : a matrix of uwreal\nb::Matrix{uwreal} : a matrix of uwreal, optional\niter=30: optional flag to set the iterations of the qr algorithm used to solve the eigenvalue problem\n\nIt returns:\n\nres = Vector{uwreal}: a vector where each elements is an eigenvalue \n\na = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\nb = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\n\nres = uweigvals(a) ##eigenvalues\nres1 = uweigvals(a,b) ## generalised eigenvalues\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.uweigvecs","page":"Linear Algebra","title":"juobs.uweigvecs","text":"uweigvecs(a::Matrix{uwreal}; iter = 30)\n\nuweigvecs(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)\n\nThis 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:\n\na::Matrix{uwreal} : a matrix of uwreal\nb::Matrix{uwreal} : a matrix of uwreal, optional\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nres = Matrix{uwreal}: a matrix where each column is an eigenvector \n\nExamples:\n\na = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\nb = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\n\nres = uweigvecs(a) ##eigenvectors in column \nres1 = uweigvecs(a,b) ## generalised eigenvectors in column \n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.uweigen","page":"Linear Algebra","title":"juobs.uweigen","text":"uweigen(a::Matrix{uwreal}; iter = 30)\n\nuweigen(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)\n\nThis 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:\n\na::Matrix{uwreal} : a matrix of uwreal\nb::Matrix{uwreal} : a matrix of uwreal, optional\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nevals = Vector{uwreal}: a vector where each elements is an eigenvalue \nevecs = Matrix{uwreal}: a matrix where the i-th column is the eigenvector of the i-th eigenvalue\n\nExamples:\n\na = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\nb = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\n\neval, evec = uweigen(a) \neval1, evec1 = uweigvecs(a,b) \n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.get_matrix","page":"Linear Algebra","title":"juobs.get_matrix","text":"get_matrix(diag::Vector{Array}, upper::Vector{Array} )\n\nThis function allows the user to build an array of matrices, where each matrix is a symmetric matrix of correlators at a given timeslice. \n\nIt takes as input:\n\ndiag::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*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\n\nThe method returns an array of symmetric matrices of dimension n for each timeslice \n\nExamples:\n\n## load data \npp_data = read_mesons(path, \"G5\", \"G5\")\npa_data = read_mesons(path, \"G5\", \"G0G5\")\naa_data = read_mesons(path, \"G0G5\", \"G0G5\")\n\n## create Corr struct\ncorr_pp = corr_obs.(pp_data)\ncorr_pa = corr_obs.(pa_data)\ncorr_aa = corr_obs.(aa_data) # array of correlators for different \\mu_q combinations\n\n## set up matrices\ncorr_diag = [corr_pp[1], corr_aa[1]] \ncorr_upper = [corr_pa[1]]\n\nmatrices = [corr_diag, corr_upper]\n\nJulia> matrices[i]\n pp[i] ap[i]\n pa[i] aa[i]\n\n## where i runs over the timeslices\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.energies","page":"Linear Algebra","title":"juobs.energies","text":"energies(evals::Vector{Array}; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nThis method computes the energy level from the eigenvalues according to:\n\nE_i(t) = log(λ(t) λ(t+1))\n\nwhere 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 \n\nExamples:\n\n## load data\npp_data = read_mesons(path, \"G5\", \"G5\")\npa_data = read_mesons(path, \"G5\", \"G0G5\")\naa_data = read_mesons(path, \"G0G5\", \"G0G5\")\n\n## create Corr struct\ncorr_pp = corr_obs.(pp_data)\ncorr_pa = corr_obs.(pa_data)\ncorr_aa = corr_obs.(aa_data) # array of correlators for different \\mu_q combinations\n\n## set up matrices \ncorr_diag = [corr_pp[1], corr_aa[1]] \ncorr_upper = [corr_pa[1]]\n\nmatrices = [corr_diag, corr_upper]\n\n## solve the GEVP\nevals = getall_eigvals(matrices, 5) #where t_0=5\nen = energies(evals)\n\nJulia> en[i] # i-th state energy at each timeslice\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.getall_eigvals","page":"Linear Algebra","title":"juobs.getall_eigvals","text":"getall_eigvals(a::Vector{Matrix}, t0; iter=30 )\n\nThis function solves a GEVP problem, returning the eigenvalues, for a list of matrices, taking as generalised matrix the one at index t0, i.e:\n\nC(t_i)v_i = λ_i C(t_0) v_i, with i=1:lenght(a)\n\nIt takes as input:\n\na::Vector{Matrix} : a vector of matrices\nt0::Int64 : idex value at which the fixed matrix is taken\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nres = Vector{Vector{uwreal}}\n\nwhere res[i] are the generalised eigenvalues of the i-th matrix of the input array. \n\nExamples:\n\n## load data\npp_data = read_mesons(path, \"G5\", \"G5\")\npa_data = read_mesons(path, \"G5\", \"G0G5\")\naa_data = read_mesons(path, \"G0G5\", \"G0G5\")\n\n## create Corr struct\ncorr_pp = corr_obs.(pp_data)\ncorr_pa = corr_obs.(pa_data)\ncorr_aa = corr_obs.(aa_data) # array of correlators for different \\mu_q combinations\n\n## set up matrices \ncorr_diag = [corr_pp[1], corr_aa[1]] \ncorr_upper = [corr_pa[1]]\n\nmatrices = [corr_diag, corr_upper]\n\n## solve the GEVP\nevals = getall_eigvals(matrices, 5) #where t_0=5\n\n\nJulia>\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.getall_eigvecs","page":"Linear Algebra","title":"juobs.getall_eigvecs","text":"getall_eigvecs(a::Vector{Matrix}, delta_t; iter=30 )\n\nThis function solves a GEVP problem, returning the eigenvectors, for a list of matrices.\n\nC(t_i)v_i = λ_i C(t_i-delta_t) v_i, with i=1:lenght(a)\n\nHere delta_t is the time shift within the two matrices of the problem, and is kept fixed. It takes as input:\n\na::Vector{Matrix} : a vector of matrices\ndelta_t::Int64 : the fixed time shift t-t_0\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nres = Vector{Matrix{uwreal}}\n\nwhere each res[i] is a matrix with the eigenvectors as columns Examples:\n\nmat_array = get_matrix(diag, upper_diag)\nevecs = getall_eigvecs(mat_array, 5)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#Observables","page":"Observables","title":"Observables","text":"","category":"section"},{"location":"obs.html","page":"Observables","title":"Observables","text":"meff\ndec_const_pcvc\ncomp_t0","category":"page"},{"location":"obs.html#juobs.meff","page":"Observables","title":"juobs.meff","text":"meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false ) \n\nmeff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false)\n\nComputes effective mass for a given correlator corr at a given plateau plat. Correlator can be passed as an Corr struct or Vector{uwreal}.\n\nThe flags pl and data allow to show the plots and return data as an extra result.\n\ndata = read_mesons(path, \"G5\", \"G5\")\ncorr_pp = corr_obs.(data)\nm = meff(corr_pp[1], [50, 60], pl=false)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#juobs.dec_const_pcvc","page":"Observables","title":"juobs.dec_const_pcvc","text":"dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu::Vector{Float64}, y0::Int64 ; pl::Bool=true, data::Bool=false)meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false)\n\ndec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false)\n\nComputes decay constant using the PCVC relation for twisted mass fermions. The decay constant is computed in the plateau plat. Correlator can be passed as an Corr struct or Vector{uwreal}. If it is passed as a uwreal vector, vector of twisted masses mu and source position y0 must be specified.\n\nThe flags pl and data allow to show the plots and return data as an extra result.\n\ndata = read_mesons(path, \"G5\", \"G5\")\ncorr_pp = corr_obs.(data)\nm = meff(corr_pp[1], [50, 60], pl=false)\nf = dec_const_pcvc(corr_pp[1], [50, 60], m, pl=false)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#juobs.comp_t0","page":"Observables","title":"juobs.comp_t0","text":"comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2)\n\ncomp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Vector{Matrix{Float64}}, Nothing}=nothing, npol::Int64=2)\n\nComputes t0 using the energy density of the action Ysl(Yang-Mills action). t0 is computed in the plateau plat. A polynomial interpolation in t is performed to find t0, where npol is the degree of the polynomial (linear fit by default)\n\nThe flag pl allows to show the plot.\n\n#Single replica\nY = read_ms(path)\nrw = read_ms(path_rw)\n\nt0 = comp_t0(Y, [38, 58], L=32)\nt0_r = comp_t0(Y, [38, 58], L=32, rw=rw)\n\n#Two replicas\nY1 = read_ms(path1)\nY2 = read_ms(path2)\nrw1 = read_ms(path_rw1)\nrw2 = read_ms(path_rw2)\n\nt0 = comp_t0([Y1, Y2], [38, 58], L=32, pl=true)\nt0_r = comp_t0(Y, [38, 58], L=32, rw=[rw1, rw2], pl=true)\n\n\n\n\n\n\n","category":"function"},{"location":"index.html#DOCUMENTATION","page":"Home","title":"DOCUMENTATION","text":"","category":"section"},{"location":"index.html#Contents","page":"Home","title":"Contents","text":"","category":"section"},{"location":"index.html","page":"Home","title":"Home","text":"Pages = [\"reader.md\", \"tools.md\", \"obs.md\", \"linalg.md\"]\nDepth = 3","category":"page"},{"location":"tools.html#Tools","page":"Tools","title":"Tools","text":"","category":"section"},{"location":"tools.html","page":"Tools","title":"Tools","text":"corr_obs\nmd_sea\nmd_val\nlin_fit\nfit_routine","category":"page"},{"location":"tools.html#juobs.corr_obs","page":"Tools","title":"juobs.corr_obs","text":"corr_obs(cdata::CData; real::Bool=true, rw::Union{Array{Float64, 2}, Nothing}=nothing, L::Int64=1)\n\ncorr_obs(cdata::Array{CData, 1}; real::Bool=true, rw::Union{Array{Array{Float64, 2}, 1}, Nothing}=nothing, L::Int64=1)\n\nCreates a Corr struct with the given CData struct cdata (read_mesons) for a single replica. An array of CData can be passed as argument for multiple replicas.\n\nThe flag real select the real or imaginary part of the correlator. If rw is specified, the method applies reweighting. rw is passed as a matrix of Float64 (read_ms1) The correlator can be normalized with the volume factor if L is fixed.\n\n#Single replica\ndata = read_mesons(path, \"G5\", \"G5\")\nrw = read_ms1(path_rw)\ncorr_pp = corr_obs.(data)\ncorr_pp_r = corr_obs.(data, rw=rw)\n\n#Two replicas\ndata = read_mesons([path_r1, path_r2], \"G5\", \"G5\")\nrw1 = read_ms1(path_rw1)\nrw2 = read_ms1(path_rw2)\n\ncorr_pp = corr_obs.(data)\ncorr_pp_r = corr_obs.(data, rw=[rw1, rw2])\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.md_sea","page":"Tools","title":"juobs.md_sea","text":"md_sea(a::uwreal, md::Vector{Matrix{Float64}}, ws::ADerrors.wspace=ADerrors.wsg)\n\nComputes the derivative of an observable A with respect to the sea quark masses.\n\nfracd Adm(sea) = sum_i fracpartial Apartial O_i fracd O_id m(sea)\n\nfracd O_idm(sea) = O_i fracpartial Spartial m - O_i fracpartial Spartial m = - (O_i - O_i) (fracpartial Spartial m - fracpartial Spartial m)\n\nwhere O_i are primary observables \n\nmd is a vector that contains the derivative of the action S with respect to the sea quark masses for each replica. md[irep][irw, icfg]\n\nmd_sea returns a tuple of uwreal observables (dAdm_l dAdm_s)_sea, where m_l and m_s are the light and strange quark masses.\n\n#Single replica\ndata = read_mesons(path, \"G5\", \"G5\")\nmd = read_md(path_md)\nrw = read_ms1(path_rw)\n\ncorr_pp = corr_obs.(data, rw=rw)\nm = meff(corr_pp[1], plat)\nm_mdl, m_mds = md_sea(m, [md], ADerrors.wsg)\nm_shifted = m + 2 * dml * m_mdl + dms * m_mds\n\n#Two replicas\ndata = read_mesons([path_r1, path_r2], \"G5\", \"G5\")\nmd1 = read_md(path_md1)\nmd2 = read_md(path_md2)\n\ncorr_pp = corr_obs.(data)\nm = meff(corr_pp[1], plat)\nm_mdl, m_mds = md_sea(m, [md1, md2], ADerrors.wsg)\nm_shifted = m + 2 * dml * m_mdl + dms * m_mds\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.md_val","page":"Tools","title":"juobs.md_val","text":"md_val(a::uwreal, obs::Corr, derm::Vector{Corr})\n\nComputes the derivative of an observable A with respect to the valence quark masses.\n\nfracd Adm(val) = sum_i fracpartial Apartial O_i fracd O_id m(val)\n\nfracd O_idm(val) = fracpartial O_ipartial m(val)\n\nwhere O_i are primary observables \n\nmd is a vector that contains the derivative of the action S with respect to the sea quark masses for each replica. md[irep][irw, icfg]\n\nmd_val returns a tuple of uwreal observables (dAdm_1 dAdm_2)_val, where m_1 and m_2 are the correlator masses.\n\ndata = read_mesons(path, \"G5\", \"G5\", legacy=true)\ndata_d1 = read_mesons(path, \"G5_d1\", \"G5_d1\", legacy=true)\ndata_d2 = read_mesons(path, \"G5_d2\", \"G5_d2\", legacy=true)\n\nrw = read_ms1(path_rw)\n\ncorr_pp = corr_obs.(data, rw=rw)\ncorr_pp_d1 = corr_obs.(data_d1, rw=rw)\ncorr_pp_d2 = corr_obs.(data_d2, rw=rw)\nderm = [[corr_pp_d1[k], corr_pp_d2[k]] for k = 1:length(pp_d1)]\n\nm = meff(corr_pp[1], plat)\nm_md1, m_md2 = md_val(m, corr_pp[1], derm[1])\nm_shifted = m + 2 * dm1 * m_md1 + dm2 * m_md2\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.lin_fit","page":"Tools","title":"juobs.lin_fit","text":"lin_fit(x::Vector{<:Real}, y::Vector{uwreal})\n\nComputes a linear fit of uwreal data points y. This method return uwreal fit parameters and chisqexpected.\n\nfitp, csqexp = lin_fit(phi2, m2)\nm2_phys = fitp[1] + fitp[2] * phi2_phys\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.fit_routine","page":"Tools","title":"juobs.fit_routine","text":"fit_routine(model::Function, xdata::Array{<:Real}, ydata::Array{uwreal}, param::Int64=3; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nfit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}, param::Int64=3; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing, covar::Bool=false)\n\nGiven a model function with a number param of parameters and an array of uwreal, this function fit ydata with the given model and print fit information The method return an array upar with the best fit parameters with their errors. The flag wpm is an optional array of Float64 of lenght 4. The first three paramenters specify the criteria to determine the summation windows:\n\nvp[1]: The autocorrelation function is summed up to t = round(vp1).\nvp[2]: The sumation window is determined using U. Wolff poposal with S_tau = wpm2\nvp[3]: The autocorrelation function Gamma(t) is summed up a point where its error deltaGamma(t) is a factor vp[3] times larger than the signal.\n\nAn additional fourth parameter vp[4], tells ADerrors to add a tail to the error with tau_exp = wpm4. Negative values of wpm[1:4] are ignored and only one component of wpm[1:3] needs to be positive. If the flag covaris set to true, fit_routine takes into account covariances between x and y for each data point.\n\n@. model(x,p) = p[1] + p[2] * exp(-(p[3]-p[1])*x)\n@. model2(x,p) = p[1] + p[2] * x[:, 1] + (p[3] + p[4] * x[:, 1]) * x[:, 2]\nfit_routine(model, xdata, ydata, param=3)\nfit_routine(model, xdata, ydata, param=3, covar=true)\n\n\n\n\n\n","category":"function"}]
[{"location":"reader.html#Reader","page":"Reader","title":"Reader","text":"","category":"section"},{"location":"reader.html","page":"Reader","title":"Reader","text":"read_mesons\nread_ms\nread_ms1\nread_md\ntruncate_data!","category":"page"},{"location":"reader.html#juobs.read_mesons","page":"Reader","title":"juobs.read_mesons","text":"read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)\n\nread_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)\n\nThis function read a mesons dat file at a given path and returns a vector of CData structures for different masses and Dirac structures. Dirac structures g1 and/or g2 can be passed as string arguments in order to filter correaltors. ADerrors id can be specified as argument. If is not specified, the id is fixed according to the ensemble name (example: \"H400\"-> id = \"H400\")\n\n*For the old version (without smearing, distance preconditioning and theta) set legacy=true.\n\nExamples:\n\nread_mesons(path)\nread_mesons(path, \"G5\")\nread_mesons(path, nothing, \"G5\")\nread_mesons(path, \"G5\", \"G5\")\nread_mesons(path, \"G5\", \"G5\", id=\"H100\")\nread_mesons(path, \"G5_d2\", \"G5_d2\", legacy=true)\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.read_ms","page":"Reader","title":"juobs.read_ms","text":"read_ms(path::String; id::Union{String, Nothing}=nothing, dtr::Int64=1, obs::String=\"Y\")\n\nReads openQCD ms dat files at a given path. This method return YData: \n\nt(t): flow time values\nobs(icfg, x0, t): the time-slice sums of the densities of the observable (Wsl, Ysl or Qsl)\nvtr: vector that contains trajectory number\nid: ensmble id\n\ndtr = dtr_cnfg / dtr_ms, where dtr_cnfg is the number of trajectories computed before saving the configuration. dtr_ms is the same but applied to the ms.dat file.\n\nExamples:\n\nY = read_ms(path)\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.read_ms1","page":"Reader","title":"juobs.read_ms1","text":"read_ms1(path::String; v::String=\"1.2\")\n\nReads openQCD ms1 dat files at a given path. This method returns a matrix W[irw, icfg] that contains the reweighting factors, where irw is the rwf index and icfg the configuration number. The function is compatible with the output files of openQCD v=1.2, 1.4 and 1.6. Version can be specified as argument.\n\nExamples:\n\nread_ms1(path)\nread_ms1(path, v=\"1.4\")\nread_ms1(path, v=\"1.6\")\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.read_md","page":"Reader","title":"juobs.read_md","text":"read_md(path::String)\n\nReads openQCD pbp.dat files at a given path. This method returns a matrix md[irw, icfg] that contains the derivatives dSdm, where mdirw=1 = dSdm_l and mdirw=2 = dSdm_s\n\nSeff = -tr(log(D+m))\n\ndSeff dm = -tr((D+m)^-1)\n\nExamples:\n\nmd = read_md(path)\n\n\n\n\n\n","category":"function"},{"location":"reader.html#juobs.truncate_data!","page":"Reader","title":"juobs.truncate_data!","text":"truncate_data!(data::YData, nc::Int64)\n\ntruncate_data!(data::Vector{YData}, nc::Vector{Int64})\n\ntruncate_data!(data::Vector{CData}, nc::Int64)\n\ntruncate_data!(data::Vector{Vector{CData}}, nc::Vector{Int64})\n\nTruncates the output of read_mesons and read_ms taking the first nc configurations.\n\nExamples:\n\n#Single replica\ndat = read_mesons(path, \"G5\", \"G5\")\nY = read_ms(path)\ntruncate_data!(dat, nc)\ntruncate_data!(Y, nc)\n\n#Two replicas\ndat = read_mesons([path1, path2], \"G5\", \"G5\")\nY = read_ms.([path1_ms, path2_ms])\ntruncate_data!(dat, [nc1, nc2])\ntruncate_data!(Y, [nc1, nc2])\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#Linear-Algebra","page":"Linear Algebra","title":"Linear Algebra","text":"","category":"section"},{"location":"linalg.html","page":"Linear Algebra","title":"Linear Algebra","text":"uweigvals\nuweigvecs\nuweigen\nget_matrix\nenergies\ngetall_eigvals\ngetall_eigvecs","category":"page"},{"location":"linalg.html#juobs.uweigvals","page":"Linear Algebra","title":"juobs.uweigvals","text":"uweigvals(a::Matrix{uwreal}; iter = 30)\n\nuweigvals(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)\n\nThis 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:\n\na::Matrix{uwreal} : a matrix of uwreal\nb::Matrix{uwreal} : a matrix of uwreal, optional\niter=30: optional flag to set the iterations of the qr algorithm used to solve the eigenvalue problem\n\nIt returns:\n\nres = Vector{uwreal}: a vector where each elements is an eigenvalue \n\na = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\nb = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\n\nres = uweigvals(a) ##eigenvalues\nres1 = uweigvals(a,b) ## generalised eigenvalues\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.uweigvecs","page":"Linear Algebra","title":"juobs.uweigvecs","text":"uweigvecs(a::Matrix{uwreal}; iter = 30)\n\nuweigvecs(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)\n\nThis 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:\n\na::Matrix{uwreal} : a matrix of uwreal\nb::Matrix{uwreal} : a matrix of uwreal, optional\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nres = Matrix{uwreal}: a matrix where each column is an eigenvector \n\nExamples:\n\na = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\nb = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\n\nres = uweigvecs(a) ##eigenvectors in column \nres1 = uweigvecs(a,b) ## generalised eigenvectors in column \n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.uweigen","page":"Linear Algebra","title":"juobs.uweigen","text":"uweigen(a::Matrix{uwreal}; iter = 30)\n\nuweigen(a::Matrix{uwreal}, b::Matrix{uwreal}; iter = 30)\n\nThis 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:\n\na::Matrix{uwreal} : a matrix of uwreal\nb::Matrix{uwreal} : a matrix of uwreal, optional\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nevals = Vector{uwreal}: a vector where each elements is an eigenvalue \nevecs = Matrix{uwreal}: a matrix where the i-th column is the eigenvector of the i-th eigenvalue\n\nExamples:\n\na = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\nb = Matrix{uwreal}(nothing, n,n) ## n*n matrix of uwreal with nothing entries\n\neval, evec = uweigen(a) \neval1, evec1 = uweigvecs(a,b) \n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.get_matrix","page":"Linear Algebra","title":"juobs.get_matrix","text":"get_matrix(diag::Vector{Array}, upper::Vector{Array} )\n\nThis function allows the user to build an array of matrices, where each matrix is a symmetric matrix of correlators at a given timeslice. \n\nIt takes as input:\n\ndiag::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*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\n\nThe method returns an array of symmetric matrices of dimension n for each timeslice \n\nExamples:\n\n## load data \npp_data = read_mesons(path, \"G5\", \"G5\")\npa_data = read_mesons(path, \"G5\", \"G0G5\")\naa_data = read_mesons(path, \"G0G5\", \"G0G5\")\n\n## create Corr struct\ncorr_pp = corr_obs.(pp_data)\ncorr_pa = corr_obs.(pa_data)\ncorr_aa = corr_obs.(aa_data) # array of correlators for different \\mu_q combinations\n\n## set up matrices\ncorr_diag = [corr_pp[1], corr_aa[1]] \ncorr_upper = [corr_pa[1]]\n\nmatrices = [corr_diag, corr_upper]\n\nJulia> matrices[i]\n pp[i] ap[i]\n pa[i] aa[i]\n\n## where i runs over the timeslices\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.energies","page":"Linear Algebra","title":"juobs.energies","text":"energies(evals::Vector{Array}; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nThis method computes the energy level from the eigenvalues according to:\n\nE_i(t) = log(λ(t) λ(t+1))\n\nwhere 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 \n\nExamples:\n\n## load data\npp_data = read_mesons(path, \"G5\", \"G5\")\npa_data = read_mesons(path, \"G5\", \"G0G5\")\naa_data = read_mesons(path, \"G0G5\", \"G0G5\")\n\n## create Corr struct\ncorr_pp = corr_obs.(pp_data)\ncorr_pa = corr_obs.(pa_data)\ncorr_aa = corr_obs.(aa_data) # array of correlators for different \\mu_q combinations\n\n## set up matrices \ncorr_diag = [corr_pp[1], corr_aa[1]] \ncorr_upper = [corr_pa[1]]\n\nmatrices = [corr_diag, corr_upper]\n\n## solve the GEVP\nevals = getall_eigvals(matrices, 5) #where t_0=5\nen = energies(evals)\n\nJulia> en[i] # i-th state energy at each timeslice\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.getall_eigvals","page":"Linear Algebra","title":"juobs.getall_eigvals","text":"getall_eigvals(a::Vector{Matrix}, t0; iter=30 )\n\nThis function solves a GEVP problem, returning the eigenvalues, for a list of matrices, taking as generalised matrix the one at index t0, i.e:\n\nC(t_i)v_i = λ_i C(t_0) v_i, with i=1:lenght(a)\n\nIt takes as input:\n\na::Vector{Matrix} : a vector of matrices\nt0::Int64 : idex value at which the fixed matrix is taken\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nres = Vector{Vector{uwreal}}\n\nwhere res[i] are the generalised eigenvalues of the i-th matrix of the input array. \n\nExamples:\n\n## load data\npp_data = read_mesons(path, \"G5\", \"G5\")\npa_data = read_mesons(path, \"G5\", \"G0G5\")\naa_data = read_mesons(path, \"G0G5\", \"G0G5\")\n\n## create Corr struct\ncorr_pp = corr_obs.(pp_data)\ncorr_pa = corr_obs.(pa_data)\ncorr_aa = corr_obs.(aa_data) # array of correlators for different \\mu_q combinations\n\n## set up matrices \ncorr_diag = [corr_pp[1], corr_aa[1]] \ncorr_upper = [corr_pa[1]]\n\nmatrices = [corr_diag, corr_upper]\n\n## solve the GEVP\nevals = getall_eigvals(matrices, 5) #where t_0=5\n\n\nJulia>\n\n\n\n\n\n","category":"function"},{"location":"linalg.html#juobs.getall_eigvecs","page":"Linear Algebra","title":"juobs.getall_eigvecs","text":"getall_eigvecs(a::Vector{Matrix}, delta_t; iter=30 )\n\nThis function solves a GEVP problem, returning the eigenvectors, for a list of matrices.\n\nC(t_i)v_i = λ_i C(t_i-delta_t) v_i, with i=1:lenght(a)\n\nHere delta_t is the time shift within the two matrices of the problem, and is kept fixed. It takes as input:\n\na::Vector{Matrix} : a vector of matrices\ndelta_t::Int64 : the fixed time shift t-t_0\niter=30 : the number of iterations of the qr algorithm used to extract the eigenvalues \n\nIt returns:\n\nres = Vector{Matrix{uwreal}}\n\nwhere each res[i] is a matrix with the eigenvectors as columns Examples:\n\nmat_array = get_matrix(diag, upper_diag)\nevecs = getall_eigvecs(mat_array, 5)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#Observables","page":"Observables","title":"Observables","text":"","category":"section"},{"location":"obs.html","page":"Observables","title":"Observables","text":"meff\nmpcac\ndec_const\ndec_const_pcvc\ncomp_t0","category":"page"},{"location":"obs.html#juobs.meff","page":"Observables","title":"juobs.meff","text":"meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing) \n\nmeff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nComputes effective mass for a given correlator corr at a given plateau plat. Correlator can be passed as an Corr struct or Vector{uwreal}.\n\nThe flags pl and data allow to show the plots and return data as an extra result.\n\ndata = read_mesons(path, \"G5\", \"G5\")\ncorr_pp = corr_obs.(data)\nm = meff(corr_pp[1], [50, 60], pl=false)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#juobs.mpcac","page":"Observables","title":"juobs.mpcac","text":"mpcac(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nmpcac(a0p::Corr, pp::Corr, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nComputes the bare PCAC mass for a given correlator a0p and pp at a given plateau plat. Correlator can be passed as an Corr struct or Vector{uwreal}.\n\nThe flags pl and data allow to show the plots and return data as an extra result. The ca variable allows to compute mpcac using the improved axial current.\n\ndata_pp = read_mesons(path, \"G5\", \"G5\")\ndata_a0p = read_mesons(path, \"G5\", \"G0G5\")\ncorr_pp = corr_obs.(data_pp)\ncorr_a0p = corr_obs.(data_a0p)\nm12 = mpcac(corr_a0p, corr_pp, [50, 60], pl=false)\n\np0 = 9.2056\np1 = -13.9847\ng2 = 1.73410\nca = -0.006033 * g2 *( 1 + exp(p0 + p1/g2))\n\nm12 = mpcac(corr_a0p, corr_pp, [50, 60], pl=false, ca=ca)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#juobs.dec_const","page":"Observables","title":"juobs.dec_const","text":"dec_const(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, y0::Int64; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\ndec_const(a0p::Corr, pp::Corr, plat::Vector{Int64}, m::uwreal; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nComputes the bare decay constant using A_0P and PP correlators . The decay constant is computed in the plateau plat. Correlator can be passed as an Corr struct or Vector{uwreal}. If it is passed as a uwreal vector, effective mass m and source position y0 must be specified.\n\nThe flags pl and data allow to show the plots and return data as an extra result. The ca variable allows to compute dec_const using the improved axial current.\n\nThe method assumes that the source is close to the boundary. It takes the following ratio to cancel boundary effects. R = fracf_A(x_0 y_0)sqrtf_P(T-y_0 y_0) * e^m (x_0 - T2)\n\ndata_pp = read_mesons(path, \"G5\", \"G5\", legacy=true)\ndata_a0p = read_mesons(path, \"G5\", \"G0G5\", legacy=true)\n\ncorr_pp = corr_obs.(data_pp, L=32)\ncorr_a0p = corr_obs.(data_a0p, L=32)\n\nm = meff(corr_pp[1], [50, 60], pl=false)\n\nbeta = 3.46\np0 = 9.2056\np1 = -13.9847\ng2 = 6 / beta\nca = -0.006033 * g2 *( 1 + exp(p0 + p1/g2))\n\nf = dec_const(corr_a0p[1], corr_pp[1], [50, 60], m, pl=true, ca=ca)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#juobs.dec_const_pcvc","page":"Observables","title":"juobs.dec_const_pcvc","text":"dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu::Vector{Float64}, y0::Int64 ; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\ndec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nComputes decay constant using the PCVC relation for twisted mass fermions. The decay constant is computed in the plateau plat. Correlator can be passed as an Corr struct or Vector{uwreal}. If it is passed as a uwreal vector, vector of twisted masses mu and source position y0 must be specified.\n\nThe flags pl and data allow to show the plots and return data as an extra result.\n\nThe method assumes that the source is in the bulk.\n\ndata = read_mesons(path, \"G5\", \"G5\")\ncorr_pp = corr_obs.(data, L=32)\nm = meff(corr_pp[1], [50, 60], pl=false)\nf = dec_const_pcvc(corr_pp[1], [50, 60], m, pl=false)\n\n\n\n\n\n","category":"function"},{"location":"obs.html#juobs.comp_t0","page":"Observables","title":"juobs.comp_t0","text":"comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\ncomp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Vector{Matrix{Float64}}, Nothing}=nothing, npol::Int64=2, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nComputes t0 using the energy density of the action Ysl(Yang-Mills action). t0 is computed in the plateau plat. A polynomial interpolation in t is performed to find t0, where npol is the degree of the polynomial (linear fit by default)\n\nThe flag pl allows to show the plot.\n\n#Single replica\nY = read_ms(path)\nrw = read_ms(path_rw)\n\nt0 = comp_t0(Y, [38, 58], L=32)\nt0_r = comp_t0(Y, [38, 58], L=32, rw=rw)\n\n#Two replicas\nY1 = read_ms(path1)\nY2 = read_ms(path2)\nrw1 = read_ms(path_rw1)\nrw2 = read_ms(path_rw2)\n\nt0 = comp_t0([Y1, Y2], [38, 58], L=32, pl=true)\nt0_r = comp_t0(Y, [38, 58], L=32, rw=[rw1, rw2], pl=true)\n\n\n\n\n\n\n","category":"function"},{"location":"index.html#DOCUMENTATION","page":"Home","title":"DOCUMENTATION","text":"","category":"section"},{"location":"index.html#Contents","page":"Home","title":"Contents","text":"","category":"section"},{"location":"index.html","page":"Home","title":"Home","text":"Pages = [\"reader.md\", \"tools.md\", \"obs.md\", \"linalg.md\"]\nDepth = 3","category":"page"},{"location":"tools.html#Tools","page":"Tools","title":"Tools","text":"","category":"section"},{"location":"tools.html","page":"Tools","title":"Tools","text":"corr_obs\nmd_sea\nmd_val\nlin_fit\nfit_routine","category":"page"},{"location":"tools.html#juobs.corr_obs","page":"Tools","title":"juobs.corr_obs","text":"corr_obs(cdata::CData; real::Bool=true, rw::Union{Array{Float64, 2}, Nothing}=nothing, L::Int64=1)\n\ncorr_obs(cdata::Array{CData, 1}; real::Bool=true, rw::Union{Array{Array{Float64, 2}, 1}, Nothing}=nothing, L::Int64=1)\n\nCreates a Corr struct with the given CData struct cdata (read_mesons) for a single replica. An array of CData can be passed as argument for multiple replicas.\n\nThe flag real select the real or imaginary part of the correlator. If rw is specified, the method applies reweighting. rw is passed as a matrix of Float64 (read_ms1) The correlator can be normalized with the volume factor if L is fixed.\n\n#Single replica\ndata = read_mesons(path, \"G5\", \"G5\")\nrw = read_ms1(path_rw)\ncorr_pp = corr_obs.(data)\ncorr_pp_r = corr_obs.(data, rw=rw)\n\n#Two replicas\ndata = read_mesons([path_r1, path_r2], \"G5\", \"G5\")\nrw1 = read_ms1(path_rw1)\nrw2 = read_ms1(path_rw2)\n\ncorr_pp = corr_obs.(data)\ncorr_pp_r = corr_obs.(data, rw=[rw1, rw2])\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.md_sea","page":"Tools","title":"juobs.md_sea","text":"md_sea(a::uwreal, md::Vector{Matrix{Float64}}, ws::ADerrors.wspace=ADerrors.wsg)\n\nComputes the derivative of an observable A with respect to the sea quark masses.\n\nfracd Adm(sea) = sum_i fracpartial Apartial O_i fracd O_id m(sea)\n\nfracd O_idm(sea) = O_i fracpartial Spartial m - O_i fracpartial Spartial m = - (O_i - O_i) (fracpartial Spartial m - fracpartial Spartial m)\n\nwhere O_i are primary observables \n\nmd is a vector that contains the derivative of the action S with respect to the sea quark masses for each replica. md[irep][irw, icfg]\n\nmd_sea returns a tuple of uwreal observables (dAdm_l dAdm_s)_sea, where m_l and m_s are the light and strange quark masses.\n\n#Single replica\ndata = read_mesons(path, \"G5\", \"G5\")\nmd = read_md(path_md)\nrw = read_ms1(path_rw)\n\ncorr_pp = corr_obs.(data, rw=rw)\nm = meff(corr_pp[1], plat)\nm_mdl, m_mds = md_sea(m, [md], ADerrors.wsg)\nm_shifted = m + 2 * dml * m_mdl + dms * m_mds\n\n#Two replicas\ndata = read_mesons([path_r1, path_r2], \"G5\", \"G5\")\nmd1 = read_md(path_md1)\nmd2 = read_md(path_md2)\n\ncorr_pp = corr_obs.(data)\nm = meff(corr_pp[1], plat)\nm_mdl, m_mds = md_sea(m, [md1, md2], ADerrors.wsg)\nm_shifted = m + 2 * dml * m_mdl + dms * m_mds\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.md_val","page":"Tools","title":"juobs.md_val","text":"md_val(a::uwreal, obs::Corr, derm::Vector{Corr})\n\nComputes the derivative of an observable A with respect to the valence quark masses.\n\nfracd Adm(val) = sum_i fracpartial Apartial O_i fracd O_id m(val)\n\nfracd O_idm(val) = fracpartial O_ipartial m(val)\n\nwhere O_i are primary observables \n\nmd is a vector that contains the derivative of the action S with respect to the sea quark masses for each replica. md[irep][irw, icfg]\n\nmd_val returns a tuple of uwreal observables (dAdm_1 dAdm_2)_val, where m_1 and m_2 are the correlator masses.\n\ndata = read_mesons(path, \"G5\", \"G5\", legacy=true)\ndata_d1 = read_mesons(path, \"G5_d1\", \"G5_d1\", legacy=true)\ndata_d2 = read_mesons(path, \"G5_d2\", \"G5_d2\", legacy=true)\n\nrw = read_ms1(path_rw)\n\ncorr_pp = corr_obs.(data, rw=rw)\ncorr_pp_d1 = corr_obs.(data_d1, rw=rw)\ncorr_pp_d2 = corr_obs.(data_d2, rw=rw)\nderm = [[corr_pp_d1[k], corr_pp_d2[k]] for k = 1:length(pp_d1)]\n\nm = meff(corr_pp[1], plat)\nm_md1, m_md2 = md_val(m, corr_pp[1], derm[1])\nm_shifted = m + 2 * dm1 * m_md1 + dm2 * m_md2\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.lin_fit","page":"Tools","title":"juobs.lin_fit","text":"lin_fit(x::Vector{<:Real}, y::Vector{uwreal})\n\nComputes a linear fit of uwreal data points y. This method return uwreal fit parameters and chisqexpected.\n\nfitp, csqexp = lin_fit(phi2, m2)\nm2_phys = fitp[1] + fitp[2] * phi2_phys\n\n\n\n\n\n","category":"function"},{"location":"tools.html#juobs.fit_routine","page":"Tools","title":"juobs.fit_routine","text":"fit_routine(model::Function, xdata::Array{<:Real}, ydata::Array{uwreal}, param::Int64=3; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)\n\nfit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}, param::Int64=3; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing, covar::Bool=false)\n\nGiven a model function with a number param of parameters and an array of uwreal, this function fit ydata with the given model and print fit information The method return an array upar with the best fit parameters with their errors. The flag wpm is an optional array of Float64 of lenght 4. The first three paramenters specify the criteria to determine the summation windows:\n\nvp[1]: The autocorrelation function is summed up to t = round(vp1).\nvp[2]: The sumation window is determined using U. Wolff poposal with S_tau = wpm2\nvp[3]: The autocorrelation function Gamma(t) is summed up a point where its error deltaGamma(t) is a factor vp[3] times larger than the signal.\n\nAn additional fourth parameter vp[4], tells ADerrors to add a tail to the error with tau_exp = wpm4. Negative values of wpm[1:4] are ignored and only one component of wpm[1:3] needs to be positive. If the flag covaris set to true, fit_routine takes into account covariances between x and y for each data point.\n\n@. model(x,p) = p[1] + p[2] * exp(-(p[3]-p[1])*x)\n@. model2(x,p) = p[1] + p[2] * x[:, 1] + (p[3] + p[4] * x[:, 1]) * x[:, 2]\nfit_routine(model, xdata, ydata, param=3)\nfit_routine(model, xdata, ydata, param=3, covar=true)\n\n\n\n\n\n","category":"function"}]
}
......@@ -50,4 +50,4 @@ m2_phys = fitp[1] + fitp[2] * phi2_phys</code></pre></div><a class="docs-sourcel
fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}, param::Int64=3; wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing, covar::Bool=false)</code></pre><p>Given a model function with a number param of parameters and an array of <code>uwreal</code>, this function fit ydata with the given <code>model</code> and print fit information The method return an array <code>upar</code> with the best fit parameters with their errors. The flag <code>wpm</code> is an optional array of Float64 of lenght 4. The first three paramenters specify the criteria to determine the summation windows:</p><ul><li><p><code>vp[1]</code>: The autocorrelation function is summed up to <span>$t = round(vp[1])$</span>.</p></li><li><p><code>vp[2]</code>: The sumation window is determined using U. Wolff poposal with <span>$S_\tau = wpm[2]$</span></p></li><li><p><code>vp[3]</code>: The autocorrelation function <span>$\Gamma(t)$</span> is summed up a point where its error <span>$\delta\Gamma(t)$</span> is a factor <code>vp[3]</code> times larger than the signal.</p></li></ul><p>An additional fourth parameter <code>vp[4]</code>, tells ADerrors to add a tail to the error with <span>$\tau_{exp} = wpm[4]$</span>. Negative values of <code>wpm[1:4]</code> are ignored and only one component of <code>wpm[1:3]</code> needs to be positive. If the flag <code>covar</code>is set to true, <code>fit_routine</code> takes into account covariances between x and y for each data point.</p><pre><code class="language-">@. model(x,p) = p[1] + p[2] * exp(-(p[3]-p[1])*x)
@. model2(x,p) = p[1] + p[2] * x[:, 1] + (p[3] + p[4] * x[:, 1]) * x[:, 2]
fit_routine(model, xdata, ydata, param=3)
fit_routine(model, xdata, ydata, param=3, covar=true)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="reader.html">« Reader</a><a class="docs-footer-nextpage" href="obs.html">Observables »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Wednesday 3 March 2021 11:34">Wednesday 3 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
fit_routine(model, xdata, ydata, param=3, covar=true)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://gitlab.ift.uam-csic.es/jugarrio/juobs">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="reader.html">« Reader</a><a class="docs-footer-nextpage" href="obs.html">Observables »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Tuesday 23 March 2021 11:06">Tuesday 23 March 2021</span>. Using Julia version 1.5.0.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>
# Observables
```@docs
meff
mpcac
dec_const
dec_const_pcvc
comp_t0
```
......
module juobs
using ADerrors, PyPlot, LaTeXStrings, LinearAlgebra, LsqFit, Optim
using ADerrors, PyPlot, LaTeXStrings, LinearAlgebra, LsqFit, LeastSquaresOptim
import Statistics: mean
include("juobs_types.jl")
......@@ -10,7 +10,7 @@ include("juobs_obs.jl")
export read_mesons, read_ms1, read_ms, read_md, truncate_data!
export get_matrix, energies, uwdot, uweigvals, uweigvecs, uweigen, invert, getall_eigvals, getall_eigvecs
export corr_obs, md_sea, md_val, plat_av, lin_fit, x_lin_fit, y_lin_fit, fit_routine
export meff, dec_const_pcvc, comp_t0
export corr_obs, corr_sym, md_sea, md_val, plat_av, lin_fit, x_lin_fit, y_lin_fit, fit_routine
export meff, mpcac, dec_const, dec_const_pcvc, comp_t0
end # module
@doc raw"""
meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false )
meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false)
meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
Computes effective mass for a given correlator corr at a given plateau `plat`.
Correlator can be passed as an `Corr` struct or `Vector{uwreal}`.
......@@ -14,12 +14,13 @@ corr_pp = corr_obs.(data)
m = meff(corr_pp[1], [50, 60], pl=false)
```
"""
function meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false, mu::Union{Vector{Float64}, Nothing}=nothing,
function meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bool=false,
mu::Union{Vector{Float64}, Nothing}=nothing, kappa::Union{Vector{Float64}, Nothing}=nothing,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dim = length(corr)
aux = 0.5 .* log.((corr[2:dim-2] ./ corr[3:dim-1]).^2)
mass = plat_av(aux, plat, wpm)
if pl == true
if pl
isnothing(wpm) ? uwerr(mass) : uwerr(mass, wpm)
x = 1:length(aux)
y = value.(aux)
......@@ -33,25 +34,284 @@ function meff(corr::Vector{uwreal}, plat::Vector{Int64}; pl::Bool=true, data::Bo
ylabel(L"$m_\mathrm{eff}$")
xlabel(L"$x_0$")
if !isnothing(kappa)
title(string(L"$\kappa_1 = $", kappa[1], L" $\kappa_2 = $", kappa[2]))
end
if !isnothing(mu)
title(string(L"$\mu_1 = $", mu[1], L" $\mu_2 = $", mu[2]))
end
display(gcf())
end
if !data
return mass
else
return (mass, aux)
end
end
function meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
if corr.mu == [0.0, 0.0]
return meff(corr.obs, plat, pl=pl, data=data, kappa=corr.kappa, mu=nothing, wpm=wpm)
else
return meff(corr.obs, plat, pl=pl, data=data, mu=corr.mu, kappa=nothing, wpm=wpm)
end
end
@doc raw"""
mpcac(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
mpcac(a0p::Corr, pp::Corr, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
Computes the bare PCAC mass for a given correlator `a0p` and `pp` at a given plateau `plat`.
Correlator can be passed as an `Corr` struct or `Vector{uwreal}`.
The flags `pl` and `data` allow to show the plots and return data as an extra result. The `ca` variable allows to compute `mpcac` using
the improved axial current.
```@example
data_pp = read_mesons(path, "G5", "G5")
data_a0p = read_mesons(path, "G5", "G0G5")
corr_pp = corr_obs.(data_pp)
corr_a0p = corr_obs.(data_a0p)
m12 = mpcac(corr_a0p, corr_pp, [50, 60], pl=false)
p0 = 9.2056
p1 = -13.9847
g2 = 1.73410
ca = -0.006033 * g2 *( 1 + exp(p0 + p1/g2))
m12 = mpcac(corr_a0p, corr_pp, [50, 60], pl=false, ca=ca)
```
"""
function mpcac(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false,
kappa::Union{Vector{Float64}, Nothing}=nothing, mu::Union{Vector{Float64}, Nothing}=nothing,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
corr_a0p = -a0p[2:end-1]
corr_pp = pp[2:end-1]
der_a0p = (corr_a0p[3:end] .- corr_a0p[1:end-2]) / 2
if ca != 0.0
der2_pp = corr_pp[1:end-2] + corr_pp[3:end] - 2 * corr_pp[2:end-1]
der_a0p = der_a0p + ca * der2_pp
end
aux = der_a0p ./ (2 .* corr_pp[2:end-1])
mass = plat_av(aux, plat, wpm)
if pl
isnothing(wpm) ? uwerr(mass) : uwerr(mass, wpm)
x = 1:length(aux)
y = value.(aux)
dy = err.(aux)
v = value(mass)
e = err(mass)
figure()
fill_between(plat[1]:plat[2], v-e, v+e, color="green", alpha=0.75)
errorbar(x, y, dy, fmt="x", color="black")
ylabel(L"$m_\mathrm{PCAC}$")
xlabel(L"$x_0$")
if !isnothing(kappa)
title(string(L"$\kappa_1 = $", kappa[1], L" $\kappa_2 = $", kappa[2]))
end
if !isnothing(mu)
title(string(L"$\mu_1 = $", mu[1], L" $\mu_2 = $", mu[2]))
end
display(gcf())
end
if data == false
if !data
return mass
else
return (mass, aux)
end
end
meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false) =
meff(corr.obs, plat, pl=pl, data=data, mu=corr.mu)
## Decay constants
function mpcac(a0p::Corr, pp::Corr, plat::Vector{Int64}; ca::Float64=0.0, pl::Bool=true, data::Bool=false,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
if a0p.kappa == pp.kappa || a0p.mu == pp.mu
if a0p.mu == [0.0, 0.0]
return mpcac(a0p.obs, pp.obs, plat, ca=ca, pl=pl, data=data, kappa=a0p.kappa, mu=nothing, wpm=wpm)
else
return mpcac(a0p.obs, pp.obs, plat, ca=ca, pl=pl, data=data, mu=a0p.mu, kappa=nothing, wpm=wpm)
end
else
error("mu or kappa values does not match")
end
end
## Decay constants
@doc raw"""
dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu::Vector{Float64}, y0::Int64 ; pl::Bool=true, data::Bool=false)meff(corr::Corr, plat::Vector{Int64}; pl::Bool=true, data::Bool=false)
dec_const(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, y0::Int64; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const(a0p::Corr, pp::Corr, plat::Vector{Int64}, m::uwreal; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const(a0pL::Vector{uwreal}, a0pR::Vector{uwreal}, ppL::Vector{uwreal}, ppR::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, y0::Int64; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const(a0pL::Corr, a0pR::Corr, ppL::Corr, ppR::Corr, plat::Vector{Int64}, m::uwreal; ca::Float64=0.0, pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
Computes the bare decay constant using ``A_0P`` and ``PP`` correlators . The decay constant is computed in the plateau `plat`.
Correlator can be passed as an `Corr` struct or `Vector{uwreal}`. If it is passed as a uwreal vector, effective mass `m` and source position `y0`
must be specified.
The flags `pl` and `data` allow to show the plots and return data as an extra result. The `ca` variable allows to compute `dec_const` using
the improved axial current.
**The method assumes that the source is close to the boundary.** It takes the following ratio to cancel boundary effects.
``R = \frac{f_A(x_0, y_0)}{\sqrt{f_P(T-y_0, y_0)}} * e^{m (x_0 - T/2)}``
**If left and right correlators are included in the input. The result is computed with the following ratio**
`` R = \sqrt{f_A(x_0, y_0) * f_A(x_0, T - 1 - y_0) / f_P(T - 1 - y_0, y_0)}``
```@example
data_pp = read_mesons(path, "G5", "G5", legacy=true)
data_a0p = read_mesons(path, "G5", "G0G5", legacy=true)
corr_pp = corr_obs.(data_pp, L=32)
corr_a0p = corr_obs.(data_a0p, L=32)
dec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false)
corr_a0pL, corr_a0pR = [corr_a0p[1], corr_a0p[2]]
corr_ppL, corr_ppR = [corr_pp[1], corr_pp[2]]
m = meff(corr_pp[1], [50, 60], pl=false)
beta = 3.46
p0 = 9.2056
p1 = -13.9847
g2 = 6 / beta
ca = -0.006033 * g2 *( 1 + exp(p0 + p1/g2))
f = dec_const(corr_a0p[1], corr_pp[1], [50, 60], m, pl=true, ca=ca)
f_ratio = dec_const(corr_a0pL, corr_a0pR, corr_ppL, corr_ppR, [50, 60], m, pl=true, ca=ca)
```
"""
function dec_const(a0p::Vector{uwreal}, pp::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, y0::Int64; ca::Float64=0.0, pl::Bool=true, data::Bool=false,
kappa::Union{Vector{Float64}, Nothing}=nothing, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
corr_a0p = -a0p
corr_pp = pp
T = length(corr_a0p)
if ca != 0.0
der_pp = (corr_pp[3:end] - corr_pp[1:end-2]) / 2
corr_a0p = corr_a0p[2:end-1] + ca * der_pp
aux = exp.((collect(1:T-2) .- (T-1)/2) .* [m for k = 1:T-2])
else
aux = exp.((collect(0:T-1) .- T/2) .* [m for k = 1:T])
end
R = corr_a0p .* aux ./ [sqrt(corr_pp[T-y0]) for k = 1:length(corr_a0p)]
R_av = plat_av(R, plat, wpm)
f = sqrt(2) * sqrt(R_av^2) / sqrt(m)
if pl
isnothing(wpm) ? uwerr(R_av) : uwerr(R_av, wpm)
isnothing(wpm) ? uwerr(f) : uwerr(f, wpm)
x = 1:length(R)
y = value.(R)
dy = err.(R)
v = value(R_av)
e = err(R_av)
figure()
lbl = string(L"$af = $", sprint(show, f))
fill_between(plat[1]:plat[2], v-e, v+e, color="green", alpha=0.75, label=L"$R$")
errorbar(x, y, dy, fmt="x", color="black", label=lbl)
legend()
ylabel(L"$R_\mathrm{av}$")
xlabel(L"$x_0$")
if !isnothing(kappa)
title(string(L"$\kappa_1 = $", kappa[1], L" $\kappa_2 = $", kappa[2]))
end
display(gcf())
end
if !data
return f
else
return (f, R)
end
end
function dec_const(a0p::Corr, pp::Corr, plat::Vector{Int64}, m::uwreal; ca::Float64=0.0, pl::Bool=true, data::Bool=false,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
if (a0p.y0 == pp.y0) && (a0p.kappa == pp.kappa)
return dec_const(a0p.obs, pp.obs, plat, m, a0p.y0, ca=ca, kappa=a0p.kappa, pl=pl, data=data, wpm=wpm)
else
error("y0 or kappa values does not match")
end
end
function dec_const(a0pL::Vector{uwreal}, a0pR::Vector{uwreal}, ppL::Vector{uwreal}, ppR::Vector{uwreal}, plat::Vector{Int64},
m::uwreal, y0::Int64; ca::Float64=0.0, pl::Bool=true, data::Bool=false, kappa::Union{Vector{Float64}, Nothing}=nothing,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
corr_pp = (ppL + ppR[end:-1:1]) / 2
T = length(corr_pp)
if ca != 0.0
der_ppL = (ppL[3:end] - ppL[1:end-2]) / 2
der_ppR = (ppR[3:end] - ppR[1:end-2]) / 2
corr_a0pL = -a0pL[2:end-1] + ca * der_ppL
corr_a0pR = -a0pR[2:end-1] + ca * der_ppR
else
corr_a0pL = -a0pL[2:end-1]
corr_a0pR = -a0pR[2:end-1]
end
f1 = [corr_pp[T - y0] for k = 1:length(corr_a0pL)]
R = ((corr_a0pL .* corr_a0pR ./ f1).^2).^(1/4)
R_av = plat_av(R, plat, wpm)
f = sqrt(2) * sqrt(R_av^2) / sqrt(m)
if pl
isnothing(wpm) ? uwerr(R_av) : uwerr(R_av, wpm)
isnothing(wpm) ? uwerr(f) : uwerr(f, wpm)
x = 1:length(R)
y = value.(R)
dy = err.(R)
v = value(R_av)
e = err(R_av)
figure()
lbl = string(L"$af = $", sprint(show, f))
fill_between(plat[1]:plat[2], v-e, v+e, color="green", alpha=0.75, label=L"$R$")
errorbar(x, y, dy, fmt="x", color="black", label=lbl)
legend()
ylabel(L"$R_\mathrm{av}$")
xlabel(L"$x_0$")
if !isnothing(kappa)
title(string(L"$\kappa_1 = $", kappa[1], L" $\kappa_2 = $", kappa[2]))
end
display(gcf())
end
if !data
return f
else
return (f, R)
end
end
function dec_const(a0pL::Corr, a0pR::Corr, ppL::Corr, ppR::Corr, plat::Vector{Int64}, m::uwreal; ca::Float64=0.0, pl::Bool=true, data::Bool=false,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
T = length(a0pL.obs)
if (a0pL.y0 == ppL.y0) && (a0pR.y0 == ppR.y0) && (a0pL.kappa == ppL.kappa) && (a0pR.kappa == ppR.kappa) && (a0pL.y0 == T - 1 - a0pR.y0)
return dec_const(a0pL.obs, a0pR.obs, ppL.obs, ppR.obs, plat, m, a0pL.y0, ca=ca, kappa=a0pL.kappa, pl=pl, data=data, wpm=wpm)
else
error("y0 or kappa values does not match")
end
end
@doc raw"""
dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu::Vector{Float64}, y0::Int64 ; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
dec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
Computes decay constant using the PCVC relation for twisted mass fermions. The decay constant is computed in the plateau `plat`.
Correlator can be passed as an `Corr` struct or `Vector{uwreal}`. If it is passed as a uwreal vector, vector of twisted masses `mu` and source position `y0`
......@@ -59,9 +319,10 @@ must be specified.
The flags `pl` and `data` allow to show the plots and return data as an extra result.
**The method assumes that the source is in the bulk.**
```@example
data = read_mesons(path, "G5", "G5")
corr_pp = corr_obs.(data)
corr_pp = corr_obs.(data, L=32)
m = meff(corr_pp[1], [50, 60], pl=false)
f = dec_const_pcvc(corr_pp[1], [50, 60], m, pl=false)
```
......@@ -78,7 +339,7 @@ function dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu
R = ((aux .* corr_pp).^2).^0.25
R_av = plat_av(R, plat, wpm)
f = sqrt(2) * (mu[1] + mu[2]) *R_av / m^1.5
if pl == true
if pl
if isnothing(wpm)
uwerr(f)
uwerr(R_av)
......@@ -99,14 +360,18 @@ function dec_const_pcvc(corr::Vector{uwreal}, plat::Vector{Int64}, m::uwreal, mu
title(string(L"$\mu_1 = $", mu[1], L" $\mu_2 = $", mu[2]))
display(gcf())
end
if data == false
if !data
return f
else
return (f, uwdot(sqrt(2) * (mu[1] + mu[2]) / m^1.5, R))
end
end
dec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false) =
dec_const_pcvc(corr.obs, plat, m, corr.mu, corr.y0, pl=pl, data=data)
function dec_const_pcvc(corr::Corr, plat::Vector{Int64}, m::uwreal; pl::Bool=true, data::Bool=false,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
return dec_const_pcvc(corr.obs, plat, m, corr.mu, corr.y0, pl=pl, data=data, wpm=wpm)
end
#t0
function get_model(x, p, n)
......@@ -117,9 +382,9 @@ function get_model(x, p, n)
return s
end
@doc raw"""
comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2)
comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Vector{Matrix{Float64}}, Nothing}=nothing, npol::Int64=2)
comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false, rw::Union{Vector{Matrix{Float64}}, Nothing}=nothing, npol::Int64=2, wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
Computes `t0` using the energy density of the action `Ysl`(Yang-Mills action).
`t0` is computed in the plateau `plat`.
......@@ -150,13 +415,13 @@ function comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false,
rw::Union{Matrix{Float64}, Nothing}=nothing, npol::Int64=2, ws::ADerrors.wspace=ADerrors.wsg,
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
Ysl = Y.Ysl
Ysl = Y.obs
t = Y.t
id = Y.id
replica = size.([Ysl], 1)
#Truncation
n_ws = findfirst(x-> x == id, ws.map_nob)
n_ws = findfirst(x-> x == ws.str2id[id], ws.map_nob)
if !isnothing(n_ws)
ivrep_ws = ws.fluc[n_ws].ivrep
......@@ -230,7 +495,7 @@ function comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false
wpm::Union{Dict{Int64,Vector{Float64}},Dict{String,Vector{Float64}}, Nothing}=nothing)
nr = length(Y)
Ysl = getfield.(Y, :Ysl)
Ysl = getfield.(Y, :obs)
t = getfield.(Y, :t)
t = t[1]
id = getfield.(Y, :id)
......@@ -240,7 +505,7 @@ function comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false
error("IDs are not equal")
end
#Truncation
n_ws = findfirst(x-> x == id[1], ws.map_nob)
n_ws = findfirst(x-> x == ws.str2id[id[1]], ws.map_nob)
if !isnothing(n_ws)
ivrep_ws = ws.fluc[n_ws].ivrep
......@@ -256,6 +521,7 @@ function comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false
error("Automatic truncation failed. R = ", replica[k], "\nTry using truncate_data!")
end
end
replica = size.(Ysl, 1)
end
Ysl = isnothing(rw) ? Ysl : apply_rw(Ysl, rw)
......@@ -312,6 +578,7 @@ function comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false
end
return t0
end
function t0_guess(t::Vector{Float64}, Ysl::Array{Float64, 3}, plat::Vector{Int64}, L::Int64)
t2E_ax = t.^2 .* mean(mean(Ysl[:, plat[1]:plat[2], :], dims=2), dims=1)[1, 1, :] / L^3
t0_aux = minimum(abs.(t2E_ax .- 0.3))
......
......@@ -81,13 +81,13 @@ function read_CHeader(path::String; legacy::Bool=false)
end
@doc raw"""
read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)
read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)
read_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)
read_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)
This function read a mesons dat file at a given path and returns a vector of `CData` structures for different masses and Dirac structures.
Dirac structures `g1` and/or `g2` can be passed as string arguments in order to filter correaltors.
ADerrors id can be specified as argument. If is not specified, the `id` is fixed according to the ensemble name (example: "H400"-> id = 400)
ADerrors id can be specified as argument. If is not specified, the `id` is fixed according to the ensemble name (example: "H400"-> id = "H400")
*For the old version (without smearing, distance preconditioning and theta) set legacy=true.
......@@ -97,17 +97,18 @@ read_mesons(path)
read_mesons(path, "G5")
read_mesons(path, nothing, "G5")
read_mesons(path, "G5", "G5")
read_mesons(path, "G5", "G5", id=1)
read_mesons(path, "G5", "G5", id="H100")
read_mesons(path, "G5_d2", "G5_d2", legacy=true)
```
"""
function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)
function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)
t1 = isnothing(g1) ? nothing : findfirst(x-> x==g1, gamma_name) - 1
t2 = isnothing(g2) ? nothing : findfirst(x-> x==g2, gamma_name) - 1
if isnothing(id)
bname = basename(path)
m = findfirst(r"[A-Z][0-9]{3}r[0-9]{3}", bname)
id = parse(Int64, bname[m[2:4]])
id = bname[m[1:4]]
#id = parse(Int64, bname[m[2:4]])
end
data = open(path, "r")
......@@ -120,8 +121,8 @@ function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union
fsize = filesize(path)
datsize = 4 + sum(c.dsize for c in c_header)*tvals*nnoise #datasize / ncnfg
ncfg = Int32((fsize-g_header.hsize-sum(c.hsize for c in c_header)) / datsize)
datsize = 4 + sum(getfield.(c_header, :dsize)) * tvals * nnoise #data_size / ncnfg
ncfg = div(fsize - g_header.hsize - sum(getfield.(c_header, :hsize)), datsize) #(total size - header_size) / data_size
corr_match = findall(x-> (x.type1==t1 || isnothing(t1)) && (x.type2==t2 || isnothing(t2)), c_header)
......@@ -171,7 +172,7 @@ function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union
return res
end
function read_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)
function read_mesons(path::Vector{String}, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{String, Nothing}=nothing, legacy::Bool=false)
res = read_mesons.(path, g1, g2, id=id, legacy=legacy)
nrep = length(res)
ncorr = length(res[1])
......@@ -274,28 +275,32 @@ function read_md(path::String)
end
@doc raw"""
read_ms(path::String; id::Union{Int64, Nothing}=nothing, dtr::Int64=1)
read_ms(path::String; id::Union{String, Nothing}=nothing, dtr::Int64=1, obs::String="Y")
Reads openQCD ms dat files at a given path. This method return YData:
- `t(t)`: flow time values
- `Ysl(icfg, x0, t)`: the time-slice sums of the densities of the Yang-Mills action
- `obs(icfg, x0, t)`: the time-slice sums of the densities of the observable (Wsl, Ysl or Qsl)
- `vtr`: vector that contains trajectory number
- `id`: ensmble id
`dtr` = `dtr_cnfg` / `dtr_ms`, where `dtr_cnfg` is the number of trajectories computed before saving the configuration. `dtr_ms`
is the same but applied to the ms.dat file.
Examples:
```@example
Y = read_ms(path)
```
"""
function read_ms(path::String; id::Union{Int64, Nothing}=nothing, dtr::Int64=1)
function read_ms(path::String; id::Union{String, Nothing}=nothing, dtr::Int64=1 , obs::String="Y")
if isnothing(id)
bname = basename(path)
m = findfirst(r"[A-Z][0-9]{3}r[0-9]{3}", bname)
id = parse(Int64, bname[m[2:4]])
id = bname[m[1:4]]
#id = parse(Int64, bname[m[2:4]])
end
data = open(path, "r")
......@@ -346,7 +351,17 @@ function read_ms(path::String; id::Union{Int64, Nothing}=nothing, dtr::Int64=1)
end
close(data)
t = Float64.(0:nn) .* dn .* eps
return YData(vntr, t, Ysl, id)
if obs == "W"
return YData(vntr, t, Wsl, id)
elseif obs == "Y"
return YData(vntr, t, Ysl, id)
elseif obs == "Q"
return YData(vntr, t, Qsl, id)
else
println("obs = ", obs," is not valid")
return nothing
end
end
@doc raw"""
......@@ -377,7 +392,7 @@ truncate_data!(Y, [nc1, nc2])
"""
function truncate_data!(data::YData, nc::Int64)
data.vtr = data.vtr[1:nc]
data.Ysl = data.Ysl[1:nc, :, :]
data.obs = data.obs[1:nc, :, :]
return nothing
end
function truncate_data!(data::Vector{YData}, nc::Vector{Int64})
......
......@@ -29,21 +29,6 @@ function apply_rw(data::Vector{<:Array{Float64}}, W::Vector{Matrix{Float64}})
return data_r
end
function dobsdp(a::uwreal, p::uwreal) # Compute da / dp
if count(p.prop .== true) != 1
error("I do not know how to compute this")
end
for i = 1:min(length(a.der), length(p.der))
if (p.prop[i] && a.prop[i])
return a.der[i] / p.der[i]
end
end
return 0.0
end
function check_corr_der(obs::Corr, derm::Vector{Corr})
g1 = Vector{String}(undef, 0)
g2 = Vector{String}(undef, 0)
......@@ -140,8 +125,36 @@ function corr_obs(cdata::Array{CData, 1}; real::Bool=true, rw::Union{Array{Array
return Corr(obs, cdata)
end
@doc raw"""
corr_sym(corrL::Corr, corrR::Corr, parity::Int64=1)
Computes the symmetrized correlator using the left correlador `corrL` and the right correlator `corrR`. The source position
of `corrR` must be `T - 1 - y0`, where `y0` is the source position of `corrL`.
```@example
pp_sym = corr_sym(ppL, ppR, +1)
a0p_sym = corr_sym(a0pL, a0pR, -1)
```
"""
function corr_sym(corrL::Corr, corrR::Corr, parity::Int64=1)
T = length(corrL.obs)
sym = [:kappa, :mu, :gamma]
if corrL.y0 != T - 1 - corrR.y0
error("Corr: Parameter mismatch")
end
for s in sym
if getfield(corrL, s) != getfield(corrR, s)
error("Corr: Parameter mismatch")
end
end
if abs(parity) != 1
error("incorrect value of parity (+- 1)")
end
#TODO: VECTORIZE
res = (corrL.obs[1:end] + parity * corrR.obs[end:-1:1]) / 2
return Corr(res, corrL.kappa, corrL.mu, corrL.gamma, corrL.y0)
end
#TODO: VECTORIZE, uwreal?
@doc raw"""
md_sea(a::uwreal, md::Vector{Matrix{Float64}}, ws::ADerrors.wspace=ADerrors.wsg)
......@@ -195,8 +208,8 @@ function md_sea(a::uwreal, md::Vector{Matrix{Float64}}, ws::ADerrors.wspace=ADer
if !all(id .== id[1])
error("ids do not match")
end
id = id[1]
id = ws.id2str[id[1]]
ivrep = getfield.(ws.fluc[p], :ivrep)
ivrep1 = fill(ivrep[1], length(ivrep))
if !all(ivrep .== ivrep1)
......@@ -208,6 +221,7 @@ function md_sea(a::uwreal, md::Vector{Matrix{Float64}}, ws::ADerrors.wspace=ADer
error("Nr obs != Nr md")
end
#md_aux as a Matrix + Automatic truncation
md_aux = md[1][:, 1:ivrep[1]]
for k = 2:length(md)
md_aux = cat(md_aux, md[k][:, 1:ivrep[k]], dims=2)
......@@ -215,27 +229,21 @@ function md_sea(a::uwreal, md::Vector{Matrix{Float64}}, ws::ADerrors.wspace=ADer
fluc_obs = getfield.(ws.fluc[p], :delta)
fluc_md = md_aux .- mean(md_aux, dims=2)
uwerr(a)
fluc_obs = mchist(a, id)
nrw = size(fluc_md, 1)
d = a.der[p]
nobs = sum(a.prop)
if nrw == 1
der1 = uwreal(0)
for k = 1:nobs
der1 = der1 - d[k] * uwreal(fluc_md[1, :] .* fluc_obs[k], id, ivrep)
end
der1 = uwreal(-fluc_md[1, :] .* fluc_obs, id, ivrep)
return (der1, der1)
elseif nrw == 2
der1 = uwreal(0)
der2 = uwreal(0)
for k = 1:nobs
der1 = der1 - d[k] * uwreal(fluc_md[1, :] .* fluc_obs[k], id, ivrep)
der2 = der2 - d[k] * uwreal(fluc_md[2, :] .* fluc_obs[k], id, ivrep)
end
der1 = uwreal(-fluc_md[1, :] .* fluc_obs, id, ivrep)
der2 = uwreal(-fluc_md[2, :] .* fluc_obs, id, ivrep)
return (der1, der2)
else
return nothing
end
end
@doc raw"""
......@@ -286,7 +294,7 @@ function md_val(a::uwreal, obs::Corr, derm::Vector{Corr})
corr = getfield(obs, :obs)
der = [dobsdp(a, corr[k]) for k = 1:length(corr)]
der = [derivative(a, corr[k]) for k = 1:length(corr)]
derm1, derm2 = derm
return (sum(der .* derm1.obs), sum(der .* derm2.obs))
end
......@@ -442,14 +450,18 @@ function fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}
C = isnothing(wpm) ? [ADerrors.cov(aux[k]) for k = 1:Ndata] : [ADerrors.cov(aux[k], wpm) for k = 1:Ndata]
chisq_full_cov(p, d) = get_chi2_cov(model, d, C, p, Nalpha)
min_fun_cov(t) = chisq_full_cov(t, dat)
sol = optimize(min_fun_cov, vcat(fit.param, dat[1:Nalpha*Ndata]), method=LBFGS())
sol = optimize(min_fun_cov, vcat(fit.param, dat[1:Nalpha*Ndata]), LevenbergMarquardt())
(upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full_cov, Optim.minimizer(sol), data) : fit_error(chisq_full_cov, Optim.minimizer(sol), data, wpm)
(upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full_cov, sol.minimizer, data) : fit_error(chisq_full_cov, sol.minimizer, data, wpm)
println("Chisq / chiexp: ", min_fun_cov(sol.minimizer), " / ", chi2_exp, " (dof: ", length(ydata) - param,")")
else
chisq_full(p, d) = get_chi2(model, d, ddat, p, Nalpha)
min_fun(t) = chisq_full(t, dat)
sol = optimize(min_fun, vcat(fit.param, dat[1:Nalpha*Ndata]), method=LBFGS())
(upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full, Optim.minimizer(sol), data) : fit_error(chisq_full, Optim.minimizer(sol), data, wpm)
sol = optimize(min_fun, vcat(fit.param, dat[1:Nalpha*Ndata]), LevenbergMarquardt())
(upar, chi2_exp) = isnothing(wpm) ? fit_error(chisq_full, sol.minimizer, data) : fit_error(chisq_full, sol.minimizer, data, wpm)
println("Chisq / chiexp: ", min_fun(sol.minimizer), " / ", chi2_exp, " (dof: ", length(ydata) - param,")")
end
#### chisq_full, min_fun out of conditional ->
......@@ -461,7 +473,6 @@ function fit_routine(model::Function, xdata::Array{uwreal}, ydata::Array{uwreal}
print("\n Fit parameter: ", i, ": ")
details(upar[i])
end
println("Chisq / chiexp: ", sol.minimum, " / ", chi2_exp, " (dof: ", length(ydata) - param,")")
return upar
end
......
......@@ -133,7 +133,7 @@ mutable struct CData
vcfg::Array{Int32}
re_data::Array{Float64}
im_data::Array{Float64}
id::Int64
id::String
CData(a, b, c, d, e) = new(a, b, c, d, e)
end
......@@ -167,13 +167,14 @@ mutable struct Corr
y0 = Int64(h[1].x0)
return new(a, kappa, mu, gamma, y0)
end
Corr(o::Vector{uwreal}, k::Vector{Float64}, m::Vector{Float64}, g::Vector{String}, y::Int64) = new(o, k, m, g, y)
end
mutable struct YData
vtr::Vector{Int32}
t::Vector{Float64}
Ysl::Array{Float64, 3}
id::Int64
obs::Array{Float64, 3}
id::String
YData(a, b, c, d) = new(a, b, c, d)
end
......@@ -214,4 +215,14 @@ end
function Base.show(io::IO, a::CData)
print(io, a.header)
end
function Base.show(io::IO, a::Corr)
fnames = fieldnames(Corr)
for k in fnames
f = getfield(a, k)
if k != :obs
print(io, "$k = $f\t")
end
end
end
\ No newline at end of file
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment