pyerrors/docs/pyerrors.html

1010 lines
100 KiB
HTML
Raw Normal View History

2021-11-07 20:53:18 +00:00
<!doctype html>
<html lang="en">
<head>
<meta charset="utf-8">
<meta name="viewport" content="width=device-width, initial-scale=1">
2022-02-15 12:25:34 +00:00
<meta name="generator" content="pdoc 10.0.1"/>
2021-11-07 20:53:18 +00:00
<title>pyerrors API documentation</title>
2022-02-14 14:08:21 +00:00
<style>/*! * Bootstrap Reboot v5.0.0 (https://getbootstrap.com/) * Copyright 2011-2021 The Bootstrap Authors * Copyright 2011-2021 Twitter, Inc. * Licensed under MIT (https://github.com/twbs/bootstrap/blob/main/LICENSE) * Forked from Normalize.css, licensed MIT (https://github.com/necolas/normalize.css/blob/master/LICENSE.md) */*,::after,::before{box-sizing:border-box}@media (prefers-reduced-motion:no-preference){:root{scroll-behavior:smooth}}body{margin:0;font-family:system-ui,-apple-system,"Segoe UI",Roboto,"Helvetica Neue",Arial,"Noto Sans","Liberation Sans",sans-serif,"Apple Color Emoji","Segoe UI Emoji","Segoe UI Symbol","Noto Color Emoji";font-size:1rem;font-weight:400;line-height:1.5;color:#212529;background-color:#fff;-webkit-text-size-adjust:100%;-webkit-tap-highlight-color:transparent}hr{margin:1rem 0;color:inherit;background-color:currentColor;border:0;opacity:.25}hr:not([size]){height:1px}h1,h2,h3,h4,h5,h6{margin-top:0;margin-bottom:.5rem;font-weight:500;line-height:1.2}h1{font-size:calc(1.375rem + 1.5vw)}@media (min-width:1200px){h1{font-size:2.5rem}}h2{font-size:calc(1.325rem + .9vw)}@media (min-width:1200px){h2{font-size:2rem}}h3{font-size:calc(1.3rem + .6vw)}@media (min-width:1200px){h3{font-size:1.75rem}}h4{font-size:calc(1.275rem + .3vw)}@media (min-width:1200px){h4{font-size:1.5rem}}h5{font-size:1.25rem}h6{font-size:1rem}p{margin-top:0;margin-bottom:1rem}abbr[data-bs-original-title],abbr[title]{-webkit-text-decoration:underline dotted;text-decoration:underline dotted;cursor:help;-webkit-text-decoration-skip-ink:none;text-decoration-skip-ink:none}address{margin-bottom:1rem;font-style:normal;line-height:inherit}ol,ul{padding-left:2rem}dl,ol,ul{margin-top:0;margin-bottom:1rem}ol ol,ol ul,ul ol,ul ul{margin-bottom:0}dt{font-weight:700}dd{margin-bottom:.5rem;margin-left:0}blockquote{margin:0 0 1rem}b,strong{font-weight:bolder}small{font-size:.875em}mark{padding:.2em;background-color:#fcf8e3}sub,sup{position:relative;font-size:.75em;line-height:0;vertical-align:baseline}sub{bottom:-.25em}sup{top:-.5em}a{color:#0d6efd;text-decoration:underline}a:hover{color:#0a58ca}a:not([href]):not([class]),a:not([href]):not([class]):hover{color:inherit;text-decoration:none}code,kbd,pre,samp{font-family:SFMono-Regular,Menlo,Monaco,Consolas,"Liberation Mono","Courier New",monospace;font-size:1em;direction:ltr;unicode-bidi:bidi-override}pre{display:block;margin-top:0;margin-bottom:1rem;overflow:auto;font-size:.875em}pre code{font-size:inherit;color:inherit;word-break:normal}code{font-size:.875em;color:#d63384;word-wrap:break-word}a>code{color:inherit}kbd{padding:.2rem .4rem;font-size:.875em;color:#fff;background-color:#212529;border-radius:.2rem}kbd kbd{padding:0;font-size:1em;font-weight:700}figure{margin:0 0 1rem}img,svg{vertical-align:middle}table{caption-side:bottom;border-collapse:collapse}caption{padding-top:.5rem;padding-bottom:.5rem;color:#6c757d;text-align:left}th{text-align:inherit;text-align:-webkit-match-parent}tbody,td,tfoot,th,thead,tr{border-color:inherit;border-style:solid;border-width:0}label{display:inline-block}button{border-radius:0}button:focus:not(:focus-visible){outline:0}button,input,optgroup,select,textarea{margin:0;font-family:inherit;font-size:inherit;line-height:inherit}button,select{text-transform:none}[role=button]{cursor:pointer}select{word-wrap:normal}select:disabled{opacity:1}[list]::-webkit-calendar-picker-indicator{display:none}[type=button],[type=reset],[type=submit],button{-webkit-appearance:button}[type=button]:not(:disabled),[type=reset]:not(:disabled),[type=submit]:not(:disabled),button:not(:disabled){cursor:pointer}::-moz-focus-inner{padding:0;border-style:none}textarea{resize:vertical}fieldset{min-width:0;padding:0;margin:0;border:0}legend{float:left;width:100%;padding:0;margin-bottom:.5rem;font-size:calc(1.275rem + .3vw);line-height:inherit}@media (min-width:1200px){legend{font-size:1.5rem}}legend+*{clear:left}::-webkit-datetime-edit-day-field,::-webkit-datetime-edit-fields-wrapper,::-webkit-datetime-edit-hour-field,::-webkit-datetime-edit-minute,::-webkit-datetime-edit-month-fie
<style>/*! syntax-highlighting.css */pre{line-height:125%;}span.linenos{color:inherit; background-color:transparent; padding-left:5px; padding-right:20px;}.pdoc-code .hll{background-color:#ffffcc}.pdoc-code{background:#f8f8f8;}.pdoc-code .c{color:#3D7B7B; font-style:italic}.pdoc-code .err{border:1px solid #FF0000}.pdoc-code .k{color:#008000; font-weight:bold}.pdoc-code .o{color:#666666}.pdoc-code .ch{color:#3D7B7B; font-style:italic}.pdoc-code .cm{color:#3D7B7B; font-style:italic}.pdoc-code .cp{color:#9C6500}.pdoc-code .cpf{color:#3D7B7B; font-style:italic}.pdoc-code .c1{color:#3D7B7B; font-style:italic}.pdoc-code .cs{color:#3D7B7B; font-style:italic}.pdoc-code .gd{color:#A00000}.pdoc-code .ge{font-style:italic}.pdoc-code .gr{color:#E40000}.pdoc-code .gh{color:#000080; font-weight:bold}.pdoc-code .gi{color:#008400}.pdoc-code .go{color:#717171}.pdoc-code .gp{color:#000080; font-weight:bold}.pdoc-code .gs{font-weight:bold}.pdoc-code .gu{color:#800080; font-weight:bold}.pdoc-code .gt{color:#0044DD}.pdoc-code .kc{color:#008000; font-weight:bold}.pdoc-code .kd{color:#008000; font-weight:bold}.pdoc-code .kn{color:#008000; font-weight:bold}.pdoc-code .kp{color:#008000}.pdoc-code .kr{color:#008000; font-weight:bold}.pdoc-code .kt{color:#B00040}.pdoc-code .m{color:#666666}.pdoc-code .s{color:#BA2121}.pdoc-code .na{color:#687822}.pdoc-code .nb{color:#008000}.pdoc-code .nc{color:#0000FF; font-weight:bold}.pdoc-code .no{color:#880000}.pdoc-code .nd{color:#AA22FF}.pdoc-code .ni{color:#717171; font-weight:bold}.pdoc-code .ne{color:#CB3F38; font-weight:bold}.pdoc-code .nf{color:#0000FF}.pdoc-code .nl{color:#767600}.pdoc-code .nn{color:#0000FF; font-weight:bold}.pdoc-code .nt{color:#008000; font-weight:bold}.pdoc-code .nv{color:#19177C}.pdoc-code .ow{color:#AA22FF; font-weight:bold}.pdoc-code .w{color:#bbbbbb}.pdoc-code .mb{color:#666666}.pdoc-code .mf{color:#666666}.pdoc-code .mh{color:#666666}.pdoc-code .mi{color:#666666}.pdoc-code .mo{color:#666666}.pdoc-code .sa{color:#BA2121}.pdoc-code .sb{color:#BA2121}.pdoc-code .sc{color:#BA2121}.pdoc-code .dl{color:#BA2121}.pdoc-code .sd{color:#BA2121; font-style:italic}.pdoc-code .s2{color:#BA2121}.pdoc-code .se{color:#AA5D1F; font-weight:bold}.pdoc-code .sh{color:#BA2121}.pdoc-code .si{color:#A45A77; font-weight:bold}.pdoc-code .sx{color:#008000}.pdoc-code .sr{color:#A45A77}.pdoc-code .s1{color:#BA2121}.pdoc-code .ss{color:#19177C}.pdoc-code .bp{color:#008000}.pdoc-code .fm{color:#0000FF}.pdoc-code .vc{color:#19177C}.pdoc-code .vg{color:#19177C}.pdoc-code .vi{color:#19177C}.pdoc-code .vm{color:#19177C}.pdoc-code .il{color:#666666}</style>
<style>/*! theme.css */:root{--pdoc-background:#fff;}.pdoc{--text:#212529;--muted:#6c757d;--link:#3660a5;--link-hover:#1659c5;--code:#f8f8f8;--active:#fff598;--accent:#eee;--accent2:#c1c1c1;--nav-hover:rgba(255, 255, 255, 0.5);--name:#0066BB;--def:#008800;--annotation:#007020;}</style>
<style>/*! layout.css */html, body{width:100%;height:100%;}html, main{scroll-behavior:smooth;}body{background-color:var(--pdoc-background);}@media (max-width:769px){#navtoggle{cursor:pointer;position:absolute;width:50px;height:40px;top:1rem;right:1rem;border-color:var(--text);color:var(--text);display:flex;opacity:0.8;}#navtoggle:hover{opacity:1;}#togglestate + div{display:none;}#togglestate:checked + div{display:inherit;}main, header{padding:2rem 3vw;}header + main{margin-top:-3rem;}.git-button{display:none !important;}nav input[type="search"]{max-width:77%;}nav input[type="search"]:first-child{margin-top:-6px;}nav input[type="search"]:valid ~ *{display:none !important;}}@media (min-width:770px){:root{--sidebar-width:clamp(12.5rem, 28vw, 22rem);}nav{position:fixed;overflow:auto;height:100vh;width:var(--sidebar-width);}main, header{padding:3rem 2rem 3rem calc(var(--sidebar-width) + 3rem);width:calc(54rem + var(--sidebar-width));max-width:100%;}header + main{margin-top:-4rem;}#navtoggle{display:none;}}#togglestate{position:absolute;height:0;}nav.pdoc{--pad:1.75rem;--indent:1.5rem;background-color:var(--accent);border-right:1px solid var(--accent2);box-shadow:0 0 20px rgba(50, 50, 50, .2) inset;padding:0 0 0 var(--pad);overflow-wrap:anywhere;scrollbar-width:thin; scrollbar-color:var(--accent2) transparent }nav.pdoc::-webkit-scrollbar{width:.4rem; }nav.pdoc::-webkit-scrollbar-thumb{background-color:var(--accent2); }nav.pdoc > div{padding:var(--pad) 0;}nav.pdoc .module-list-button{display:inline-flex;align-items:center;color:var(--text);border-color:var(--muted);margin-bottom:1rem;}nav.pdoc .module-list-button:hover{border-color:var(--text);}nav.pdoc input[type=search]{display:block;outline-offset:0;width:calc(100% - var(--pad));}nav.pdoc .logo{max-width:calc(100% - var(--pad));max-height:35vh;display:block;margin:0 auto 1rem;transform:translate(calc(-.5 * var(--pad)), 0);}nav.pdoc ul{list-style:none;padding-left:0;}nav.pdoc li{display:block;margin:0;padding:.2rem 0 .2rem var(--indent);transition:all 100ms;}nav.pdoc > div > ul > li{padding-left:0;}nav.pdoc li:hover{background-color:var(--nav-hover);}nav.pdoc a, nav.pdoc a:hover{color:var(--text);}nav.pdoc a{display:block;}nav.pdoc > h2:first-of-type{margin-top:1.5rem;}nav.pdoc .class:before{content:"class ";color:var(--muted);}nav.pdoc .function:after{content:"()";color:var(--muted);}nav.pdoc footer:before{content:"";display:block;width:calc(100% - var(--pad));border-top:solid var(--accent2) 1px;margin-top:1.5rem;padding-top:.5rem;}nav.pdoc footer{font-size:small;}</style>
2022-02-15 12:25:34 +00:00
<style>/*! content.css */.pdoc{color:var(--text);box-sizing:border-box;line-height:1.5;background:none;}.pdoc .pdoc-button{display:inline-block;border:solid black 1px;border-radius:2px;font-size:.75rem;padding:calc(0.5em - 1px) 1em;transition:100ms all;}.pdoc .visually-hidden{position:absolute !important;width:1px !important;height:1px !important;padding:0 !important;margin:-1px !important;overflow:hidden !important;clip:rect(0, 0, 0, 0) !important;white-space:nowrap !important;border:0 !important;}.pdoc h1, .pdoc h2, .pdoc h3{font-weight:300;margin:.3em 0;padding:.2em 0;}.pdoc a{text-decoration:none;color:var(--link);}.pdoc a:hover{color:var(--link-hover);}.pdoc blockquote{margin-left:2rem;}.pdoc pre{border-top:1px solid var(--accent2);border-bottom:1px solid var(--accent2);margin-top:0;margin-bottom:1em;padding:.5rem 0 .5rem .5rem;overflow-x:auto;background-color:var(--code);}.pdoc code{color:var(--text);padding:.2em .4em;margin:0;font-size:85%;background-color:var(--code);border-radius:6px;}.pdoc a > code{color:inherit;}.pdoc pre > code{display:inline-block;font-size:inherit;background:none;border:none;padding:0;}.pdoc .modulename{margin-top:0;font-weight:bold;}.pdoc .modulename a{color:var(--link);transition:100ms all;}.pdoc .git-button{float:right;border:solid var(--link) 1px;}.pdoc .git-button:hover{background-color:var(--link);color:var(--pdoc-background);}.pdoc details{filter:opacity(1);}.pdoc details:not([open]){height:0;}.pdoc details > summary{position:absolute;top:-35px;right:0;font-size:.75rem;color:var(--muted);padding:0 .7em;user-select:none;cursor:pointer;}.pdoc details > summary:focus{outline:0;}.pdoc .docstring{margin-bottom:1.5rem;}.pdoc > section:first-of-type > .docstring{margin-bottom:2.5rem;}.pdoc .docstring .pdoc-code{margin-left:1em;margin-right:1em;}.pdoc h1:target,.pdoc h2:target,.pdoc h3:target,.pdoc h4:target,.pdoc h5:target,.pdoc h6:target{background-color:var(--active);box-shadow:-1rem 0 0 0 var(--active);}.pdoc div:target > .attr,.pdoc section:target > .attr,.pdoc dd:target > a{background-color:var(--active);}.pdoc .attr:hover{filter:contrast(0.95);}.pdoc .headerlink{position:absolute;width:0;margin-left:-1.5rem;line-height:1.4rem;font-size:1.5rem;font-weight:normal;transition:all 100ms ease-in-out;opacity:0;user-select:none;}.pdoc .attr > .headerlink{margin-left:-2.5rem;}.pdoc *:hover > .headerlink,.pdoc *:target > .attr > .headerlink{opacity:1;}.pdoc .attr{display:block;color:var(--text);margin:.5rem 0 .5rem;padding:.4rem 5rem .4rem 1rem;background-color:var(--accent);}.pdoc .classattr{margin-left:2rem;}.pdoc .name{color:var(--name);font-weight:bold;}.pdoc .def{color:var(--def);font-weight:bold;}.pdoc .signature{white-space:pre-wrap;}.pdoc .annotation{color:var(--annotation);}.pdoc .inherited{margin-left:2rem;}.pdoc .inherited dt{font-weight:700;}.pdoc .inherited dt, .pdoc .inherited dd{display:inline;margin-left:0;margin-bottom:.5rem;}.pdoc .inherited dd:not(:last-child):after{content:", ";}.pdoc .inherited .class:before{content:"class ";}.pdoc .inherited .function a:after{content:"()";}.pdoc .search-result .docstring{overflow:auto;max-height:25vh;}.pdoc .search-result.focused > .attr{background-color:var(--active);}.pdoc .attribution{margin-top:2rem;display:block;opacity:0.5;transition:all 200ms;filter:grayscale(100%);}.pdoc .attribution:hover{opacity:1;filter:grayscale(0%);}.pdoc .attribution img{margin-left:5px;height:35px;vertical-align:middle;width:70px;transition:all 200ms;}.pdoc table{display:block;width:max-content;max-width:100%;overflow:auto;margin-bottom:1rem;}.pdoc table th{font-weight:600;}.pdoc table th, .pdoc table td{padding:6px 13px;border:1px solid var(--accent2);}</style>
<style>/*! custom.css */</style><script>
2021-11-07 20:53:18 +00:00
window.MathJax = {
tex: {
inlineMath: [['$', '$'], ['\\(', '\\)']]
}
};
</script>
<script src="https://polyfill.io/v3/polyfill.min.js?features=es6"></script>
<script id="MathJax-script" async src="https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-mml-chtml.js"></script>
<script>
/* Re-invoke MathJax when DOM content changes, for example during search. */
document.addEventListener("DOMContentLoaded", () => {
new MutationObserver(() => MathJax.typeset()).observe(
document.querySelector("main.pdoc").parentNode,
{childList: true}
);
})
2022-02-14 14:08:21 +00:00
</script></head>
<body>
<nav class="pdoc">
<label id="navtoggle" for="togglestate" class="pdoc-button"><svg xmlns='http://www.w3.org/2000/svg' viewBox='0 0 30 30'><path stroke-linecap='round' stroke="currentColor" stroke-miterlimit='10' stroke-width='2' d='M4 7h22M4 15h22M4 23h22'/></svg></label>
<input id="togglestate" type="checkbox" aria-hidden="true" tabindex="-1">
<div>
<input type="search" placeholder="Search..." role="searchbox" aria-label="search"
pattern=".+" required>
<h2>Contents</h2>
<ul>
2021-11-07 20:53:18 +00:00
<li><a href="#what-is-pyerrors">What is pyerrors?</a>
<ul>
2021-11-16 11:46:28 +00:00
<li><a href="#basic-example">Basic example</a></li>
2021-11-07 20:53:18 +00:00
</ul></li>
<li><a href="#the-obs-class">The <code>Obs</code> class</a>
<ul>
2021-11-09 11:51:27 +00:00
<li><a href="#error-propagation">Error propagation</a></li>
2022-01-16 15:44:50 +00:00
<li><a href="#error-estimation">Error estimation</a></li>
<li><a href="#multiple-ensemblesreplica">Multiple ensembles/replica</a></li>
2021-11-07 20:53:18 +00:00
<li><a href="#irregular-monte-carlo-chains">Irregular Monte Carlo chains</a></li>
</ul></li>
2021-11-07 21:09:48 +00:00
<li><a href="#correlators">Correlators</a></li>
2022-02-14 14:08:21 +00:00
<li><a href="#complex-valued-observables">Complex valued observables</a></li>
<li><a href="#error-propagation-in-iterative-algorithms">Error propagation in iterative algorithms</a>
<ul>
<li><a href="#least-squares-fits">Least squares fits</a></li>
<li><a href="#total-least-squares-fits">Total least squares fits</a></li>
</ul></li>
2021-11-07 20:53:18 +00:00
<li><a href="#matrix-operations">Matrix operations</a></li>
2021-11-15 14:30:41 +00:00
<li><a href="#export-data">Export data</a>
<ul>
<li><a href="#jackknife-samples">Jackknife samples</a></li>
</ul></li>
2022-02-14 14:08:21 +00:00
<li><a href="#citing">Citing</a></li>
2021-11-07 20:53:18 +00:00
</ul>
2022-02-14 14:08:21 +00:00
<h2>Submodules</h2>
<ul>
<li><a href="pyerrors/correlators.html">pyerrors.correlators</a></li>
<li><a href="pyerrors/covobs.html">pyerrors.covobs</a></li>
<li><a href="pyerrors/dirac.html">pyerrors.dirac</a></li>
<li><a href="pyerrors/fits.html">pyerrors.fits</a></li>
<li><a href="pyerrors/input.html">pyerrors.input</a></li>
<li><a href="pyerrors/linalg.html">pyerrors.linalg</a></li>
<li><a href="pyerrors/misc.html">pyerrors.misc</a></li>
<li><a href="pyerrors/mpm.html">pyerrors.mpm</a></li>
<li><a href="pyerrors/obs.html">pyerrors.obs</a></li>
<li><a href="pyerrors/roots.html">pyerrors.roots</a></li>
<li><a href="pyerrors/version.html">pyerrors.version</a></li>
</ul>
<a class="attribution" title="pdoc: Python API documentation generator" href="https://pdoc.dev">
built with <span class="visually-hidden">pdoc</span><img
alt="pdoc logo"
src="data:image/svg+xml,%3Csvg%20xmlns%3D%22http%3A//www.w3.org/2000/svg%22%20role%3D%22img%22%20aria-label%3D%22pdoc%20logo%22%20width%3D%22300%22%20height%3D%22150%22%20viewBox%3D%22-1%200%2060%2030%22%3E%3Ctitle%3Epdoc%3C/title%3E%3Cpath%20d%3D%22M29.621%2021.293c-.011-.273-.214-.475-.511-.481a.5.5%200%200%200-.489.503l-.044%201.393c-.097.551-.695%201.215-1.566%201.704-.577.428-1.306.486-2.193.182-1.426-.617-2.467-1.654-3.304-2.487l-.173-.172a3.43%203.43%200%200%200-.365-.306.49.49%200%200%200-.286-.196c-1.718-1.06-4.931-1.47-7.353.191l-.219.15c-1.707%201.187-3.413%202.131-4.328%201.03-.02-.027-.49-.685-.141-1.763.233-.721.546-2.408.772-4.076.042-.09.067-.187.046-.288.166-1.347.277-2.625.241-3.351%201.378-1.008%202.271-2.586%202.271-4.362%200-.976-.272-1.935-.788-2.774-.057-.094-.122-.18-.184-.268.033-.167.052-.339.052-.516%200-1.477-1.202-2.679-2.679-2.679-.791%200-1.496.352-1.987.9a6.3%206.3%200%200%200-1.001.029c-.492-.564-1.207-.929-2.012-.929-1.477%200-2.679%201.202-2.679%202.679A2.65%202.65%200%200%200%20.97%206.554c-.383.747-.595%201.572-.595%202.41%200%202.311%201.507%204.29%203.635%205.107-.037.699-.147%202.27-.423%203.294l-.137.461c-.622%202.042-2.515%208.257%201.727%2010.643%201.614.908%203.06%201.248%204.317%201.248%202.665%200%204.492-1.524%205.322-2.401%201.476-1.559%202.886-1.854%206.491.82%201.877%201.393%203.514%201.753%204.861%201.068%202.223-1.713%202.811-3.867%203.399-6.374.077-.846.056-1.469.054-1.537zm-4.835%204.313c-.054.305-.156.586-.242.629-.034-.007-.131-.022-.307-.157-.145-.111-.314-.478-.456-.908.221.121.432.25.675.355.115.039.219.051.33.081zm-2.251-1.238c-.05.33-.158.648-.252.694-.022.001-.125-.018-.307-.157-.217-.166-.488-.906-.639-1.573.358.344.754.693%201.198%201.036zm-3.887-2.337c-.006-.116-.018-.231-.041-.342.635.145%201.189.368%201.599.625.097.231.166.481.174.642-.03.049-.055.101-.067.158-.046.013-.128.026-.298.004-.278-.037-.901-.57-1.367-1.087zm-1.127-.497c.116.306.176.625.12.71-.019.014-.117.045-.345.016-.206-.027-.604-.332-.986-.695.41-.051.816-.056%201.211-.031zm-4.535%201.535c.209.22.379.47.358.598-.006.041-.088.138-.351.234-.144.055-.539-.063-.979-.259a11.66%2011.66%200%200%200%20.972-.573zm.983-.664c.359-.237.738-.418%201.126-.554.25.237.479.548.457.694-.006.042-.087.138-.351.235-.174.064-.694-.105-1.232-.375zm-3.381%201.794c-.022.145-.061.29-.149.401-.133.166-.358.248-.69.251h-.002c-.133%200-.306-.26-.45-.621.417.091.854.07%201.291-.031zm-2.066-8.077a4.78%204.78%200%200%201-.775-.584c.172-.115.505-.254.88-.378l-.105.962zm-.331%202.302a10.32%2010.32%200%200%201-.828-.502c.202-.143.576-.328.984-.49l-.156.992zm-.45%202.157l-.701-.403c.214-.115.536-.249.891-.376a11.57%2011.57%200%200%201-.19.779zm-.181%201.716c.064.398.194.702.298.893-.194-.051-.435-.162-.736-.398.061-.119.224-.3.438-.495zM8.87%204.141c0%20.152-.123.276-.276.276s-.275-.124-.275-.276.123-.276.276-.276.275.124.275.276zm-.735-.389a1.15%201.15%200%200%200-.314.783%201.16%201.16%200%200%200%201.162%201.162c.457%200%20.842-.27%201.032-.653.026.117.042.238.042.362a1.68%201.68%200%200%201-1.679%201.679%201.68%201.68%200%200%201-1.679-1.679c0-.843.626-1.535%201.436-1.654zM5.059%205.406A1.68%201.68%200%200%201%203.38%207.085a1.68%201.68%200%200%201-1.679-1.679c0-.037.009-.072.011-.109.21.3.541.508.935.508a1.16%201.16%200%200%200%201.162-1.162%201.14%201.14%200%200%200-.474-.912c.015%200%20.03-.005.045-.005.926.001%201.679.754%201.679%201.68zM3.198%204.141c0%20.152-.123.276-.276.276s-.275-.124-.275-.276.123-.276.276-.276.275.124.275.276zM1.375%208.964c0-.52.103-1.035.288-1.52.466.394%201.06.64%201.717.64%201.144%200%202.116-.725%202.499-1.738.383%201.012%201.355%201.738%202.499%201.738.867%200%201.631-.421%202.121-1.062.307.605.478%201.267.478%201.942%200%202.486-2.153%204.51-4.801%204.51s-4.801-2.023-4.801-4.51zm24.342%2019.349c-.985.498-2.267.168-3.813-.979-3.073-2.281-5.453-3.199-7.813-.705-1.315%201.391-4.163%203.365-8.423.97-3.174-1.786-2.239-6.266-1.261-9.479l.146-.492c.276-1.02.395-2.457.444-3.268a6.11%206.11%200%200%200%201.18.115%206.01%206.01%200%200%200%202.536-.562l-.006.175c-.802.215-1.848.612
</a>
</div>
</nav>
2021-11-07 20:53:18 +00:00
<main class="pdoc">
<section>
<h1 class="modulename">
pyerrors </h1>
<div class="docstring"><h1 id="what-is-pyerrors">What is pyerrors?</h1>
2021-11-07 21:04:06 +00:00
<p><code><a href="">pyerrors</a></code> is a python package for error computation and propagation of Markov chain Monte Carlo data.
2021-12-11 23:00:30 +00:00
It is based on the gamma method <a href="https://arxiv.org/abs/hep-lat/0306017">arXiv:hep-lat/0306017</a>. Some of its features are:</p>
2021-11-07 21:04:06 +00:00
<ul>
2021-12-11 23:00:30 +00:00
<li>automatic differentiation for exact liner error propagation as suggested in <a href="https://arxiv.org/abs/1809.01289">arXiv:1809.01289</a> (partly based on the <a href="https://github.com/HIPS/autograd">autograd</a> package).</li>
<li>treatment of slow modes in the simulation as suggested in <a href="https://arxiv.org/abs/1009.5228">arXiv:1009.5228</a>.</li>
<li>coherent error propagation for data from different Markov chains.</li>
<li>non-linear fits with x- and y-errors and exact linear error propagation based on automatic differentiation as introduced in <a href="https://arxiv.org/abs/1809.01289">arXiv:1809.01289</a>.</li>
<li>real and complex matrix operations and their error propagation based on automatic differentiation (Matrix inverse, Cholesky decomposition, calculation of eigenvalues and eigenvectors, singular value decomposition...).</li>
2021-11-07 21:04:06 +00:00
</ul>
2021-11-07 20:53:18 +00:00
2021-12-11 23:00:30 +00:00
<p>There exist similar publicly available implementations of gamma method error analysis suites in <a href="https://gitlab.ift.uam-csic.es/alberto/aderrors">Fortran</a>, <a href="https://gitlab.ift.uam-csic.es/alberto/aderrors.jl">Julia</a> and <a href="https://github.com/mbruno46/pyobs">Python</a>.</p>
2021-11-30 14:58:46 +00:00
2021-11-16 11:46:28 +00:00
<h2 id="basic-example">Basic example</h2>
2021-11-07 20:53:18 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="kn">import</span> <span class="nn">numpy</span> <span class="k">as</span> <span class="nn">np</span>
2021-11-07 20:53:18 +00:00
<span class="kn">import</span> <span class="nn">pyerrors</span> <span class="k">as</span> <span class="nn">pe</span>
2021-11-16 11:50:48 +00:00
<span class="n">my_obs</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble_name&#39;</span><span class="p">])</span> <span class="c1"># Initialize an Obs object</span>
2021-11-16 11:46:28 +00:00
<span class="n">my_new_obs</span> <span class="o">=</span> <span class="mi">2</span> <span class="o">*</span> <span class="n">np</span><span class="o">.</span><span class="n">log</span><span class="p">(</span><span class="n">my_obs</span><span class="p">)</span> <span class="o">/</span> <span class="n">my_obs</span> <span class="o">**</span> <span class="mi">2</span> <span class="c1"># Construct derived Obs object</span>
2021-11-16 11:50:48 +00:00
<span class="n">my_new_obs</span><span class="o">.</span><span class="n">gamma_method</span><span class="p">()</span> <span class="c1"># Estimate the statistical error</span>
<span class="nb">print</span><span class="p">(</span><span class="n">my_new_obs</span><span class="p">)</span> <span class="c1"># Print the result to stdout</span>
2021-11-15 13:13:15 +00:00
<span class="o">&gt;</span> <span class="mf">0.31498</span><span class="p">(</span><span class="mi">72</span><span class="p">)</span>
2021-11-07 20:53:18 +00:00
</code></pre></div>
<h1 id="the-obs-class">The <code>Obs</code> class</h1>
2021-11-09 11:51:27 +00:00
<p><code><a href="">pyerrors</a></code> introduces a new datatype, <code>Obs</code>, which simplifies error propagation and estimation for auto- and cross-correlated data.
2022-02-21 14:52:31 +00:00
An <code>Obs</code> object can be initialized with two arguments, the first is a list containing the samples for an observable from a Monte Carlo chain.
2021-11-09 11:51:27 +00:00
The samples can either be provided as python list or as numpy array.
The second argument is a list containing the names of the respective Monte Carlo chains as strings. These strings uniquely identify a Monte Carlo chain/ensemble.</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="kn">import</span> <span class="nn">pyerrors</span> <span class="k">as</span> <span class="nn">pe</span>
2021-11-07 20:53:18 +00:00
<span class="n">my_obs</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble_name&#39;</span><span class="p">])</span>
</code></pre></div>
2021-11-09 11:51:27 +00:00
<h2 id="error-propagation">Error propagation</h2>
2022-02-21 14:52:31 +00:00
<p>When performing mathematical operations on <code>Obs</code> objects the correct error propagation is intrinsically taken care of using a first order Taylor expansion
2021-11-16 11:46:28 +00:00
$$\delta_f^i=\sum_\alpha \bar{f}_\alpha \delta_\alpha^i\,,\quad \delta_\alpha^i=a_\alpha^i-\bar{a}_\alpha\,,$$
2021-11-16 11:50:48 +00:00
as introduced in <a href="https://arxiv.org/abs/hep-lat/0306017">arXiv:hep-lat/0306017</a>.
The required derivatives $\bar{f}_\alpha$ are evaluated up to machine precision via automatic differentiation as suggested in <a href="https://arxiv.org/abs/1809.01289">arXiv:1809.01289</a>.</p>
2021-11-09 11:51:27 +00:00
<p>The <code>Obs</code> class is designed such that mathematical numpy functions can be used on <code>Obs</code> just as for regular floats.</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="kn">import</span> <span class="nn">numpy</span> <span class="k">as</span> <span class="nn">np</span>
2021-11-09 11:51:27 +00:00
<span class="kn">import</span> <span class="nn">pyerrors</span> <span class="k">as</span> <span class="nn">pe</span>
<span class="n">my_obs1</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples1</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble_name&#39;</span><span class="p">])</span>
<span class="n">my_obs2</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples2</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble_name&#39;</span><span class="p">])</span>
<span class="n">my_sum</span> <span class="o">=</span> <span class="n">my_obs1</span> <span class="o">+</span> <span class="n">my_obs2</span>
<span class="n">my_m_eff</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">log</span><span class="p">(</span><span class="n">my_obs1</span> <span class="o">/</span> <span class="n">my_obs2</span><span class="p">)</span>
2021-11-16 11:46:28 +00:00
<span class="n">iamzero</span> <span class="o">=</span> <span class="n">my_m_eff</span> <span class="o">-</span> <span class="n">my_m_eff</span>
<span class="c1"># Check that value and fluctuations are zero within machine precision</span>
<span class="nb">print</span><span class="p">(</span><span class="n">iamzero</span> <span class="o">==</span> <span class="mf">0.0</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="kc">True</span>
2021-11-09 11:51:27 +00:00
</code></pre></div>
<h2 id="error-estimation">Error estimation</h2>
2021-11-15 15:17:23 +00:00
<p>The error estimation within <code><a href="">pyerrors</a></code> is based on the gamma method introduced in <a href="https://arxiv.org/abs/hep-lat/0306017">arXiv:hep-lat/0306017</a>.
After having arrived at the derived quantity of interest the <code>gamma_method</code> can be called as detailed in the following example.</p>
2021-11-15 10:12:16 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_sum</span><span class="o">.</span><span class="n">gamma_method</span><span class="p">()</span>
2021-11-15 15:17:23 +00:00
<span class="nb">print</span><span class="p">(</span><span class="n">my_sum</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mf">1.70</span><span class="p">(</span><span class="mi">57</span><span class="p">)</span>
2021-11-15 10:12:16 +00:00
<span class="n">my_sum</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
2021-11-15 15:17:23 +00:00
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">1.70000000e+00</span> <span class="o">+/-</span> <span class="mf">5.72046658e-01</span> <span class="o">+/-</span> <span class="mf">7.56746598e-02</span> <span class="p">(</span><span class="mf">33.650</span><span class="o">%</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="n">t_int</span> <span class="mf">2.71422900e+00</span> <span class="o">+/-</span> <span class="mf">6.40320983e-01</span> <span class="n">S</span> <span class="o">=</span> <span class="mf">2.00</span>
2021-11-15 10:43:47 +00:00
<span class="o">&gt;</span> <span class="mi">1000</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble_name&#39;</span> <span class="p">:</span> <span class="mi">1000</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">1000</span><span class="p">)</span>
2021-11-15 10:12:16 +00:00
</code></pre></div>
2021-11-15 15:50:10 +00:00
<p>We use the following definition of the integrated autocorrelation time established in <a href="https://link.springer.com/article/10.1007/BF01022990">Madras &amp; Sokal 1988</a>
2021-11-16 11:46:28 +00:00
$$\tau_\mathrm{int}=\frac{1}{2}+\sum_{t=1}^{W}\rho(t)\geq \frac{1}{2}\,.$$
2021-12-11 23:07:13 +00:00
The window $W$ is determined via the automatic windowing procedure described in <a href="https://arxiv.org/abs/hep-lat/0306017">arXiv:hep-lat/0306017</a>.
2021-11-15 15:50:10 +00:00
The standard value for the parameter $S$ of this automatic windowing procedure is $S=2$. Other values for $S$ can be passed to the <code>gamma_method</code> as parameter.</p>
2021-11-15 10:12:16 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_sum</span><span class="o">.</span><span class="n">gamma_method</span><span class="p">(</span><span class="n">S</span><span class="o">=</span><span class="mf">3.0</span><span class="p">)</span>
2021-11-15 10:12:16 +00:00
<span class="n">my_sum</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
2021-11-15 15:17:23 +00:00
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">1.70000000e+00</span> <span class="o">+/-</span> <span class="mf">6.30675201e-01</span> <span class="o">+/-</span> <span class="mf">1.04585650e-01</span> <span class="p">(</span><span class="mf">37.099</span><span class="o">%</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="n">t_int</span> <span class="mf">3.29909703e+00</span> <span class="o">+/-</span> <span class="mf">9.77310102e-01</span> <span class="n">S</span> <span class="o">=</span> <span class="mf">3.00</span>
2021-11-15 10:43:47 +00:00
<span class="o">&gt;</span> <span class="mi">1000</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble_name&#39;</span> <span class="p">:</span> <span class="mi">1000</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">1000</span><span class="p">)</span>
2021-11-15 10:12:16 +00:00
</code></pre></div>
2021-11-15 17:41:36 +00:00
<p>The integrated autocorrelation time $\tau_\mathrm{int}$ and the autocorrelation function $\rho(W)$ can be monitored via the methods <code><a href="pyerrors/obs.html#Obs.plot_tauint">pyerrors.obs.Obs.plot_tauint</a></code> and <code><a href="pyerrors/obs.html#Obs.plot_tauint">pyerrors.obs.Obs.plot_tauint</a></code>.</p>
2021-11-15 10:12:16 +00:00
2022-02-21 14:52:31 +00:00
<p>If the parameter $S$ is set to zero it is assumed that the dataset does not exhibit any autocorrelation and the windowsize is chosen to be zero.
2021-11-17 13:43:40 +00:00
In this case the error estimate is identical to the sample standard error.</p>
2021-11-09 11:51:27 +00:00
<h3 id="exponential-tails">Exponential tails</h3>
2021-11-15 10:57:31 +00:00
<p>Slow modes in the Monte Carlo history can be accounted for by attaching an exponential tail to the autocorrelation function $\rho$ as suggested in <a href="https://arxiv.org/abs/1009.5228">arXiv:1009.5228</a>. The longest autocorrelation time in the history, $\tau_\mathrm{exp}$, can be passed to the <code>gamma_method</code> as parameter. In this case the automatic windowing procedure is vacated and the parameter $S$ does not affect the error estimate.</p>
2021-11-15 10:12:16 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_sum</span><span class="o">.</span><span class="n">gamma_method</span><span class="p">(</span><span class="n">tau_exp</span><span class="o">=</span><span class="mf">7.2</span><span class="p">)</span>
2021-11-15 10:12:16 +00:00
<span class="n">my_sum</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
2021-11-15 15:17:23 +00:00
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">1.70000000e+00</span> <span class="o">+/-</span> <span class="mf">6.28097762e-01</span> <span class="o">+/-</span> <span class="mf">5.79077524e-02</span> <span class="p">(</span><span class="mf">36.947</span><span class="o">%</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="n">t_int</span> <span class="mf">3.27218667e+00</span> <span class="o">+/-</span> <span class="mf">7.99583654e-01</span> <span class="n">tau_exp</span> <span class="o">=</span> <span class="mf">7.20</span><span class="p">,</span> <span class="n">N_sigma</span> <span class="o">=</span> <span class="mi">1</span>
2021-11-15 10:43:47 +00:00
<span class="o">&gt;</span> <span class="mi">1000</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble_name&#39;</span> <span class="p">:</span> <span class="mi">1000</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">1000</span><span class="p">)</span>
2021-11-15 10:12:16 +00:00
</code></pre></div>
2021-12-11 23:07:13 +00:00
<p>For the full API see <code><a href="pyerrors/obs.html#Obs.gamma_method">pyerrors.obs.Obs.gamma_method</a></code>.</p>
2021-11-15 10:12:16 +00:00
2021-11-07 20:53:18 +00:00
<h2 id="multiple-ensemblesreplica">Multiple ensembles/replica</h2>
2021-11-15 10:57:31 +00:00
<p>Error propagation for multiple ensembles (Markov chains with different simulation parameters) is handled automatically. Ensembles are uniquely identified by their <code>name</code>.</p>
2021-11-08 14:53:27 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">obs1</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples1</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">])</span>
2021-11-09 12:09:45 +00:00
<span class="n">obs2</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples2</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble2&#39;</span><span class="p">])</span>
2021-11-08 14:53:27 +00:00
<span class="n">my_sum</span> <span class="o">=</span> <span class="n">obs1</span> <span class="o">+</span> <span class="n">obs2</span>
<span class="n">my_sum</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
2021-11-16 12:00:28 +00:00
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">2.00697958e+00</span>
2021-11-08 14:53:27 +00:00
<span class="o">&gt;</span> <span class="mi">1500</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">2</span> <span class="n">ensembles</span><span class="p">:</span>
2021-11-15 09:55:51 +00:00
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble1&#39;</span> <span class="p">:</span> <span class="mi">1000</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">1000</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble2&#39;</span> <span class="p">:</span> <span class="mi">500</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">500</span><span class="p">)</span>
2021-11-08 14:53:27 +00:00
</code></pre></div>
2021-11-15 14:30:41 +00:00
<p><code><a href="">pyerrors</a></code> identifies multiple replica (independent Markov chains with identical simulation parameters) by the vertical bar <code>|</code> in the name of the data set.</p>
2021-11-08 14:53:27 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">obs1</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples1</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1|r01&#39;</span><span class="p">])</span>
2021-11-09 12:09:45 +00:00
<span class="n">obs2</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples2</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1|r02&#39;</span><span class="p">])</span>
2021-11-08 14:53:27 +00:00
2021-11-15 09:55:51 +00:00
<span class="o">&gt;</span> <span class="n">my_sum</span> <span class="o">=</span> <span class="n">obs1</span> <span class="o">+</span> <span class="n">obs2</span>
<span class="o">&gt;</span> <span class="n">my_sum</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
2021-11-16 12:00:28 +00:00
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">2.00697958e+00</span>
2021-11-08 14:53:27 +00:00
<span class="o">&gt;</span> <span class="mi">1500</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
2021-11-15 09:55:51 +00:00
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble1&#39;</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Replicum</span> <span class="s1">&#39;r01&#39;</span> <span class="p">:</span> <span class="mi">1000</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">1000</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Replicum</span> <span class="s1">&#39;r02&#39;</span> <span class="p">:</span> <span class="mi">500</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">1</span> <span class="n">to</span> <span class="mi">500</span><span class="p">)</span>
2021-11-08 14:53:27 +00:00
</code></pre></div>
2021-11-15 10:12:16 +00:00
<h3 id="error-estimation-for-multiple-ensembles">Error estimation for multiple ensembles</h3>
2021-11-15 10:57:31 +00:00
<p>In order to keep track of different error analysis parameters for different ensembles one can make use of global dictionaries as detailed in the following example.</p>
2021-11-15 10:12:16 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="o">.</span><span class="n">S_dict</span><span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">]</span> <span class="o">=</span> <span class="mf">2.5</span>
2021-11-15 10:12:16 +00:00
<span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="o">.</span><span class="n">tau_exp_dict</span><span class="p">[</span><span class="s1">&#39;ensemble2&#39;</span><span class="p">]</span> <span class="o">=</span> <span class="mf">8.0</span>
<span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="o">.</span><span class="n">tau_exp_dict</span><span class="p">[</span><span class="s1">&#39;ensemble3&#39;</span><span class="p">]</span> <span class="o">=</span> <span class="mf">2.0</span>
</code></pre></div>
2021-11-15 10:43:47 +00:00
<p>In case the <code>gamma_method</code> is called without any parameters it will use the values specified in the dictionaries for the respective ensembles.
Passing arguments to the <code>gamma_method</code> still dominates over the dictionaries.</p>
2021-11-07 20:53:18 +00:00
<h2 id="irregular-monte-carlo-chains">Irregular Monte Carlo chains</h2>
2022-01-21 14:55:18 +00:00
<p><code>Obs</code> objects defined on irregular Monte Carlo chains can be initialized with the parameter <code>idl</code>.</p>
2021-11-08 14:53:27 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="c1"># Observable defined on configurations 20 to 519</span>
2021-11-08 15:04:55 +00:00
<span class="n">obs1</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples1</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">],</span> <span class="n">idl</span><span class="o">=</span><span class="p">[</span><span class="nb">range</span><span class="p">(</span><span class="mi">20</span><span class="p">,</span> <span class="mi">520</span><span class="p">)])</span>
2021-11-16 12:00:28 +00:00
<span class="n">obs1</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">9.98319881e-01</span>
<span class="o">&gt;</span> <span class="mi">500</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble1&#39;</span> <span class="p">:</span> <span class="mi">500</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">20</span> <span class="n">to</span> <span class="mi">519</span><span class="p">)</span>
2021-11-08 15:10:26 +00:00
2021-11-08 15:04:55 +00:00
<span class="c1"># Observable defined on every second configuration between 5 and 1003</span>
<span class="n">obs2</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples2</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">],</span> <span class="n">idl</span><span class="o">=</span><span class="p">[</span><span class="nb">range</span><span class="p">(</span><span class="mi">5</span><span class="p">,</span> <span class="mi">1005</span><span class="p">,</span> <span class="mi">2</span><span class="p">)])</span>
2021-11-16 12:00:28 +00:00
<span class="n">obs2</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">9.99100712e-01</span>
<span class="o">&gt;</span> <span class="mi">500</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble1&#39;</span> <span class="p">:</span> <span class="mi">500</span> <span class="n">configurations</span> <span class="p">(</span><span class="kn">from</span> <span class="mi">5</span> <span class="n">to</span> <span class="mi">1003</span> <span class="ow">in</span> <span class="n">steps</span> <span class="n">of</span> <span class="mi">2</span><span class="p">)</span>
2021-11-08 15:10:26 +00:00
2021-11-08 15:04:55 +00:00
<span class="c1"># Observable defined on configurations 2, 9, 28, 29 and 501</span>
<span class="n">obs3</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples3</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">],</span> <span class="n">idl</span><span class="o">=</span><span class="p">[[</span><span class="mi">2</span><span class="p">,</span> <span class="mi">9</span><span class="p">,</span> <span class="mi">28</span><span class="p">,</span> <span class="mi">29</span><span class="p">,</span> <span class="mi">501</span><span class="p">]])</span>
2021-11-16 12:00:28 +00:00
<span class="n">obs3</span><span class="o">.</span><span class="n">details</span><span class="p">()</span>
<span class="o">&gt;</span> <span class="n">Result</span> <span class="mf">1.01718064e+00</span>
<span class="o">&gt;</span> <span class="mi">5</span> <span class="n">samples</span> <span class="ow">in</span> <span class="mi">1</span> <span class="n">ensemble</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="err">·</span> <span class="n">Ensemble</span> <span class="s1">&#39;ensemble1&#39;</span> <span class="p">:</span> <span class="mi">5</span> <span class="n">configurations</span> <span class="p">(</span><span class="n">irregular</span> <span class="nb">range</span><span class="p">)</span>
2021-11-08 14:53:27 +00:00
</code></pre></div>
2022-01-21 14:55:18 +00:00
<p><code>Obs</code> objects defined on regular and irregular histories of the same ensemble can be computed with each other and the correct error propagation and estimation is automatically taken care of.</p>
2021-11-08 15:04:55 +00:00
<p><strong>Warning:</strong> Irregular Monte Carlo chains can result in odd patterns in the autocorrelation functions.
2021-11-15 10:57:31 +00:00
Make sure to check the autocorrelation time with e.g. <code><a href="pyerrors/obs.html#Obs.plot_rho">pyerrors.obs.Obs.plot_rho</a></code> or <code><a href="pyerrors/obs.html#Obs.plot_tauint">pyerrors.obs.Obs.plot_tauint</a></code>.</p>
2021-11-08 15:04:55 +00:00
2021-12-11 23:07:13 +00:00
<p>For the full API see <code><a href="pyerrors/obs.html#Obs">pyerrors.obs.Obs</a></code>.</p>
2021-11-07 20:53:18 +00:00
2021-11-09 11:51:27 +00:00
<h1 id="correlators">Correlators</h1>
2021-11-07 20:53:18 +00:00
2022-02-15 12:25:34 +00:00
<p>When one is not interested in single observables but correlation functions, <code><a href="">pyerrors</a></code> offers the <code>Corr</code> class which simplifies the corresponding error propagation and provides the user with a set of standard methods. In order to initialize a <code>Corr</code> objects one needs to arrange the data as a list of <code>Obs</code></p>
2021-12-23 15:43:56 +00:00
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_corr</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Corr</span><span class="p">([</span><span class="n">obs_0</span><span class="p">,</span> <span class="n">obs_1</span><span class="p">,</span> <span class="n">obs_2</span><span class="p">,</span> <span class="n">obs_3</span><span class="p">])</span>
2021-12-23 15:43:56 +00:00
<span class="nb">print</span><span class="p">(</span><span class="n">my_corr</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="n">x0</span><span class="o">/</span><span class="n">a</span> <span class="n">Corr</span><span class="p">(</span><span class="n">x0</span><span class="o">/</span><span class="n">a</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="o">------------------</span>
<span class="o">&gt;</span> <span class="mi">0</span> <span class="mf">0.7957</span><span class="p">(</span><span class="mi">80</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">1</span> <span class="mf">0.5156</span><span class="p">(</span><span class="mi">51</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">2</span> <span class="mf">0.3227</span><span class="p">(</span><span class="mi">33</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">3</span> <span class="mf">0.2041</span><span class="p">(</span><span class="mi">21</span><span class="p">)</span>
</code></pre></div>
<p>In case the correlation functions are not defined on the outermost timeslices, for example because of fixed boundary conditions, a padding can be introduced.</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_corr</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Corr</span><span class="p">([</span><span class="n">obs_0</span><span class="p">,</span> <span class="n">obs_1</span><span class="p">,</span> <span class="n">obs_2</span><span class="p">,</span> <span class="n">obs_3</span><span class="p">],</span> <span class="n">padding</span><span class="o">=</span><span class="p">[</span><span class="mi">1</span><span class="p">,</span> <span class="mi">1</span><span class="p">])</span>
2021-12-23 15:43:56 +00:00
<span class="nb">print</span><span class="p">(</span><span class="n">my_corr</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="n">x0</span><span class="o">/</span><span class="n">a</span> <span class="n">Corr</span><span class="p">(</span><span class="n">x0</span><span class="o">/</span><span class="n">a</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="o">------------------</span>
<span class="o">&gt;</span> <span class="mi">0</span>
<span class="o">&gt;</span> <span class="mi">1</span> <span class="mf">0.7957</span><span class="p">(</span><span class="mi">80</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">2</span> <span class="mf">0.5156</span><span class="p">(</span><span class="mi">51</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">3</span> <span class="mf">0.3227</span><span class="p">(</span><span class="mi">33</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">4</span> <span class="mf">0.2041</span><span class="p">(</span><span class="mi">21</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">5</span>
</code></pre></div>
<p>The individual entries of a correlator can be accessed via slicing</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="nb">print</span><span class="p">(</span><span class="n">my_corr</span><span class="p">[</span><span class="mi">3</span><span class="p">])</span>
2021-12-23 15:43:56 +00:00
<span class="o">&gt;</span> <span class="mf">0.3227</span><span class="p">(</span><span class="mi">33</span><span class="p">)</span>
</code></pre></div>
<p>Error propagation with the <code>Corr</code> class works very similar to <code>Obs</code> objects. Mathematical operations are overloaded and <code>Corr</code> objects can be computed together with other <code>Corr</code> objects, <code>Obs</code> objects or real numbers and integers.</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_new_corr</span> <span class="o">=</span> <span class="mf">0.3</span> <span class="o">*</span> <span class="n">my_corr</span><span class="p">[</span><span class="mi">2</span><span class="p">]</span> <span class="o">*</span> <span class="n">my_corr</span> <span class="o">*</span> <span class="n">my_corr</span> <span class="o">+</span> <span class="mi">12</span> <span class="o">/</span> <span class="n">my_corr</span>
2021-12-23 15:43:56 +00:00
</code></pre></div>
<p><code><a href="">pyerrors</a></code> provides the user with a set of regularly used methods for the manipulation of correlator objects:</p>
<ul>
<li><code>Corr.gamma_method</code> applies the gamma method to all entries of the correlator.</li>
<li><code>Corr.m_eff</code> to construct effective masses. Various variants for periodic and fixed temporal boundary conditions are available.</li>
<li><code>Corr.deriv</code> returns the first derivative of the correlator as <code>Corr</code>. Different discretizations of the numerical derivative are available.</li>
<li><code>Corr.second_deriv</code> returns the second derivative of the correlator as <code>Corr</code>. Different discretizations of the numerical derivative are available.</li>
<li><code>Corr.symmetric</code> symmetrizes parity even correlations functions, assuming periodic boundary conditions.</li>
<li><code>Corr.anti_symmetric</code> anti-symmetrizes parity odd correlations functions, assuming periodic boundary conditions.</li>
<li><code>Corr.T_symmetry</code> averages a correlator with its time symmetry partner, assuming fixed boundary conditions.</li>
<li><code>Corr.plateau</code> extracts a plateau value from the correlator in a given range.</li>
<li><code>Corr.roll</code> periodically shifts the correlator.</li>
<li><code>Corr.reverse</code> reverses the time ordering of the correlator.</li>
<li><code>Corr.correlate</code> constructs a disconnected correlation function from the correlator and another <code>Corr</code> or <code>Obs</code> object.</li>
<li><code>Corr.reweight</code> reweights the correlator.</li>
</ul>
2022-01-21 14:55:18 +00:00
<p><code><a href="">pyerrors</a></code> can also handle matrices of correlation functions and extract energy states from these matrices via a generalized eigenvalue problem (see <code><a href="pyerrors/correlators.html#Corr.GEVP">pyerrors.correlators.Corr.GEVP</a></code>).</p>
2021-12-23 15:43:56 +00:00
2021-12-11 23:07:13 +00:00
<p>For the full API see <code><a href="pyerrors/correlators.html#Corr">pyerrors.correlators.Corr</a></code>.</p>
2021-11-07 20:53:18 +00:00
2022-02-14 14:08:21 +00:00
<h1 id="complex-valued-observables">Complex valued observables</h1>
2021-11-07 21:09:48 +00:00
2021-11-15 14:59:01 +00:00
<p><code><a href="">pyerrors</a></code> can handle complex valued observables via the class <code><a href="pyerrors/obs.html#CObs">pyerrors.obs.CObs</a></code>.
<code>CObs</code> are initialized with a real and an imaginary part which both can be <code>Obs</code> valued.</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_real_part</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples1</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">])</span>
2021-11-15 14:59:01 +00:00
<span class="n">my_imag_part</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Obs</span><span class="p">([</span><span class="n">samples2</span><span class="p">],</span> <span class="p">[</span><span class="s1">&#39;ensemble1&#39;</span><span class="p">])</span>
<span class="n">my_cobs</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">CObs</span><span class="p">(</span><span class="n">my_real_part</span><span class="p">,</span> <span class="n">my_imag_part</span><span class="p">)</span>
<span class="n">my_cobs</span><span class="o">.</span><span class="n">gamma_method</span><span class="p">()</span>
<span class="nb">print</span><span class="p">(</span><span class="n">my_cobs</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="p">(</span><span class="mf">0.9959</span><span class="p">(</span><span class="mi">91</span><span class="p">)</span><span class="o">+</span><span class="mf">0.659</span><span class="p">(</span><span class="mi">28</span><span class="p">)</span><span class="n">j</span><span class="p">)</span>
</code></pre></div>
<p>Elementary mathematical operations are overloaded and samples are properly propagated as for the <code>Obs</code> class.</p>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_derived_cobs</span> <span class="o">=</span> <span class="p">(</span><span class="n">my_cobs</span> <span class="o">+</span> <span class="n">my_cobs</span><span class="o">.</span><span class="n">conjugate</span><span class="p">())</span> <span class="o">/</span> <span class="n">np</span><span class="o">.</span><span class="n">abs</span><span class="p">(</span><span class="n">my_cobs</span><span class="p">)</span>
2021-11-15 14:59:01 +00:00
<span class="n">my_derived_cobs</span><span class="o">.</span><span class="n">gamma_method</span><span class="p">()</span>
<span class="nb">print</span><span class="p">(</span><span class="n">my_derived_cobs</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="p">(</span><span class="mf">1.668</span><span class="p">(</span><span class="mi">23</span><span class="p">)</span><span class="o">+</span><span class="mf">0.0</span><span class="n">j</span><span class="p">)</span>
</code></pre></div>
2021-11-07 21:09:48 +00:00
2022-02-14 14:08:21 +00:00
<h1 id="error-propagation-in-iterative-algorithms">Error propagation in iterative algorithms</h1>
<p><code><a href="">pyerrors</a></code> supports exact linear error propagation for iterative algorithms like various variants of non-linear least sqaures fits or root finding. The derivatives required for the error propagation are calculated as described in <a href="https://arxiv.org/abs/1809.01289">arXiv:1809.01289</a>.</p>
<h2 id="least-squares-fits">Least squares fits</h2>
<p>Standard non-linear least square fits with errors on the dependent but not the independent variables can be performed with <code><a href="pyerrors/fits.html#least_squares">pyerrors.fits.least_squares</a></code>. As default solver the Levenberg-Marquardt algorithm implemented in <a href="https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.least_squares.html">scipy</a> is used.</p>
<p>Fit functions have to be of the following form</p>
<div class="pdoc-code codehilite"><pre><span></span><code><span class="kn">import</span> <span class="nn">autograd.numpy</span> <span class="k">as</span> <span class="nn">anp</span>
<span class="k">def</span> <span class="nf">func</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">x</span><span class="p">):</span>
<span class="k">return</span> <span class="n">a</span><span class="p">[</span><span class="mi">1</span><span class="p">]</span> <span class="o">*</span> <span class="n">anp</span><span class="o">.</span><span class="n">exp</span><span class="p">(</span><span class="o">-</span><span class="n">a</span><span class="p">[</span><span class="mi">0</span><span class="p">]</span> <span class="o">*</span> <span class="n">x</span><span class="p">)</span>
</code></pre></div>
2022-02-21 14:52:31 +00:00
<p><strong>It is important that numerical functions refer to <code>autograd.numpy</code> instead of <code>numpy</code> for the automatic differentiation in iterative algorithms to work properly.</strong></p>
2022-02-14 14:08:21 +00:00
<p>Fits can then be performed via</p>
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">fit_result</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">fits</span><span class="o">.</span><span class="n">least_squares</span><span class="p">(</span><span class="n">x</span><span class="p">,</span> <span class="n">y</span><span class="p">,</span> <span class="n">func</span><span class="p">)</span>
<span class="nb">print</span><span class="p">(</span><span class="s2">&quot;</span><span class="se">\n</span><span class="s2">&quot;</span><span class="p">,</span> <span class="n">fit_result</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="n">Fit</span> <span class="k">with</span> <span class="mi">2</span> <span class="n">parameters</span>
<span class="o">&gt;</span> <span class="n">Method</span><span class="p">:</span> <span class="n">Levenberg</span><span class="o">-</span><span class="n">Marquardt</span>
<span class="o">&gt;</span> <span class="err">`</span><span class="n">ftol</span><span class="err">`</span> <span class="n">termination</span> <span class="n">condition</span> <span class="ow">is</span> <span class="n">satisfied</span><span class="o">.</span>
<span class="o">&gt;</span> <span class="n">chisquare</span><span class="o">/</span><span class="n">d</span><span class="o">.</span><span class="n">o</span><span class="o">.</span><span class="n">f</span><span class="o">.</span><span class="p">:</span> <span class="mf">0.9593035785160936</span>
<span class="o">&gt;</span> <span class="n">Goodness</span> <span class="n">of</span> <span class="n">fit</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="n">χ</span><span class="err">²</span><span class="o">/</span><span class="n">d</span><span class="o">.</span><span class="n">o</span><span class="o">.</span><span class="n">f</span><span class="o">.</span> <span class="o">=</span> <span class="mf">0.959304</span>
<span class="o">&gt;</span> <span class="n">p</span><span class="o">-</span><span class="n">value</span> <span class="o">=</span> <span class="mf">0.5673</span>
<span class="o">&gt;</span> <span class="n">Fit</span> <span class="n">parameters</span><span class="p">:</span>
<span class="o">&gt;</span> <span class="mi">0</span> <span class="mf">0.0548</span><span class="p">(</span><span class="mi">28</span><span class="p">)</span>
<span class="o">&gt;</span> <span class="mi">1</span> <span class="mf">1.933</span><span class="p">(</span><span class="mi">64</span><span class="p">)</span>
</code></pre></div>
<p>where x is a <code>list</code> or <code>numpy.array</code> of <code>floats</code> and y is a <code>list</code> or <code>numpy.array</code> of <code>Obs</code>.</p>
<p>Data stored in <code>Corr</code> objects can be fitted directly using the <code>Corr.fit</code> method.</p>
<div class="pdoc-code codehilite"><pre><span></span><code><span class="n">my_corr</span> <span class="o">=</span> <span class="n">pe</span><span class="o">.</span><span class="n">Corr</span><span class="p">(</span><span class="n">y</span><span class="p">)</span>
<span class="n">fit_result</span> <span class="o">=</span> <span class="n">my_corr</span><span class="o">.</span><span class="n">fit</span><span class="p">(</span><span class="n">func</span><span class="p">,</span> <span class="n">fitrange</span><span class="o">=</span><span class="p">[</span><span class="mi">12</span><span class="p">,</span> <span class="mi">25</span><span class="p">])</span>
</code></pre></div>
<p>this can simplify working with absolute fit ranges and takes care of gaps in the data automatically.</p>
<p>For fit functions with multiple independent variables the fit function can be of the form</p>
<div class="pdoc-code codehilite"><pre><span></span><code><span class="k">def</span> <span class="nf">func</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">x</span><span class="p">):</span>
<span class="p">(</span><span class="n">x1</span><span class="p">,</span> <span class="n">x2</span><span class="p">)</span> <span class="o">=</span> <span class="n">x</span>
<span class="k">return</span> <span class="n">a</span><span class="p">[</span><span class="mi">0</span><span class="p">]</span> <span class="o">*</span> <span class="n">x1</span> <span class="o">**</span> <span class="mi">2</span> <span class="o">+</span> <span class="n">a</span><span class="p">[</span><span class="mi">1</span><span class="p">]</span> <span class="o">*</span> <span class="n">x2</span>
</code></pre></div>
2021-11-07 20:53:18 +00:00
2022-02-14 14:08:21 +00:00
<h2 id="total-least-squares-fits">Total least squares fits</h2>
<p><code><a href="">pyerrors</a></code> can also fit data with errors on both the dependent and independent variables using the total least squares method also referred to orthogonal distance regression as implemented in <a href="https://docs.scipy.org/doc/scipy/reference/odr.html">scipy</a>, see <code><a href="pyerrors/fits.html#least_squares">pyerrors.fits.least_squares</a></code>. The syntax is identical to the standard least squares case, the only diffrence being that <code>x</code> also has to be a <code>list</code> or <code>numpy.array</code> of <code>Obs</code>.</p>
<p>For the full API see <code><a href="pyerrors/fits.html">pyerrors.fits</a></code> for fits and <code><a href="pyerrors/roots.html">pyerrors.roots</a></code> for finding roots of functions.</p>
2021-11-07 21:09:48 +00:00
2021-11-07 20:53:18 +00:00
<h1 id="matrix-operations">Matrix operations</h1>
2022-02-14 14:26:48 +00:00
<p><code><a href="">pyerrors</a></code> provides wrappers for <code>Obs</code>- and <code>CObs</code>-valued matrix operations based on <code>numpy.linalg</code>. The supported functions include:</p>
2022-01-21 14:55:18 +00:00
<ul>
<li><code>inv</code> for the matrix inverse.</li>
<li><code>cholseky</code> for the Cholesky decomposition.</li>
<li><code>det</code> for the matrix determinant.</li>
<li><code>eigh</code> for eigenvalues and eigenvectors of hermitean matrices.</li>
<li><code>eig</code> for eigenvalues of general matrices.</li>
<li><code>pinv</code> for the Moore-Penrose pseudoinverse.</li>
<li><code>svd</code> for the singular-value-decomposition.</li>
</ul>
<p>For the full API see <code><a href="pyerrors/linalg.html">pyerrors.linalg</a></code>.</p>
2021-11-07 21:09:48 +00:00
2021-11-15 14:30:41 +00:00
<h1 id="export-data">Export data</h1>
2022-02-21 14:52:31 +00:00
<p>The preferred exported file format within <code><a href="">pyerrors</a></code> is json.gz. The exact specifications of this format will be listed here soon.</p>
2021-11-15 14:30:41 +00:00
<h2 id="jackknife-samples">Jackknife samples</h2>
2021-12-11 23:07:13 +00:00
<p>For comparison with other analysis workflows <code><a href="">pyerrors</a></code> can generate jackknife samples from an <code>Obs</code> object or import jackknife samples into an <code>Obs</code> object.
See <code><a href="pyerrors/obs.html#Obs.export_jackknife">pyerrors.obs.Obs.export_jackknife</a></code> and <code><a href="pyerrors/obs.html#import_jackknife">pyerrors.obs.import_jackknife</a></code> for details.</p>
2021-11-15 14:30:41 +00:00
2022-02-14 14:08:21 +00:00
<h1 id="citing">Citing</h1>
<p>If you use <code><a href="">pyerrors</a></code> for research that leads to a publication please consider citing:</p>
<ul>
2022-02-14 14:19:54 +00:00
<li>Ulli Wolff, <em>Monte Carlo errors with less errors</em>. Comput.Phys.Commun. 156 (2004) 143-153, Comput.Phys.Commun. 176 (2007) 383 (erratum).</li>
<li>Stefan Schaefer, Rainer Sommer, Francesco Virotta, <em>Critical slowing down and error analysis in lattice QCD simulations</em>. Nucl.Phys.B 845 (2011) 93-119.</li>
<li>Alberto Ramos, <em>Automatic differentiation for error analysis of Monte Carlo data</em>. Comput.Phys.Commun. 238 (2019) 19-35.</li>
2022-02-14 14:08:21 +00:00
</ul>
2021-11-07 20:53:18 +00:00
</div>
<details>
<summary>View Source</summary>
2022-02-14 14:08:21 +00:00
<div class="pdoc-code codehilite"><pre><span></span><span class="sa">r</span><span class="sd">&#39;&#39;&#39;</span>
2021-11-07 20:53:18 +00:00
<span class="sd"># What is pyerrors?</span>
<span class="sd">`pyerrors` is a python package for error computation and propagation of Markov chain Monte Carlo data.</span>
2021-12-11 23:00:30 +00:00
<span class="sd">It is based on the gamma method [arXiv:hep-lat/0306017](https://arxiv.org/abs/hep-lat/0306017). Some of its features are:</span>
<span class="sd">- automatic differentiation for exact liner error propagation as suggested in [arXiv:1809.01289](https://arxiv.org/abs/1809.01289) (partly based on the [autograd](https://github.com/HIPS/autograd) package).</span>
<span class="sd">- treatment of slow modes in the simulation as suggested in [arXiv:1009.5228](https://arxiv.org/abs/1009.5228).</span>
<span class="sd">- coherent error propagation for data from different Markov chains.</span>
<span class="sd">- non-linear fits with x- and y-errors and exact linear error propagation based on automatic differentiation as introduced in [arXiv:1809.01289](https://arxiv.org/abs/1809.01289).</span>
<span class="sd">- real and complex matrix operations and their error propagation based on automatic differentiation (Matrix inverse, Cholesky decomposition, calculation of eigenvalues and eigenvectors, singular value decomposition...).</span>
<span class="sd">There exist similar publicly available implementations of gamma method error analysis suites in [Fortran](https://gitlab.ift.uam-csic.es/alberto/aderrors), [Julia](https://gitlab.ift.uam-csic.es/alberto/aderrors.jl) and [Python](https://github.com/mbruno46/pyobs).</span>
2021-11-30 14:58:46 +00:00
2021-11-16 11:46:28 +00:00
<span class="sd">## Basic example</span>
2021-11-07 20:53:18 +00:00
<span class="sd">```python</span>
<span class="sd">import numpy as np</span>
<span class="sd">import pyerrors as pe</span>
2021-11-16 11:50:48 +00:00
<span class="sd">my_obs = pe.Obs([samples], [&#39;ensemble_name&#39;]) # Initialize an Obs object</span>
2021-11-16 11:46:28 +00:00
<span class="sd">my_new_obs = 2 * np.log(my_obs) / my_obs ** 2 # Construct derived Obs object</span>
2021-11-16 11:50:48 +00:00
<span class="sd">my_new_obs.gamma_method() # Estimate the statistical error</span>
<span class="sd">print(my_new_obs) # Print the result to stdout</span>
2021-11-15 13:13:15 +00:00
<span class="sd">&gt; 0.31498(72)</span>
2021-11-07 20:53:18 +00:00
<span class="sd">```</span>
2021-11-15 13:13:15 +00:00
2021-11-07 20:53:18 +00:00
<span class="sd"># The `Obs` class</span>
2021-11-09 11:51:27 +00:00
<span class="sd">`pyerrors` introduces a new datatype, `Obs`, which simplifies error propagation and estimation for auto- and cross-correlated data.</span>
2022-02-21 14:52:31 +00:00
<span class="sd">An `Obs` object can be initialized with two arguments, the first is a list containing the samples for an observable from a Monte Carlo chain.</span>
2021-11-09 11:51:27 +00:00
<span class="sd">The samples can either be provided as python list or as numpy array.</span>
<span class="sd">The second argument is a list containing the names of the respective Monte Carlo chains as strings. These strings uniquely identify a Monte Carlo chain/ensemble.</span>
2021-11-07 20:53:18 +00:00
<span class="sd">```python</span>
<span class="sd">import pyerrors as pe</span>
<span class="sd">my_obs = pe.Obs([samples], [&#39;ensemble_name&#39;])</span>
<span class="sd">```</span>
2021-11-09 11:51:27 +00:00
<span class="sd">## Error propagation</span>
2022-02-21 14:52:31 +00:00
<span class="sd">When performing mathematical operations on `Obs` objects the correct error propagation is intrinsically taken care of using a first order Taylor expansion</span>
2021-11-16 11:46:28 +00:00
<span class="sd">$$\delta_f^i=\sum_\alpha \bar{f}_\alpha \delta_\alpha^i\,,\quad \delta_\alpha^i=a_\alpha^i-\bar{a}_\alpha\,,$$</span>
2021-11-09 11:51:27 +00:00
<span class="sd">as introduced in [arXiv:hep-lat/0306017](https://arxiv.org/abs/hep-lat/0306017).</span>
<span class="sd">The required derivatives $\bar{f}_\alpha$ are evaluated up to machine precision via automatic differentiation as suggested in [arXiv:1809.01289](https://arxiv.org/abs/1809.01289).</span>
<span class="sd">The `Obs` class is designed such that mathematical numpy functions can be used on `Obs` just as for regular floats.</span>
<span class="sd">```python</span>
<span class="sd">import numpy as np</span>
<span class="sd">import pyerrors as pe</span>
<span class="sd">my_obs1 = pe.Obs([samples1], [&#39;ensemble_name&#39;])</span>
<span class="sd">my_obs2 = pe.Obs([samples2], [&#39;ensemble_name&#39;])</span>
<span class="sd">my_sum = my_obs1 + my_obs2</span>
<span class="sd">my_m_eff = np.log(my_obs1 / my_obs2)</span>
2021-11-16 11:46:28 +00:00
<span class="sd">iamzero = my_m_eff - my_m_eff</span>
<span class="sd"># Check that value and fluctuations are zero within machine precision</span>
<span class="sd">print(iamzero == 0.0)</span>
<span class="sd">&gt; True</span>
2021-11-09 11:51:27 +00:00
<span class="sd">```</span>
<span class="sd">## Error estimation</span>
2021-11-15 10:57:31 +00:00
<span class="sd">The error estimation within `pyerrors` is based on the gamma method introduced in [arXiv:hep-lat/0306017](https://arxiv.org/abs/hep-lat/0306017).</span>
<span class="sd">After having arrived at the derived quantity of interest the `gamma_method` can be called as detailed in the following example.</span>
2021-11-15 10:12:16 +00:00
<span class="sd">```python</span>
<span class="sd">my_sum.gamma_method()</span>
2021-11-15 15:17:23 +00:00
<span class="sd">print(my_sum)</span>
<span class="sd">&gt; 1.70(57)</span>
2021-11-15 10:12:16 +00:00
<span class="sd">my_sum.details()</span>
2021-11-15 15:17:23 +00:00
<span class="sd">&gt; Result 1.70000000e+00 +/- 5.72046658e-01 +/- 7.56746598e-02 (33.650%)</span>
<span class="sd">&gt; t_int 2.71422900e+00 +/- 6.40320983e-01 S = 2.00</span>
2021-11-15 10:43:47 +00:00
<span class="sd">&gt; 1000 samples in 1 ensemble:</span>
<span class="sd">&gt; · Ensemble &#39;ensemble_name&#39; : 1000 configurations (from 1 to 1000)</span>
2021-11-15 15:17:23 +00:00
2021-11-15 10:12:16 +00:00
<span class="sd">```</span>
2021-11-15 15:50:10 +00:00
<span class="sd">We use the following definition of the integrated autocorrelation time established in [Madras &amp; Sokal 1988](https://link.springer.com/article/10.1007/BF01022990)</span>
2021-11-16 11:46:28 +00:00
<span class="sd">$$\tau_\mathrm{int}=\frac{1}{2}+\sum_{t=1}^{W}\rho(t)\geq \frac{1}{2}\,.$$</span>
2021-12-11 23:07:13 +00:00
<span class="sd">The window $W$ is determined via the automatic windowing procedure described in [arXiv:hep-lat/0306017](https://arxiv.org/abs/hep-lat/0306017).</span>
2021-11-15 15:50:10 +00:00
<span class="sd">The standard value for the parameter $S$ of this automatic windowing procedure is $S=2$. Other values for $S$ can be passed to the `gamma_method` as parameter.</span>
2021-11-15 10:12:16 +00:00
<span class="sd">```python</span>
<span class="sd">my_sum.gamma_method(S=3.0)</span>
<span class="sd">my_sum.details()</span>
2021-11-15 15:17:23 +00:00
<span class="sd">&gt; Result 1.70000000e+00 +/- 6.30675201e-01 +/- 1.04585650e-01 (37.099%)</span>
<span class="sd">&gt; t_int 3.29909703e+00 +/- 9.77310102e-01 S = 3.00</span>
2021-11-15 10:43:47 +00:00
<span class="sd">&gt; 1000 samples in 1 ensemble:</span>
<span class="sd">&gt; · Ensemble &#39;ensemble_name&#39; : 1000 configurations (from 1 to 1000)</span>
2021-11-15 10:12:16 +00:00
<span class="sd">```</span>
2021-11-15 17:41:36 +00:00
<span class="sd">The integrated autocorrelation time $\tau_\mathrm{int}$ and the autocorrelation function $\rho(W)$ can be monitored via the methods `pyerrors.obs.Obs.plot_tauint` and `pyerrors.obs.Obs.plot_tauint`.</span>
2021-11-15 10:12:16 +00:00
2022-02-21 14:52:31 +00:00
<span class="sd">If the parameter $S$ is set to zero it is assumed that the dataset does not exhibit any autocorrelation and the windowsize is chosen to be zero.</span>
2021-11-17 13:43:40 +00:00
<span class="sd">In this case the error estimate is identical to the sample standard error.</span>
2021-11-09 11:51:27 +00:00
<span class="sd">### Exponential tails</span>
2021-11-15 10:57:31 +00:00
<span class="sd">Slow modes in the Monte Carlo history can be accounted for by attaching an exponential tail to the autocorrelation function $\rho$ as suggested in [arXiv:1009.5228](https://arxiv.org/abs/1009.5228). The longest autocorrelation time in the history, $\tau_\mathrm{exp}$, can be passed to the `gamma_method` as parameter. In this case the automatic windowing procedure is vacated and the parameter $S$ does not affect the error estimate.</span>
2021-11-15 10:12:16 +00:00
<span class="sd">```python</span>
2021-11-15 10:43:47 +00:00
<span class="sd">my_sum.gamma_method(tau_exp=7.2)</span>
2021-11-15 10:12:16 +00:00
<span class="sd">my_sum.details()</span>
2021-11-15 15:17:23 +00:00
<span class="sd">&gt; Result 1.70000000e+00 +/- 6.28097762e-01 +/- 5.79077524e-02 (36.947%)</span>
<span class="sd">&gt; t_int 3.27218667e+00 +/- 7.99583654e-01 tau_exp = 7.20, N_sigma = 1</span>
2021-11-15 10:43:47 +00:00
<span class="sd">&gt; 1000 samples in 1 ensemble:</span>
<span class="sd">&gt; · Ensemble &#39;ensemble_name&#39; : 1000 configurations (from 1 to 1000)</span>
2021-11-15 10:12:16 +00:00
<span class="sd">```</span>
2021-12-11 23:07:13 +00:00
<span class="sd">For the full API see `pyerrors.obs.Obs.gamma_method`.</span>
2021-11-15 10:12:16 +00:00
2021-11-07 20:53:18 +00:00
<span class="sd">## Multiple ensembles/replica</span>
2021-11-15 10:57:31 +00:00
<span class="sd">Error propagation for multiple ensembles (Markov chains with different simulation parameters) is handled automatically. Ensembles are uniquely identified by their `name`.</span>
2021-11-08 14:53:27 +00:00
<span class="sd">```python</span>
<span class="sd">obs1 = pe.Obs([samples1], [&#39;ensemble1&#39;])</span>
2021-11-09 12:09:45 +00:00
<span class="sd">obs2 = pe.Obs([samples2], [&#39;ensemble2&#39;])</span>
2021-11-08 14:53:27 +00:00
<span class="sd">my_sum = obs1 + obs2</span>
<span class="sd">my_sum.details()</span>
2021-11-16 12:00:28 +00:00
<span class="sd">&gt; Result 2.00697958e+00</span>
2021-11-08 14:53:27 +00:00
<span class="sd">&gt; 1500 samples in 2 ensembles:</span>
2021-11-15 09:55:51 +00:00
<span class="sd">&gt; · Ensemble &#39;ensemble1&#39; : 1000 configurations (from 1 to 1000)</span>
<span class="sd">&gt; · Ensemble &#39;ensemble2&#39; : 500 configurations (from 1 to 500)</span>
2021-11-08 14:53:27 +00:00
<span class="sd">```</span>
2021-11-15 14:30:41 +00:00
<span class="sd">`pyerrors` identifies multiple replica (independent Markov chains with identical simulation parameters) by the vertical bar `|` in the name of the data set.</span>
2021-11-08 14:53:27 +00:00
<span class="sd">```python</span>
<span class="sd">obs1 = pe.Obs([samples1], [&#39;ensemble1|r01&#39;])</span>
2021-11-09 12:09:45 +00:00
<span class="sd">obs2 = pe.Obs([samples2], [&#39;ensemble1|r02&#39;])</span>
2021-11-08 14:53:27 +00:00
2021-11-15 09:55:51 +00:00
<span class="sd">&gt; my_sum = obs1 + obs2</span>
<span class="sd">&gt; my_sum.details()</span>
2021-11-16 12:00:28 +00:00
<span class="sd">&gt; Result 2.00697958e+00</span>
2021-11-08 14:53:27 +00:00
<span class="sd">&gt; 1500 samples in 1 ensemble:</span>
2021-11-15 09:55:51 +00:00
<span class="sd">&gt; · Ensemble &#39;ensemble1&#39;</span>
<span class="sd">&gt; · Replicum &#39;r01&#39; : 1000 configurations (from 1 to 1000)</span>
<span class="sd">&gt; · Replicum &#39;r02&#39; : 500 configurations (from 1 to 500)</span>
2021-11-08 14:53:27 +00:00
<span class="sd">```</span>
2021-11-15 10:12:16 +00:00
<span class="sd">### Error estimation for multiple ensembles</span>
2021-11-15 10:57:31 +00:00
<span class="sd">In order to keep track of different error analysis parameters for different ensembles one can make use of global dictionaries as detailed in the following example.</span>
2021-11-15 10:12:16 +00:00
<span class="sd">```python</span>
<span class="sd">pe.Obs.S_dict[&#39;ensemble1&#39;] = 2.5</span>
<span class="sd">pe.Obs.tau_exp_dict[&#39;ensemble2&#39;] = 8.0</span>
<span class="sd">pe.Obs.tau_exp_dict[&#39;ensemble3&#39;] = 2.0</span>
<span class="sd">```</span>
2021-11-15 10:43:47 +00:00
<span class="sd">In case the `gamma_method` is called without any parameters it will use the values specified in the dictionaries for the respective ensembles.</span>
<span class="sd">Passing arguments to the `gamma_method` still dominates over the dictionaries.</span>
2021-11-07 20:53:18 +00:00
<span class="sd">## Irregular Monte Carlo chains</span>
2022-01-21 14:55:18 +00:00
<span class="sd">`Obs` objects defined on irregular Monte Carlo chains can be initialized with the parameter `idl`.</span>
2021-11-08 15:04:55 +00:00
2021-11-08 14:53:27 +00:00
<span class="sd">```python</span>
2021-11-08 15:04:55 +00:00
<span class="sd"># Observable defined on configurations 20 to 519</span>
<span class="sd">obs1 = pe.Obs([samples1], [&#39;ensemble1&#39;], idl=[range(20, 520)])</span>
2021-11-16 12:00:28 +00:00
<span class="sd">obs1.details()</span>
<span class="sd">&gt; Result 9.98319881e-01</span>
<span class="sd">&gt; 500 samples in 1 ensemble:</span>
<span class="sd">&gt; · Ensemble &#39;ensemble1&#39; : 500 configurations (from 20 to 519)</span>
2021-11-08 15:10:26 +00:00
2021-11-08 15:04:55 +00:00
<span class="sd"># Observable defined on every second configuration between 5 and 1003</span>
<span class="sd">obs2 = pe.Obs([samples2], [&#39;ensemble1&#39;], idl=[range(5, 1005, 2)])</span>
2021-11-16 12:00:28 +00:00
<span class="sd">obs2.details()</span>
<span class="sd">&gt; Result 9.99100712e-01</span>
<span class="sd">&gt; 500 samples in 1 ensemble:</span>
<span class="sd">&gt; · Ensemble &#39;ensemble1&#39; : 500 configurations (from 5 to 1003 in steps of 2)</span>
2021-11-08 15:10:26 +00:00
2021-11-08 15:04:55 +00:00
<span class="sd"># Observable defined on configurations 2, 9, 28, 29 and 501</span>
<span class="sd">obs3 = pe.Obs([samples3], [&#39;ensemble1&#39;], idl=[[2, 9, 28, 29, 501]])</span>
2021-11-16 12:00:28 +00:00
<span class="sd">obs3.details()</span>
<span class="sd">&gt; Result 1.01718064e+00</span>
<span class="sd">&gt; 5 samples in 1 ensemble:</span>
<span class="sd">&gt; · Ensemble &#39;ensemble1&#39; : 5 configurations (irregular range)</span>
2021-11-08 15:04:55 +00:00
<span class="sd">```</span>
2021-11-08 14:53:27 +00:00
2022-01-21 14:55:18 +00:00
<span class="sd">`Obs` objects defined on regular and irregular histories of the same ensemble can be computed with each other and the correct error propagation and estimation is automatically taken care of.</span>
2021-11-08 15:04:55 +00:00
<span class="sd">**Warning:** Irregular Monte Carlo chains can result in odd patterns in the autocorrelation functions.</span>
2021-11-15 10:57:31 +00:00
<span class="sd">Make sure to check the autocorrelation time with e.g. `pyerrors.obs.Obs.plot_rho` or `pyerrors.obs.Obs.plot_tauint`.</span>
2021-11-08 14:53:27 +00:00
2021-12-11 23:07:13 +00:00
<span class="sd">For the full API see `pyerrors.obs.Obs`.</span>
2021-11-07 20:53:18 +00:00
2021-11-07 21:09:48 +00:00
<span class="sd"># Correlators</span>
2022-02-15 12:25:34 +00:00
<span class="sd">When one is not interested in single observables but correlation functions, `pyerrors` offers the `Corr` class which simplifies the corresponding error propagation and provides the user with a set of standard methods. In order to initialize a `Corr` objects one needs to arrange the data as a list of `Obs`</span>
2021-12-23 15:43:56 +00:00
<span class="sd">```python</span>
<span class="sd">my_corr = pe.Corr([obs_0, obs_1, obs_2, obs_3])</span>
<span class="sd">print(my_corr)</span>
<span class="sd">&gt; x0/a Corr(x0/a)</span>
<span class="sd">&gt; ------------------</span>
<span class="sd">&gt; 0 0.7957(80)</span>
<span class="sd">&gt; 1 0.5156(51)</span>
<span class="sd">&gt; 2 0.3227(33)</span>
<span class="sd">&gt; 3 0.2041(21)</span>
<span class="sd">```</span>
<span class="sd">In case the correlation functions are not defined on the outermost timeslices, for example because of fixed boundary conditions, a padding can be introduced.</span>
<span class="sd">```python</span>
2022-01-18 15:09:34 +00:00
<span class="sd">my_corr = pe.Corr([obs_0, obs_1, obs_2, obs_3], padding=[1, 1])</span>
2021-12-23 15:43:56 +00:00
<span class="sd">print(my_corr)</span>
<span class="sd">&gt; x0/a Corr(x0/a)</span>
<span class="sd">&gt; ------------------</span>
<span class="sd">&gt; 0</span>
<span class="sd">&gt; 1 0.7957(80)</span>
<span class="sd">&gt; 2 0.5156(51)</span>
<span class="sd">&gt; 3 0.3227(33)</span>
<span class="sd">&gt; 4 0.2041(21)</span>
<span class="sd">&gt; 5</span>
<span class="sd">```</span>
<span class="sd">The individual entries of a correlator can be accessed via slicing</span>
<span class="sd">```python</span>
<span class="sd">print(my_corr[3])</span>
<span class="sd">&gt; 0.3227(33)</span>
<span class="sd">```</span>
<span class="sd">Error propagation with the `Corr` class works very similar to `Obs` objects. Mathematical operations are overloaded and `Corr` objects can be computed together with other `Corr` objects, `Obs` objects or real numbers and integers.</span>
<span class="sd">```python</span>
<span class="sd">my_new_corr = 0.3 * my_corr[2] * my_corr * my_corr + 12 / my_corr</span>
<span class="sd">```</span>
<span class="sd">`pyerrors` provides the user with a set of regularly used methods for the manipulation of correlator objects:</span>
<span class="sd">- `Corr.gamma_method` applies the gamma method to all entries of the correlator.</span>
<span class="sd">- `Corr.m_eff` to construct effective masses. Various variants for periodic and fixed temporal boundary conditions are available.</span>
<span class="sd">- `Corr.deriv` returns the first derivative of the correlator as `Corr`. Different discretizations of the numerical derivative are available.</span>
<span class="sd">- `Corr.second_deriv` returns the second derivative of the correlator as `Corr`. Different discretizations of the numerical derivative are available.</span>
<span class="sd">- `Corr.symmetric` symmetrizes parity even correlations functions, assuming periodic boundary conditions.</span>
<span class="sd">- `Corr.anti_symmetric` anti-symmetrizes parity odd correlations functions, assuming periodic boundary conditions.</span>
<span class="sd">- `Corr.T_symmetry` averages a correlator with its time symmetry partner, assuming fixed boundary conditions.</span>
<span class="sd">- `Corr.plateau` extracts a plateau value from the correlator in a given range.</span>
<span class="sd">- `Corr.roll` periodically shifts the correlator.</span>
<span class="sd">- `Corr.reverse` reverses the time ordering of the correlator.</span>
<span class="sd">- `Corr.correlate` constructs a disconnected correlation function from the correlator and another `Corr` or `Obs` object.</span>
<span class="sd">- `Corr.reweight` reweights the correlator.</span>
2022-01-21 14:55:18 +00:00
<span class="sd">`pyerrors` can also handle matrices of correlation functions and extract energy states from these matrices via a generalized eigenvalue problem (see `pyerrors.correlators.Corr.GEVP`).</span>
2021-12-23 15:43:56 +00:00
2021-12-11 23:07:13 +00:00
<span class="sd">For the full API see `pyerrors.correlators.Corr`.</span>
2021-11-09 11:51:27 +00:00
2022-02-14 14:08:21 +00:00
<span class="sd"># Complex valued observables</span>
2021-11-15 14:59:01 +00:00
<span class="sd">`pyerrors` can handle complex valued observables via the class `pyerrors.obs.CObs`.</span>
<span class="sd">`CObs` are initialized with a real and an imaginary part which both can be `Obs` valued.</span>
<span class="sd">```python</span>
<span class="sd">my_real_part = pe.Obs([samples1], [&#39;ensemble1&#39;])</span>
<span class="sd">my_imag_part = pe.Obs([samples2], [&#39;ensemble1&#39;])</span>
<span class="sd">my_cobs = pe.CObs(my_real_part, my_imag_part)</span>
<span class="sd">my_cobs.gamma_method()</span>
<span class="sd">print(my_cobs)</span>
<span class="sd">&gt; (0.9959(91)+0.659(28)j)</span>
<span class="sd">```</span>
<span class="sd">Elementary mathematical operations are overloaded and samples are properly propagated as for the `Obs` class.</span>
<span class="sd">```python</span>
<span class="sd">my_derived_cobs = (my_cobs + my_cobs.conjugate()) / np.abs(my_cobs)</span>
<span class="sd">my_derived_cobs.gamma_method()</span>
<span class="sd">print(my_derived_cobs)</span>
<span class="sd">&gt; (1.668(23)+0.0j)</span>
<span class="sd">```</span>
2021-11-07 21:09:48 +00:00
2022-02-14 14:08:21 +00:00
<span class="sd"># Error propagation in iterative algorithms</span>
<span class="sd">`pyerrors` supports exact linear error propagation for iterative algorithms like various variants of non-linear least sqaures fits or root finding. The derivatives required for the error propagation are calculated as described in [arXiv:1809.01289](https://arxiv.org/abs/1809.01289).</span>
<span class="sd">## Least squares fits</span>
<span class="sd">Standard non-linear least square fits with errors on the dependent but not the independent variables can be performed with `pyerrors.fits.least_squares`. As default solver the Levenberg-Marquardt algorithm implemented in [scipy](https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.least_squares.html) is used.</span>
<span class="sd">Fit functions have to be of the following form</span>
<span class="sd">```python</span>
<span class="sd">import autograd.numpy as anp</span>
<span class="sd">def func(a, x):</span>
<span class="sd"> return a[1] * anp.exp(-a[0] * x)</span>
<span class="sd">```</span>
2022-02-21 14:52:31 +00:00
<span class="sd">**It is important that numerical functions refer to `autograd.numpy` instead of `numpy` for the automatic differentiation in iterative algorithms to work properly.**</span>
2022-02-14 14:08:21 +00:00
<span class="sd">Fits can then be performed via</span>
<span class="sd">```python</span>
<span class="sd">fit_result = pe.fits.least_squares(x, y, func)</span>
<span class="sd">print(&quot;\n&quot;, fit_result)</span>
<span class="sd">&gt; Fit with 2 parameters</span>
<span class="sd">&gt; Method: Levenberg-Marquardt</span>
<span class="sd">&gt; `ftol` termination condition is satisfied.</span>
<span class="sd">&gt; chisquare/d.o.f.: 0.9593035785160936</span>
<span class="sd">&gt; Goodness of fit:</span>
<span class="sd">&gt; χ²/d.o.f. = 0.959304</span>
<span class="sd">&gt; p-value = 0.5673</span>
<span class="sd">&gt; Fit parameters:</span>
<span class="sd">&gt; 0 0.0548(28)</span>
<span class="sd">&gt; 1 1.933(64)</span>
<span class="sd">```</span>
<span class="sd">where x is a `list` or `numpy.array` of `floats` and y is a `list` or `numpy.array` of `Obs`.</span>
<span class="sd">Data stored in `Corr` objects can be fitted directly using the `Corr.fit` method.</span>
<span class="sd">```python</span>
<span class="sd">my_corr = pe.Corr(y)</span>
<span class="sd">fit_result = my_corr.fit(func, fitrange=[12, 25])</span>
<span class="sd">```</span>
<span class="sd">this can simplify working with absolute fit ranges and takes care of gaps in the data automatically.</span>
<span class="sd">For fit functions with multiple independent variables the fit function can be of the form</span>
<span class="sd">```python</span>
<span class="sd">def func(a, x):</span>
<span class="sd"> (x1, x2) = x</span>
<span class="sd"> return a[0] * x1 ** 2 + a[1] * x2</span>
<span class="sd">```</span>
<span class="sd">## Total least squares fits</span>
<span class="sd">`pyerrors` can also fit data with errors on both the dependent and independent variables using the total least squares method also referred to orthogonal distance regression as implemented in [scipy](https://docs.scipy.org/doc/scipy/reference/odr.html), see `pyerrors.fits.least_squares`. The syntax is identical to the standard least squares case, the only diffrence being that `x` also has to be a `list` or `numpy.array` of `Obs`.</span>
<span class="sd">For the full API see `pyerrors.fits` for fits and `pyerrors.roots` for finding roots of functions.</span>
2021-11-07 21:09:48 +00:00
2021-11-07 20:53:18 +00:00
<span class="sd"># Matrix operations</span>
2022-02-14 14:26:48 +00:00
<span class="sd">`pyerrors` provides wrappers for `Obs`- and `CObs`-valued matrix operations based on `numpy.linalg`. The supported functions include:</span>
2022-01-21 14:55:18 +00:00
<span class="sd">- `inv` for the matrix inverse.</span>
<span class="sd">- `cholseky` for the Cholesky decomposition.</span>
<span class="sd">- `det` for the matrix determinant.</span>
<span class="sd">- `eigh` for eigenvalues and eigenvectors of hermitean matrices.</span>
<span class="sd">- `eig` for eigenvalues of general matrices.</span>
<span class="sd">- `pinv` for the Moore-Penrose pseudoinverse.</span>
<span class="sd">- `svd` for the singular-value-decomposition.</span>
<span class="sd">For the full API see `pyerrors.linalg`.</span>
2021-11-07 20:53:18 +00:00
2021-11-15 14:30:41 +00:00
<span class="sd"># Export data</span>
2022-02-14 14:19:54 +00:00
2022-02-21 14:52:31 +00:00
<span class="sd">The preferred exported file format within `pyerrors` is json.gz. The exact specifications of this format will be listed here soon.</span>
2021-11-15 14:30:41 +00:00
<span class="sd">## Jackknife samples</span>
2021-12-11 23:07:13 +00:00
<span class="sd">For comparison with other analysis workflows `pyerrors` can generate jackknife samples from an `Obs` object or import jackknife samples into an `Obs` object.</span>
<span class="sd">See `pyerrors.obs.Obs.export_jackknife` and `pyerrors.obs.import_jackknife` for details.</span>
2021-11-15 14:30:41 +00:00
2022-02-14 14:08:21 +00:00
<span class="sd"># Citing</span>
<span class="sd">If you use `pyerrors` for research that leads to a publication please consider citing:</span>
2022-02-14 14:19:54 +00:00
<span class="sd">- Ulli Wolff, *Monte Carlo errors with less errors*. Comput.Phys.Commun. 156 (2004) 143-153, Comput.Phys.Commun. 176 (2007) 383 (erratum).</span>
<span class="sd">- Stefan Schaefer, Rainer Sommer, Francesco Virotta, *Critical slowing down and error analysis in lattice QCD simulations*. Nucl.Phys.B 845 (2011) 93-119.</span>
<span class="sd">- Alberto Ramos, *Automatic differentiation for error analysis of Monte Carlo data*. Comput.Phys.Commun. 238 (2019) 19-35.</span>
2021-11-07 20:53:18 +00:00
<span class="sd">&#39;&#39;&#39;</span>
<span class="kn">from</span> <span class="nn">.obs</span> <span class="kn">import</span> <span class="o">*</span>
<span class="kn">from</span> <span class="nn">.correlators</span> <span class="kn">import</span> <span class="o">*</span>
<span class="kn">from</span> <span class="nn">.fits</span> <span class="kn">import</span> <span class="o">*</span>
2021-11-24 17:27:16 +00:00
<span class="kn">from</span> <span class="nn">.misc</span> <span class="kn">import</span> <span class="o">*</span>
2021-11-07 20:53:18 +00:00
<span class="kn">from</span> <span class="nn">.</span> <span class="kn">import</span> <span class="n">dirac</span>
2021-11-11 10:51:37 +00:00
<span class="kn">from</span> <span class="nn">.</span> <span class="kn">import</span> <span class="nb">input</span>
2021-11-07 20:53:18 +00:00
<span class="kn">from</span> <span class="nn">.</span> <span class="kn">import</span> <span class="n">linalg</span>
<span class="kn">from</span> <span class="nn">.</span> <span class="kn">import</span> <span class="n">mpm</span>
<span class="kn">from</span> <span class="nn">.</span> <span class="kn">import</span> <span class="n">roots</span>
<span class="kn">from</span> <span class="nn">.version</span> <span class="kn">import</span> <span class="n">__version__</span>
</pre></div>
</details>
</section>
</main>
<script>
function escapeHTML(html) {
return document.createElement('div').appendChild(document.createTextNode(html)).parentNode.innerHTML;
}
const originalContent = document.querySelector("main.pdoc");
let currentContent = originalContent;
function setContent(innerHTML) {
let elem;
if (innerHTML) {
elem = document.createElement("main");
elem.classList.add("pdoc");
elem.innerHTML = innerHTML;
} else {
elem = originalContent;
}
if (currentContent !== elem) {
currentContent.replaceWith(elem);
currentContent = elem;
}
}
function getSearchTerm() {
return (new URL(window.location)).searchParams.get("search");
}
const searchBox = document.querySelector(".pdoc input[type=search]");
searchBox.addEventListener("input", function () {
let url = new URL(window.location);
if (searchBox.value.trim()) {
url.hash = "";
url.searchParams.set("search", searchBox.value);
} else {
url.searchParams.delete("search");
}
history.replaceState("", "", url.toString());
onInput();
});
window.addEventListener("popstate", onInput);
let search, searchErr;
async function initialize() {
try {
search = await new Promise((resolve, reject) => {
const script = document.createElement("script");
script.type = "text/javascript";
script.async = true;
script.onload = () => resolve(window.pdocSearch);
script.onerror = (e) => reject(e);
script.src = "search.js";
document.getElementsByTagName("head")[0].appendChild(script);
});
} catch (e) {
console.error("Cannot fetch pdoc search index");
searchErr = "Cannot fetch search index.";
}
onInput();
document.querySelector("nav.pdoc").addEventListener("click", e => {
if (e.target.hash) {
searchBox.value = "";
searchBox.dispatchEvent(new Event("input"));
}
});
}
function onInput() {
setContent((() => {
const term = getSearchTerm();
if (!term) {
return null
}
if (searchErr) {
return `<h3>Error: ${searchErr}</h3>`
}
if (!search) {
return "<h3>Searching...</h3>"
}
window.scrollTo({top: 0, left: 0, behavior: 'auto'});
const results = search(term);
let html;
if (results.length === 0) {
html = `No search results for '${escapeHTML(term)}'.`
} else {
html = `<h4>${results.length} search result${results.length > 1 ? "s" : ""} for '${escapeHTML(term)}'.</h4>`;
}
for (let result of results.slice(0, 10)) {
let doc = result.doc;
let url = `${doc.modulename.replaceAll(".", "/")}.html`;
if (doc.qualname) {
url += `#${doc.qualname}`;
}
let heading;
switch (result.doc.type) {
case "function":
2022-01-16 15:44:50 +00:00
heading = `<span class="def">${doc.funcdef}</span> <span class="name">${doc.fullname}</span><span class="signature">${doc.signature}:</span>`;
2021-11-07 20:53:18 +00:00
break;
case "class":
heading = `<span class="def">class</span> <span class="name">${doc.fullname}</span>`;
2022-01-16 15:44:50 +00:00
if (doc.bases)
heading += `<wbr>(<span class="base">${doc.bases}</span>)`;
heading += `:`;
break;
case "variable":
heading = `<span class="name">${doc.fullname}</span>`;
if (doc.annotation)
heading += `<span class="annotation">${doc.annotation}</span>`;
if (doc.default_value)
heading += `<span class="default_value">${doc.default_value}</span>`;
2021-11-07 20:53:18 +00:00
break;
default:
heading = `<span class="name">${doc.fullname}</span>`;
break;
}
html += `
<section class="search-result">
<a href="${url}" class="attr ${doc.type}">${heading}</a>
<div class="docstring">${doc.doc}</div>
</section>
`;
}
return html;
})());
}
if (getSearchTerm()) {
initialize();
searchBox.value = getSearchTerm();
onInput();
} else {
searchBox.addEventListener("focus", initialize, {once: true});
}
searchBox.addEventListener("keydown", e => {
if (["ArrowDown", "ArrowUp", "Enter"].includes(e.key)) {
let focused = currentContent.querySelector(".search-result.focused");
if (!focused) {
currentContent.querySelector(".search-result").classList.add("focused");
} else if (
e.key === "ArrowDown"
&& focused.nextElementSibling
&& focused.nextElementSibling.classList.contains("search-result")
) {
focused.classList.remove("focused");
focused.nextElementSibling.classList.add("focused");
focused.nextElementSibling.scrollIntoView({
behavior: "smooth",
block: "nearest",
inline: "nearest"
});
} else if (
e.key === "ArrowUp"
&& focused.previousElementSibling
&& focused.previousElementSibling.classList.contains("search-result")
) {
focused.classList.remove("focused");
focused.previousElementSibling.classList.add("focused");
focused.previousElementSibling.scrollIntoView({
behavior: "smooth",
block: "nearest",
inline: "nearest"
});
} else if (
e.key === "Enter"
) {
focused.querySelector("a").click();
}
}
});
</script></body>
</html>