tom/personal_site
1
0
Fork 0
mirror of https://github.com/TomHodson/tomhodson.github.com.git synced 2025-06-26 10:01:18 +02:00
Code Issues Packages Projects Releases Wiki Activity
personal_site/_thesis/1.1_FK_Intro.html
Tom Hodson e9d85bfa1e updates!
2022-08-08 11:59:36 +01:00

1651 lines
90 KiB
HTML
Raw Blame History

This file contains invisible Unicode characters

This file contains invisible Unicode characters that are indistinguishable to humans but may be processed differently by a computer. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

---
title: The Falikov-Kimball Model - Introduction
excerpt: The Falikov-Kimball is about the simplest possible testbed we could have for the many electron problem.
layout: none
image:
---
<!DOCTYPE html>
<html xmlns="http://www.w3.org/1999/xhtml" lang="" xml:lang="">
<head>
<meta charset="utf-8" />
<meta name="generator" content="pandoc" />
<meta name="viewport" content="width=device-width, initial-scale=1.0, user-scalable=yes" />
<meta name="description" content="The Falikov-Kimball is about the simplest possible testbed we could have for the many electron problem." />
<title>The Falikov-Kimball Model - Introduction</title>
<!-- <style>
html {
line-height: 1.5;
font-family: Georgia, serif;
font-size: 20px;
color: #1a1a1a;
background-color: #fdfdfd;
}
body {
margin: 0 auto;
max-width: 36em;
padding-left: 50px;
padding-right: 50px;
padding-top: 50px;
padding-bottom: 50px;
hyphens: auto;
overflow-wrap: break-word;
text-rendering: optimizeLegibility;
font-kerning: normal;
}
@media (max-width: 600px) {
body {
font-size: 0.9em;
padding: 1em;
}
h1 {
font-size: 1.8em;
}
}
@media print {
body {
background-color: transparent;
color: black;
font-size: 12pt;
}
p, h2, h3 {
orphans: 3;
widows: 3;
}
h2, h3, h4 {
page-break-after: avoid;
}
}
p {
margin: 1em 0;
}
a {
color: #1a1a1a;
}
a:visited {
color: #1a1a1a;
}
img {
max-width: 100%;
}
h1, h2, h3, h4, h5, h6 {
margin-top: 1.4em;
}
h5, h6 {
font-size: 1em;
font-style: italic;
}
h6 {
font-weight: normal;
}
ol, ul {
padding-left: 1.7em;
margin-top: 1em;
}
li > ol, li > ul {
margin-top: 0;
}
blockquote {
margin: 1em 0 1em 1.7em;
padding-left: 1em;
border-left: 2px solid #e6e6e6;
color: #606060;
}
code {
font-family: Menlo, Monaco, 'Lucida Console', Consolas, monospace;
font-size: 85%;
margin: 0;
}
pre {
margin: 1em 0;
overflow: auto;
}
pre code {
padding: 0;
overflow: visible;
overflow-wrap: normal;
}
.sourceCode {
background-color: transparent;
overflow: visible;
}
hr {
background-color: #1a1a1a;
border: none;
height: 1px;
margin: 1em 0;
}
table {
margin: 1em 0;
border-collapse: collapse;
width: 100%;
overflow-x: auto;
display: block;
font-variant-numeric: lining-nums tabular-nums;
}
table caption {
margin-bottom: 0.75em;
}
tbody {
margin-top: 0.5em;
border-top: 1px solid #1a1a1a;
border-bottom: 1px solid #1a1a1a;
}
th {
border-top: 1px solid #1a1a1a;
padding: 0.25em 0.5em 0.25em 0.5em;
}
td {
padding: 0.125em 0.5em 0.25em 0.5em;
}
header {
margin-bottom: 4em;
text-align: center;
}
#TOC li {
list-style: none;
}
#TOC ul {
padding-left: 1.3em;
}
#TOC > ul {
padding-left: 0;
}
#TOC a:not(:hover) {
text-decoration: none;
}
code{white-space: pre-wrap;}
span.smallcaps{font-variant: small-caps;}
span.underline{text-decoration: underline;}
div.column{display: inline-block; vertical-align: top; width: 50%;}
div.hanging-indent{margin-left: 1.5em; text-indent: -1.5em;}
ul.task-list{list-style: none;}
pre > code.sourceCode { white-space: pre; position: relative; }
pre > code.sourceCode > span { display: inline-block; line-height: 1.25; }
pre > code.sourceCode > span:empty { height: 1.2em; }
.sourceCode { overflow: visible; }
code.sourceCode > span { color: inherit; text-decoration: inherit; }
div.sourceCode { margin: 1em 0; }
pre.sourceCode { margin: 0; }
@media screen {
div.sourceCode { overflow: auto; }
}
@media print {
pre > code.sourceCode { white-space: pre-wrap; }
pre > code.sourceCode > span { text-indent: -5em; padding-left: 5em; }
}
pre.numberSource code
{ counter-reset: source-line 0; }
pre.numberSource code > span
{ position: relative; left: -4em; counter-increment: source-line; }
pre.numberSource code > span > a:first-child::before
{ content: counter(source-line);
position: relative; left: -1em; text-align: right; vertical-align: baseline;
border: none; display: inline-block;
-webkit-touch-callout: none; -webkit-user-select: none;
-khtml-user-select: none; -moz-user-select: none;
-ms-user-select: none; user-select: none;
padding: 0 4px; width: 4em;
color: #aaaaaa;
}
pre.numberSource { margin-left: 3em; border-left: 1px solid #aaaaaa; padding-left: 4px; }
div.sourceCode
{ }
@media screen {
pre > code.sourceCode > span > a:first-child::before { text-decoration: underline; }
}
code span.al { color: #ff0000; font-weight: bold; } /* Alert */
code span.an { color: #60a0b0; font-weight: bold; font-style: italic; } /* Annotation */
code span.at { color: #7d9029; } /* Attribute */
code span.bn { color: #40a070; } /* BaseN */
code span.bu { } /* BuiltIn */
code span.cf { color: #007020; font-weight: bold; } /* ControlFlow */
code span.ch { color: #4070a0; } /* Char */
code span.cn { color: #880000; } /* Constant */
code span.co { color: #60a0b0; font-style: italic; } /* Comment */
code span.cv { color: #60a0b0; font-weight: bold; font-style: italic; } /* CommentVar */
code span.do { color: #ba2121; font-style: italic; } /* Documentation */
code span.dt { color: #902000; } /* DataType */
code span.dv { color: #40a070; } /* DecVal */
code span.er { color: #ff0000; font-weight: bold; } /* Error */
code span.ex { } /* Extension */
code span.fl { color: #40a070; } /* Float */
code span.fu { color: #06287e; } /* Function */
code span.im { } /* Import */
code span.in { color: #60a0b0; font-weight: bold; font-style: italic; } /* Information */
code span.kw { color: #007020; font-weight: bold; } /* Keyword */
code span.op { color: #666666; } /* Operator */
code span.ot { color: #007020; } /* Other */
code span.pp { color: #bc7a00; } /* Preprocessor */
code span.sc { color: #4070a0; } /* SpecialChar */
code span.ss { color: #bb6688; } /* SpecialString */
code span.st { color: #4070a0; } /* String */
code span.va { color: #19177c; } /* Variable */
code span.vs { color: #4070a0; } /* VerbatimString */
code span.wa { color: #60a0b0; font-weight: bold; font-style: italic; } /* Warning */
div.csl-bib-body { }
div.csl-entry {
clear: both;
}
.hanging div.csl-entry {
margin-left:2em;
text-indent:-2em;
}
div.csl-left-margin {
min-width:2em;
float:left;
}
div.csl-right-inline {
margin-left:2em;
padding-left:1em;
}
div.csl-indent {
margin-left: 2em;
}
</style> -->
<!-- <script
src="https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-chtml-full.js"
type="text/javascript"></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.0.1/es5/tex-mml-chtml.js"></script>
-->
<script src="/assets/mathjax/tex-mml-svg.js" id="MathJax-script" async></script>
<!--[if lt IE 9]>
<script src="//cdnjs.cloudflare.com/ajax/libs/html5shiv/3.7.3/html5shiv-printshiv.min.js"></script>
<![endif]-->
<link rel="stylesheet" href="/assets/css/styles.css">
<script src="/assets/js/index.js"></script>
</head>
<body>
{% include header.html %}
<main>
<nav id="TOC" role="doc-toc">
<ul>
<li><a href="#contributions"
id="toc-contributions">Contributions</a></li>
<li><a href="#introduction" id="toc-introduction">Introduction</a>
<ul>
<li><a href="#localisation" id="toc-localisation">Localisation</a>
<ul>
<li><a href="#the-falikov-kimball-model"
id="toc-the-falikov-kimball-model">The Falikov Kimball Model</a></li>
</ul></li>
<li><a href="#falikov-kimball-and-hubbard-models"
id="toc-falikov-kimball-and-hubbard-models">Falikov Kimball and Hubbard
models</a>
<ul>
<li><a href="#hubbard-model" id="toc-hubbard-model">Hubbard
model</a></li>
<li><a href="#falikov-kimball-model"
id="toc-falikov-kimball-model">Falikov-Kimball model</a></li>
<li><a href="#thermodynamics-of-the-fk-model"
id="toc-thermodynamics-of-the-fk-model">Thermodynamics of the FK
model</a></li>
<li><a href="#thermodynamics"
id="toc-thermodynamics">Thermodynamics</a></li>
<li><a href="#markov-chain-monte-carlo"
id="toc-markov-chain-monte-carlo">Markov Chain Monte Carlo</a></li>
</ul></li>
<li><a href="#localisation-1" id="toc-localisation-1">Localisation</a>
<ul>
<li><a href="#thermalisation"
id="toc-thermalisation">Thermalisation</a></li>
<li><a href="#anderson-localisation"
id="toc-anderson-localisation">Anderson Localisation</a></li>
<li><a href="#many-body-localisation"
id="toc-many-body-localisation">Many Body Localisation</a></li>
<li><a href="#disorder-free-localisation"
id="toc-disorder-free-localisation">Disorder Free localisation</a></li>
<li><a href="#diagnostics-of-localisation"
id="toc-diagnostics-of-localisation">Diagnostics of
Localisation</a></li>
</ul></li>
<li><a href="#numerical-methods" id="toc-numerical-methods">Numerical
Methods</a>
<ul>
<li><a href="#markov-chain-monte-carlo-1"
id="toc-markov-chain-monte-carlo-1">Markov Chain Monte Carlo}</a></li>
<li><a href="#applying-mcmc-to-the-fk-model"
id="toc-applying-mcmc-to-the-fk-model">Applying MCMC to the FK
model}</a></li>
</ul></li>
<li><a href="#markov-chain-monte-carlo-in-practice"
id="toc-markov-chain-monte-carlo-in-practice">Markov Chain Monte-Carlo
in Practice}</a>
<ul>
<li><a href="#quick-intro-to-mcmc" id="toc-quick-intro-to-mcmc">Quick
Intro to MCMC}</a></li>
<li><a href="#convergence-time" id="toc-convergence-time">Convergence
Time}</a></li>
<li><a href="#auto-correlation-time"
id="toc-auto-correlation-time">Auto-correlation Time}</a></li>
<li><a href="#the-metropolis-hastings-algorithm"
id="toc-the-metropolis-hastings-algorithm">The Metropolis-Hastings
Algorithm}</a></li>
<li><a href="#choosing-the-proposal-distribution"
id="toc-choosing-the-proposal-distribution">Choosing the proposal
distribution}</a></li>
<li><a href="#two-step-trick" id="toc-two-step-trick">Two Step
Trick</a></li>
</ul></li>
</ul></li>
</ul>
</nav>
<h1 id="contributions">Contributions</h1>
<p>This material is this chapter expands on work presented in</p>
<p><span class="citation" data-cites="citekey"><sup><a
href="#ref-citekey"
role="doc-biblioref"><strong>citekey?</strong></a></sup></span> <a
href="https://link.aps.org/doi/10.1103/PhysRevB.104.045116">One-dimensional
long-range Falikov-Kimball model: Thermal phase transition and
disorder-free localization</a>, Hodson, T. and Willsher, J. and Knolle,
J., Phys. Rev. B, <strong>104</strong>, 4, 2021,</p>
<p>Johannes had the initial idea to use a long range Ising term to
stablise order in a one dimension Falikov-Kimball model. Josef developed
a proof of concept during a summer project at Imperial. The three of us
brought the project to fruition.</p>
<div class="sourceCode" id="cb1"><pre
class="sourceCode python"><code class="sourceCode python"><span id="cb1-1"><a href="#cb1-1" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation abanin_recent_2017 <span class="kw">not</span> found abaninRecentProgressManybody2017</span>
<span id="cb1-2"><a href="#cb1-2" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation anderson_absence_1958<span class="op">-</span><span class="dv">1</span> <span class="kw">not</span> found andersonAbsenceDiffusionCertain1958</span>
<span id="cb1-3"><a href="#cb1-3" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation antipov_interaction<span class="op">-</span>tuned_2016<span class="op">-</span><span class="dv">1</span> <span class="kw">not</span> found andersonAbsenceDiffusionCertain1958</span>
<span id="cb1-4"><a href="#cb1-4" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation binder_finite_1981 <span class="kw">not</span> found binderFiniteSizeScaling1981</span>
<span id="cb1-5"><a href="#cb1-5" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation croy_anderson_2011 <span class="kw">not</span> found croyAndersonLocalization1D2011</span>
<span id="cb1-6"><a href="#cb1-6" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation dalessio_quantum_2016 <span class="kw">not</span> found dalessioQuantumChaosEigenstate2016</span>
<span id="cb1-7"><a href="#cb1-7" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation dyson_existence_1969 <span class="kw">not</span> found dysonExistencePhasetransitionOnedimensional1969</span>
<span id="cb1-8"><a href="#cb1-8" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation fig:binder <span class="kw">not</span> found</span>
<span id="cb1-9"><a href="#cb1-9" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation goldshtein_pure_1977 <span class="kw">not</span> found goldshteinPurePointSpectrum1977</span>
<span id="cb1-10"><a href="#cb1-10" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation huang_accelerated_2017 <span class="kw">not</span> found huangAcceleratedMonteCarlo2017</span>
<span id="cb1-11"><a href="#cb1-11" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation hubbard_j._electron_1963 <span class="kw">not</span> found</span>
<span id="cb1-12"><a href="#cb1-12" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation imbrie_diagonalization_2016 <span class="kw">not</span> found</span>
<span id="cb1-13"><a href="#cb1-13" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation imbrie_many<span class="op">-</span>body_2016 <span class="kw">not</span> found</span>
<span id="cb1-14"><a href="#cb1-14" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation izrailev_anomalous_2012 <span class="kw">not</span> found</span>
<span id="cb1-15"><a href="#cb1-15" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation izrailev_localization_1999 <span class="kw">not</span> found</span>
<span id="cb1-16"><a href="#cb1-16" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation kapfer_sampling_2013 <span class="kw">not</span> found</span>
<span id="cb1-17"><a href="#cb1-17" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation kennedy_itinerant_1986 <span class="kw">not</span> found</span>
<span id="cb1-18"><a href="#cb1-18" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation khatami_fluctuation<span class="op">-</span>dissipation_2013 <span class="kw">not</span> found</span>
<span id="cb1-19"><a href="#cb1-19" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation kramer_localization_1993 <span class="kw">not</span> found</span>
<span id="cb1-20"><a href="#cb1-20" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation krauth_introduction_1996 <span class="kw">not</span> found</span>
<span id="cb1-21"><a href="#cb1-21" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation lieb_absence_1968 <span class="kw">not</span> found</span>
<span id="cb1-22"><a href="#cb1-22" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation lipkin_validity_1965 <span class="kw">not</span> found</span>
<span id="cb1-23"><a href="#cb1-23" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation maska_thermodynamics_2006<span class="op">-</span><span class="dv">1</span> <span class="kw">not</span> found</span>
<span id="cb1-24"><a href="#cb1-24" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation musial_monte_2002 <span class="kw">not</span> found</span>
<span id="cb1-25"><a href="#cb1-25" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation nagaoka_ferromagnetism_1966 <span class="kw">not</span> found</span>
<span id="cb1-26"><a href="#cb1-26" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation noauthor_hubbard_2013 <span class="kw">not</span> found</span>
<span id="cb1-27"><a href="#cb1-27" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation peierls_isings_1936 <span class="kw">not</span> found</span>
<span id="cb1-28"><a href="#cb1-28" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation roberts_weak_1997 <span class="kw">not</span> found</span>
<span id="cb1-29"><a href="#cb1-29" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation ruelle_statistical_1968 <span class="kw">not</span> found</span>
<span id="cb1-30"><a href="#cb1-30" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation schiulaz_dynamics_2015 <span class="kw">not</span> found</span>
<span id="cb1-31"><a href="#cb1-31" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation smith_disorder<span class="op">-</span>free_2017 <span class="kw">not</span> found</span>
<span id="cb1-32"><a href="#cb1-32" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation smith_dynamical_2018 <span class="kw">not</span> found</span>
<span id="cb1-33"><a href="#cb1-33" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation srednicki_chaos_1994 <span class="kw">not</span> found</span>
<span id="cb1-34"><a href="#cb1-34" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation thouless_long<span class="op">-</span>range_1969 <span class="kw">not</span> found</span>
<span id="cb1-35"><a href="#cb1-35" aria-hidden="true" tabindex="-1"></a>[WARNING] Citeproc: citation yao_quasi<span class="op">-</span>many<span class="op">-</span>body_2016 <span class="kw">not</span> found</span></code></pre></div>
<h1 id="introduction">Introduction</h1>
<h2 id="localisation">Localisation</h2>
<p>The discovery of localisation in quantum systems surprising at the
time given the seeming ubiquity of extended Bloch states. Later, when
thermalisation in quantum systems gained interest, localisation
phenomena again stood out as counterexamples to the eigenstate
thermalisation hypothesis<span class="citation"
data-cites="abanin_recent_2017 srednicki_chaos_1994"><sup><a
href="#ref-abanin_recent_2017"
role="doc-biblioref"><strong>abanin_recent_2017?</strong></a>,<a
href="#ref-srednicki_chaos_1994"
role="doc-biblioref"><strong>srednicki_chaos_1994?</strong></a></sup></span>,
allowing quantum systems to avoid to retain memory of their initial
conditions in the face of thermal noise.</p>
<p>The simplest and first discovered kind is Anderson localisation,
first studied in 1958<span class="citation"
data-cites="anderson_absence_1958-1"><sup><a
href="#ref-anderson_absence_1958-1"
role="doc-biblioref"><strong>anderson_absence_1958-1?</strong></a></sup></span>
in the context of non-interacting fermions subject to a static or
quenched disorder potential <span class="math inline">\(V_j\)</span>
drawn uniformly from the interval <span
class="math inline">\([-W,W]\)</span></p>
<p><span class="math display">\[
H = -t\sum_{\langle jk \rangle} c^\daggerger_j c_k + \sum_j V_j
c_j^\daggerger c_j
\]</span></p>
<p>this model exhibits exponentially localised eigenfunctions <span
class="math inline">\(\psi(x) = f(x) e^{-x/\lambda}\)</span> which
cannot contribute to transport processes. Initially it was thought that
in one dimensional disordered models, all states would be localised,
however it was later shown that in the presence of correlated disorder,
bands of extended states can exist<span class="citation"
data-cites="izrailev_localization_1999 croy_anderson_2011 izrailev_anomalous_2012"><sup><a
href="#ref-izrailev_localization_1999"
role="doc-biblioref"><strong>izrailev_localization_1999?</strong></a>,<a
href="#ref-croy_anderson_2011"
role="doc-biblioref"><strong>croy_anderson_2011?</strong></a>,<a
href="#ref-izrailev_anomalous_2012"
role="doc-biblioref"><strong>izrailev_anomalous_2012?</strong></a></sup></span>.</p>
<p>Later localisation was found in interacting many-body systems with
quenched disorder:</p>
<p><span class="math display">\[
H = -t\sum_{\langle jk \rangle} c^\daggerger_j c_k + \sum_j V_j
c_j^\daggerger c_j + U\sum_{jk} n_j n_k
\]</span></p>
<p>where the number operators <span class="math inline">\(n_j =
c^\dagger_j c_j\)</span>. Here, in contrast to the Anderson model,
localisation phenomena can be proven robust to weak perturbations of the
Hamiltonian. This is called many-body localisation (MBL)<span
class="citation" data-cites="imbrie_many-body_2016"><sup><a
href="#ref-imbrie_many-body_2016"
role="doc-biblioref"><strong>imbrie_many-body_2016?</strong></a></sup></span>.</p>
<p>Both MBL and Anderson localisation depend crucially on the presence
of quenched disorder. This has led to ongoing interest in the
possibility of disorder-free localisation, in which the disorder
necessary to generate localisation is generated entirely from the
dynamics of the model. This contracts with typical models of disordered
systems in which disorder is explicielty introduced into the Hamilton or
the initial state.</p>
<p>The concept of disorder-free localisation was first proposed in the
context of Helium mixtures<span class="citation"
data-cites="kagan1984localization"><sup><a
href="#ref-kagan1984localization"
role="doc-biblioref">1</a></sup></span> and then extended to heavy-light
mixtures in which multiple species with large mass ratios interact. The
idea is that the heavier particles act as an effective disorder
potential for the lighter ones, inducing localisation. Two such
models<span class="citation"
data-cites="yao_quasi-many-body_2016 schiulaz_dynamics_2015"><sup><a
href="#ref-yao_quasi-many-body_2016"
role="doc-biblioref"><strong>yao_quasi-many-body_2016?</strong></a>,<a
href="#ref-schiulaz_dynamics_2015"
role="doc-biblioref"><strong>schiulaz_dynamics_2015?</strong></a></sup></span>
instead find that the models thermalise exponentially slowly in system
size, which Ref.<span class="citation"
data-cites="yao_quasi-many-body_2016"><sup><a
href="#ref-yao_quasi-many-body_2016"
role="doc-biblioref"><strong>yao_quasi-many-body_2016?</strong></a></sup></span>
dubs Quasi-MBL.</p>
<p>True disorder-free localisation does occur in exactly solveable
models with extensively many conserved quantities<span class="citation"
data-cites="smith_disorder-free_2017"><sup><a
href="#ref-smith_disorder-free_2017"
role="doc-biblioref"><strong>smith_disorder-free_2017?</strong></a></sup></span>.
As conserved quantites have no time dynamics this can be thought of as
taking the separation of timescales to the infinite limit.</p>
<h3 id="the-falikov-kimball-model">The Falikov Kimball Model</h3>
<p>In the Falikov Kimball (FK) model spinless fermions <span
class="math inline">\(c_{i\uparrow}\)</span> are coupled via a repulsive
on-site interaction to a classical degree of freedom <span
class="math inline">\(n_{i\downarrow}\)</span>.</p>
<p><span class="math display">\[\begin{aligned}
H &amp;= -t \sum_{&lt;ij&gt;} c^\daggerger_{i\uparrow}c_{j\uparrow} + U
\sum_{i} (n_{i \uparrow} - 1/2)( n_{i\downarrow} - 1/2) \\
&amp; - \mu \sum_i \left( n_{i \uparrow} + n_{i \downarrow}
\right) + \sum_{ij} V_{ij} (n_{i\downarrow} - 1/2)(n_{j\downarrow} -
1/2)
\end{aligned}\]</span> <strong>replace with hamiltonian from the
paper</strong></p>
<p>This notation emphasises that this can also be thought of as an
asymmetric Hubbard model in which the spin down electrons cannot hop and
are subject to an additional long range potential. However, to avoid the
confusion of talking about two distinct species of spinless electrons we
will use a different interpretation and refer to the classical degrees
of freedom as the “ionic sector” and the quantum degrees of freedom as
the “electronic sector”. The final term that induces interactions
between the classical particles has been added by us to stabilise the
formation of an ordered phase in 1D. The classical variables commute
with the Hamiltonian <span class="math inline">\([H, n_{i\downarrow}] =
0\)</span> so like the lattice gauge model in Ref<span class="citation"
data-cites="smith_disorder-free_2017"><sup><a
href="#ref-smith_disorder-free_2017"
role="doc-biblioref"><strong>smith_disorder-free_2017?</strong></a></sup></span>}
the FK model has extensively many conserved quantities which can act as
an effective disorder potential for the electronic sector.</p>
<p>Due to Pauli exclusion, the maximum filling occurs when one of each
species occupies each lattice site such that <span
class="math inline">\(\sum_i (n_{i\downarrow} + n_{i\uparrow} )/ N =
2\)</span>. Here we focus on the half filled case which also displays
particle-hole symmetry (see later).</p>
<h2 id="falikov-kimball-and-hubbard-models">Falikov Kimball and Hubbard
models</h2>
<p>We will first introduce the standard Hubbard and Falikov-Kimball (FK)
models then look at some of their properties. We’ll then cover why the
Falikov-Kimball model represents an interesting system in which to study
disorder free localisation.</p>
<h3 id="hubbard-model">Hubbard model</h3>
<p>The Hubbard model gives a very simple setting in which to study
interacting, itinerant electrons. It is a tight binding model of spin
half electrons with finite bandwidth <span
class="math inline">\(t\)</span> and a repulsive on-site interaction
<span class="math inline">\(U &gt; 0\)</span>.</p>
<p><span class="math display">\[
H = -\sum_{&lt;ij&gt;,\sigma} t_{\sigma}
c^\dagger_{i\sigma}c_{j\sigma} + U \sum_{i} (n_{i \uparrow} - 1/2)(
n_{i\downarrow} - 1/2) - \mu \sum_i \left( n_{i \uparrow} + n_{i
\downarrow} \right)
\]</span></p>
<p>in standard notation. The standard Hubbard model corresponds to the
case <span class="math inline">\(t_{\uparrow} = t_{\downarrow}\)</span>.
Here we have used the particle-hole symmetric version of the interaction
term, which is more often given as <span class="math inline">\(n_{i
\uparrow} n_{i\downarrow}\)</span>. The difference just amounts to a
redefinition of the chemical potential.</p>
<p>Hubbard originally used the model at half filling <span
class="math inline">\(\mu = 0\)</span> to explain the Mott
metal-insulator (MI) transition, however it has seen applications to
high-temperature superconductivity and become target for cold-atom
optical trap experiments.<span class="citation"
data-cites="noauthor_hubbard_2013"><sup><a
href="#ref-noauthor_hubbard_2013"
role="doc-biblioref"><strong>noauthor_hubbard_2013?</strong></a></sup></span>,
greiner_quantum_2002, jordens_mott_2008}. While simple, only a few
analytic results exist, namely the Bethe ansatz<span class="citation"
data-cites="lieb_absence_1968"><sup><a href="#ref-lieb_absence_1968"
role="doc-biblioref"><strong>lieb_absence_1968?</strong></a></sup></span>}
which proves the absence of even a zero temperature phase transition in
the 1D model and Nagaoka’s theorem<span class="citation"
data-cites="nagaoka_ferromagnetism_1966"><sup><a
href="#ref-nagaoka_ferromagnetism_1966"
role="doc-biblioref"><strong>nagaoka_ferromagnetism_1966?</strong></a></sup></span>}
which proves that the three dimensional model has a ferromagnetic ground
state in the vicinity of half filling.</p>
<h3 id="falikov-kimball-model">Falikov-Kimball model</h3>
<p>The Falikov-Kimball model corresponds to the case <span
class="math inline">\(t_{\downarrow} = 0\)</span>. It can be interpreted
as two coupled spinless electron bands with infinite mass ratio. An
itinerant light species with creation operator <span
class="math inline">\(c^\dagger_{i\uparrow}\)</span> coupled to an
infinitely heavy, immobile species with density operator <span
class="math inline">\(n_{i\downarrow}\)</span>. These are often called c
and f electrons or electrons and ions. The model was first introduced by
Hubbard in 1963 as a model of interacting localised and de-localised
electron bands and gained its name from Falikov and Kimball’s use of it
to study the MI transition in rare-earth materials<span class="citation"
data-cites="hubbard_j._electron_1963"><sup><a
href="#ref-hubbard_j._electron_1963"
role="doc-biblioref"><strong>hubbard_j._electron_1963?</strong></a></sup></span>,
falicov_simple_1969}.</p>
<p>Here we will use refer to the light spinless species as
<code>electrons' with creation operator $c^\dagger_{i}$ and the heavy species as</code>ions’
with density operator <span class="math inline">\(n_i\)</span>. When the
the density operator of the electrons is needed I’ll always use <span
class="math inline">\(c^\dagger_{i}c_{i}\)</span>. We also set <span
class="math inline">\(t = 1\)</span>.</p>
<p><span class="math display">\[
H_{\mathrm{FK}} = -\sum_{&lt;ij&gt;} c^\dagger_{i}c_{j} + U \sum_{i}
(c^\dagger_{i}c_{i} - 1/2)( n_i - 1/2) - \mu \sum_i
\left(c^\dagger_{i}c_{i} + n_{i}\right)
\]</span> % ### Particle-Hole Symmetry The Hubbard and FK models on a
bipartite lattice have particle-hole (PH) symmetry <span
class="math inline">\(P^\dagger H P = - H\)</span>, accordingly they
have symmetric energy spectra. The associated symmetry operator <span
class="math inline">\(P\)</span> exchanges creation and annihilation
operators along with a sign change between the two sublattices.</p>
<p><span class="math display">\[ d^\dagger_{i\sigma} = (-1)^i
c_{i\sigma}\]</span> % The entirely filled state <span
class="math inline">\(\ket{\Omega} = \sum_{j\rho} c^\dagger_{j\rho}
\ket{0}\)</span> becomes the new vacuum state: <span
class="math display">\[d_{i\sigma} \ket{\Omega} = (-1)^i
c^\dagger_{i\sigma} \sum_{j\rho} c^\dagger_{j\rho} \ket{0} = 0\]</span>
% The number operator <span class="math inline">\(n_{i\sigma} =
0,1\)</span> counts holes rather than particles: <span
class="math display">\[ d^\dagger_{i\sigma} d_{i \sigma} = c_{i\sigma}
c^\dagger_{i\sigma} = 1 - c^\dagger_{i\sigma} c_{i\sigma}\]</span> %
With the last equality following from the fermionic commutation
relations. In the case of nearest neighbour hopping on a bipartite
lattice this transformation also leaves the hopping term unchanged:
<span class="math display">\[ d^\dagger_{i\sigma} d_{j \sigma} =
(-1)^{i+j} c_{i\sigma} c^\dagger_{j\sigma} = c^\dagger_{i\sigma}
c_{j\sigma} \]</span> % Since when <span
class="math inline">\(i\)</span> and <span
class="math inline">\(j\)</span> label sites on separate sublattices,
<span class="math inline">\((-1)^{i-j} = -1\)</span> and this is
absorbed into rearranging the operators via their anti-commutator.</p>
<p>Defining the particle density <span
class="math inline">\(\rho\)</span> as the number of fermions per site:
<span class="math display">\[
\rho = \frac{1}{N} \sum_i \left( n_{i \uparrow} + n_{i \downarrow}
\right)
\]</span> % The PH symmetry maps the Hamiltonian to itself with the sign
of the chemical potential reversed and the density is inverted about
half filling: <span class="math display">\[ \text{PH} : H(t, U, \mu)
\rightarrow H(t, U, -\mu) \]</span> <span class="math display">\[ \rho
\rightarrow 2 - \rho \]</span> % The Hamiltonian is symmetric under PH
at <span class="math inline">\(\mu = 0\)</span> and so must all the
observables, hence half filling <span class="math inline">\(\rho =
1\)</span> occurs here. This symmetry and known observable acts as a
useful test for the numerical calculations.</p>
<h3 id="thermodynamics-of-the-fk-model">Thermodynamics of the FK
model</h3>
\begin{figure}
<p>} \end{figure}</p>
<p>At half filling and in dimensions greater than one, the FK model
exhibits a phase transition at some <span
class="math inline">\(U\)</span> dependent critical temperature <span
class="math inline">\(T_c(U)\)</span> to a low temperature charge
density wave state in which the ions occupy one of the two sublattices A
and B<span class="citation"
data-cites="maska_thermodynamics_2006-1"><sup><a
href="#ref-maska_thermodynamics_2006-1"
role="doc-biblioref"><strong>maska_thermodynamics_2006-1?</strong></a></sup></span>}.
The order parameter is the square of the staggered magnetisation: <span
class="math display">\[
M = \sum_{i \in A} n_i - \sum_{i \in B} n_i
\]</span> % In the disordered phase Ref.<span class="citation"
data-cites="antipov_interaction-tuned_2016-1"><sup><a
href="#ref-antipov_interaction-tuned_2016-1"
role="doc-biblioref"><strong>antipov_interaction-tuned_2016-1?</strong></a></sup></span>}
identifies an interplay between Anderson localisation at weak
interaction and a Mott insulator phase in the strongly interacting
regime.</p>
<p>In the one dimensional FK model, however, Peierls’ argument<span
class="citation" data-cites="peierls_isings_1936"><sup><a
href="#ref-peierls_isings_1936"
role="doc-biblioref"><strong>peierls_isings_1936?</strong></a></sup></span>,
kennedy_itinerant_1986} and the Bethe ansatz<span class="citation"
data-cites="lieb_absence_1968"><sup><a href="#ref-lieb_absence_1968"
role="doc-biblioref"><strong>lieb_absence_1968?</strong></a></sup></span>}
make it clear that there is no ordered CDW phase. Peierls’ argument is
that one should consider the difference in free energy <span
class="math inline">\(\Delta F = \Delta E - T\Delta S\)</span> between
an ordered state and a state with single domain wall in the order
parameter. In the Ising model this would be having the spins pointing up
in one part of the model and down in the other, for a CDW phase it means
having the ions occupy the A sublattice in one part and the B sublattice
in the other.</p>
<p>Short range interactions will produce a constant energy penalty for
such a domain wall that does not scale with system size while in 1D
there are <span class="math inline">\(L\)</span> such states so the
domain wall is associated with entropy <span class="math inline">\(S
\propto \ln L\)</span> which dominates in the thermodynamic limit. The
slow logarithmic scaling suggests we should be wary of finite size
scaling effects.</p>
<p>One dimensional systems are more amenable to numerical and
experimental study so we add long range staggered interactions to bring
back the ordered phase:</p>
<p><span class="math display">\[ H_{\textrm{int}} = 4J \sum_{ij}
\frac{(-1)^{|i-j|}}{ |i - j|^{\alpha} } (n_i - 1/2) (n_j - 1/2) = J
\sum_{ij} |i - j|^{-\alpha} \tau_i \tau_j\]</span> % at half-filling the
modified Hamiltonian is then: <span class="math display">\[
H_{\mathrm{FK}}^* &amp;= -\sum_{&lt;ij&gt;} c^\dagger_{i}c_{j} + U
\sum_{i} (c^\dagger_{i}c_{i} - 1/2)( n_i - 1/2) \\
&amp;+ 4J \sum_{ij} \frac{(-1)^{|i-j|}}{ |i - j|^{\alpha} } (n_i -
1/2) (n_j - 1/2) \\
&amp;= -\sum_{&lt;ij&gt;} c^\dagger_{i}c_{j} + 2U \sum_{i} (-1)^i
(c^\dagger_{i}c_{i} - 1/2)\tau_i + J \sum_{ij} |i - j|^{-\alpha} \tau_i
\tau_j \\
\]</span> % The form of this interaction comes from interpreting the
f-electrons as a classical Ising chain using a staggered mapping <span
class="math inline">\(\tau_i = (-1)^i (2n_i^ f - 1)\)</span> so that
ferromagnetic order in the <span class="math inline">\(\tau_i\)</span>
variables corresponds to a CDW state in the <span
class="math inline">\(n_i^f\)</span> variables. It also preserves the
particle hole symmetry because for the ions the PH transformation
corresponds to <span class="math inline">\(n_i \rightarrow 1 -
n_i\)</span>. When <span class="math inline">\(U = 0\)</span> the model
decouples into a long ranged Ising model and free fermions.</p>
<p>Our extension to the FK model could now be though of as spinless
fermions coupled to a long range Ising (LRI) model. The LRI model has
been extensively studied and its behaviour may be bear relation to the
behaviour of our modified FK model.</p>
<p><span class="math display">\[H_{\mathrm{LRI}} = \sum_{ij}
J(\abs{i-j}) \tau_i \tau_j = J \sum_{i\neq j} |i - j|^{-\alpha} \tau_i
\tau_j\]</span> % Rigorous renormalisation group arguments show that the
LRI model has an ordered phase in 1D for $1 &lt; &lt; 2 $<span
class="citation" data-cites="dyson_existence_1969"><sup><a
href="#ref-dyson_existence_1969"
role="doc-biblioref"><strong>dyson_existence_1969?</strong></a></sup></span>}.
Peierls’ argument can be extended<span class="citation"
data-cites="thouless_long-range_1969"><sup><a
href="#ref-thouless_long-range_1969"
role="doc-biblioref"><strong>thouless_long-range_1969?</strong></a></sup></span>}
to provide intuition for why this is the case. Again considering the
energy difference between the ordered state <span
class="math inline">\(\ket{\ldots\uparrow\uparrow\uparrow\uparrow\ldots}\)</span>
and a domain wall state <span
class="math inline">\(\ket{\ldots\uparrow\uparrow\downarrow\downarrow\ldots}\)</span>.
In the case of the LRI model, careful counting shows that this energy
penalty is: <span class="math display">\[\Delta E \propto
\sum_{n=1}^{\infty} n J(n)\]</span> % because each interaction between
spins separated across the domain by a bond length <span
class="math inline">\(n\)</span> can be drawn between <span
class="math inline">\(n\)</span> equivalent pairs of sites. Ruelle
proved rigorously for a very general class of 1D systems, that if <span
class="math inline">\(\Delta E\)</span> or its many-body generalisation
converges in the thermodynamic limit then the free energy is
analytic<span class="citation"
data-cites="ruelle_statistical_1968"><sup><a
href="#ref-ruelle_statistical_1968"
role="doc-biblioref"><strong>ruelle_statistical_1968?</strong></a></sup></span>}.
This rules out a finite order phase transition, though not one of the
Kosterlitz-Thouless type. Dyson also proves this though with a slightly
different condition on <span class="math inline">\(J(n)\)</span><span
class="citation" data-cites="dyson_existence_1969"><sup><a
href="#ref-dyson_existence_1969"
role="doc-biblioref"><strong>dyson_existence_1969?</strong></a></sup></span>}.</p>
<p>With a power law form for <span class="math inline">\(J(n)\)</span>,
there are three cases to consider:</p>
<ol type="1">
<li>$ = 0$ For infinite range interactions the Ising model is exactly
solveable and mean field theory is exact<span class="citation"
data-cites="lipkin_validity_1965"><sup><a
href="#ref-lipkin_validity_1965"
role="doc-biblioref"><strong>lipkin_validity_1965?</strong></a></sup></span>}.</li>
<li>$ $ For slowly decaying interactions <span
class="math inline">\(\sum_n J(n)\)</span> does not converge so the
Hamiltonian is non-extensive, a case which won’t be further considered
here.</li>
<li>$ 1 &lt; &lt; 2 $ A phase transition to an ordered state at a finite
temperature.</li>
<li>$ = 2 $ The energy of domain walls diverges logarithmically, and
this turns out to be a Kostelitz-Thouless transition<span
class="citation" data-cites="thouless_long-range_1969"><sup><a
href="#ref-thouless_long-range_1969"
role="doc-biblioref"><strong>thouless_long-range_1969?</strong></a></sup></span>}.</li>
<li>$ 2 &lt; $ For quickly decaying interactions, domain walls have a
finite energy penalty, hence Peirels’ argument holds and there is no
phase transition.</li>
</ol>
<h3 id="thermodynamics">Thermodynamics</h3>
<p>On bipartite lattices in dimensions 2 and above the FK model exhibits
a finite temperature phase transition to an ordered charge density wave
(CDW) phase<span class="citation"
data-cites="maska_thermodynamics_2006-1"><sup><a
href="#ref-maska_thermodynamics_2006-1"
role="doc-biblioref"><strong>maska_thermodynamics_2006-1?</strong></a></sup></span>.
In this phase, the ions are confined to one of the two sublattices,
breaking the <span class="math inline">\(\mathbb{Z}_2\)</span>
symmetry.</p>
<p>In 1D, however, Periel’s argument<span class="citation"
data-cites="peierls_isings_1936 kennedy_itinerant_1986"><sup><a
href="#ref-peierls_isings_1936"
role="doc-biblioref"><strong>peierls_isings_1936?</strong></a>,<a
href="#ref-kennedy_itinerant_1986"
role="doc-biblioref"><strong>kennedy_itinerant_1986?</strong></a></sup></span>
states that domain walls only introduce a constant energy penalty into
the free energy while bringing a entropic contribution logarithmic in
system size. Hence the 1D model does not have a finite temperature phase
transition. However 1D systems are much easier to study numerically and
admit simpler realisations experimentally. We therefore introduce a long
range coupling between the ions in order to stabilise a CDW phase in 1D.
This leads to a disordered system that is gaped by the CDW background
but with correlated fluctuations leading to a disorder-free correlation
induced mobility edge in one dimension.</p>
<h3 id="markov-chain-monte-carlo">Markov Chain Monte Carlo</h3>
<p>To evaluate thermodynamic averages we perform a classical Markov
Chain Monte Carlo random walk over the space of ionic configurations, at
each step diagonalising the effective electronic Hamiltonian<span
class="citation" data-cites="maska_thermodynamics_2006-1"><sup><a
href="#ref-maska_thermodynamics_2006-1"
role="doc-biblioref"><strong>maska_thermodynamics_2006-1?</strong></a></sup></span>}.
Using a binder-cumulant method<span class="citation"
data-cites="binder_finite_1981 musial_monte_2002"><sup><a
href="#ref-binder_finite_1981"
role="doc-biblioref"><strong>binder_finite_1981?</strong></a>,<a
href="#ref-musial_monte_2002"
role="doc-biblioref"><strong>musial_monte_2002?</strong></a></sup></span>,
we demonstrate the model has a finite temperature phase transition when
the interaction is sufficiently long ranged. We then estimate the
density of states and the inverse participation ratio as a function of
energy to diagnose localisation properties. We show preliminary results
that the in-gap states induced at finite temperature are localised while
the states in the unperturbed bands remain extended, evidence for a
mobility edge.</p>
<div id="fig:binder" class="fignos">
<figure>
<img src="/assets/thesis/figure_code/fk_chapter/binder.png"
style="width:100.0%" alt="Figure 1: Hello I am the figure caption!" />
<figcaption aria-hidden="true"><span>Figure 1:</span> Hello I am the
figure caption!</figcaption>
</figure>
</div>
<p>Macro definitions in this cell <span class="math display">\[
\newcommand{\expval}[1]{\langle #1 \rangle}
\newcommand{\ket}[1]{|#1\rangle}
\newcommand{\bra}[1]{\langle#1|}
\newcommand{\op}[2]{|#1\rangle \langle#2|}
\]</span></p>
<p><span class="math display">\[
\expval{O}, \op{\alpha}{\beta}, \ket{\psi}
\]</span></p>
<h2 id="localisation-1">Localisation</h2>
<h3 id="thermalisation">Thermalisation</h3>
<p>Isolated classical systems generally thermalise if they are large
enough. Since classical dynamics is the limit of some underlying quantum
dynamics, it seems reasonable to suggest that isolated quantum systems
also thermalise in some related sense. However it is not immediately
obvious what thermalisation should mean in a quantum setting where the
presence of unitary time dynamics implies full information about a
system’s initial state is always preserved.</p>
<p>A potential solution lies in the eigenstate thermalisation
hypothesis. It states that if a system thermalises, then we can isolate
small patches of a larger system, trace everyhing else out, and get a
thermal density matrix.</p>
<p>Following Ref.<span class="citation"
data-cites="abanin_recent_2017"><sup><a href="#ref-abanin_recent_2017"
role="doc-biblioref"><strong>abanin_recent_2017?</strong></a></sup></span>,
consider the time evolution of a local operator <span
class="math inline">\(\hat{O}\)</span> <span class="math display">\[
\expval{\hat{O}}{\psi(t)} = \sum_{\alpha \beta} C^*_\alpha C_\beta
e^{i(E_\alpha - E_\beta)} O_{\alpha \beta}\]</span></p>
<p>Where <span class="math inline">\(C_\alpha\)</span> are determined by
the initial state and <span class="math inline">\(O_{\alpha \beta} =
\expval{\alpha | \hat{O} | \beta}\)</span> are the matrix elements of
<span class="math inline">\(\hat{O}\)</span> with respect to the energy
eigenstates. Srednicki<span class="citation"
data-cites="srednicki_chaos_1994"><sup><a
href="#ref-srednicki_chaos_1994"
role="doc-biblioref"><strong>srednicki_chaos_1994?</strong></a></sup></span>}
introduced the ansatz that for local operators:</p>
<p><span class="math display">\[O_{\alpha \beta} =
O(E)\delta_{\alpha\beta} + e^{-S(E)/2} f(E,\omega)
R_{\alpha\beta}\]</span></p>
<p>with <span class="math inline">\(E = (E_\alpha + E_\beta)\)</span>,
<span class="math inline">\(\omega = (E_\alpha - E_\beta)\)</span> and
<span class="math inline">\(R_{\alpha\beta}\)</span> are sampled from
some distribution with zero mean and unit variance. The first term
asserts that the diagonal elements are given by the thermal expectation
value <span class="math inline">\(O(E) = Tr[e^{-\beta \hat{H}}
\hat{O}]/\mathcal{Z}\)</span> with <span
class="math inline">\(\beta\)</span> an effective temperature defined by
equating the energy to the expectation of the Hamiltonian at that
temperature <span class="math inline">\(E = Tr[H e^{-\beta
\hat{H}}/\mathcal{Z}]\)</span>.</p>
<p>The second term deals with thermodynamic fluctuations scaled by the
entropy <span class="math inline">\(S(E) = -Tr(\rho \log \rho)\)</span>
where <span class="math inline">\(\rho = e^{-\beta \hat{H}}\)</span> and
<span class="math inline">\(\mathcal{Z} = Tr[e^{-\beta
\hat{H}}]\)</span>.</p>
<p>With this ansatz the long time average of the observable becomes
equal to the thermal expectations with fluctuations suppressed by the
entropic term <span class="math inline">\(e^{-S(E)}\)</span> and the
rapidly varying phase factors <span class="math inline">\(e^{i(E_\alpha
- E_\beta)}\)</span>. This statement of the ETH has verified for the
quantum hard sphere model<span class="citation"
data-cites="srednicki_chaos_1994"><sup><a
href="#ref-srednicki_chaos_1994"
role="doc-biblioref"><strong>srednicki_chaos_1994?</strong></a></sup></span>
and numerically for other models<span class="citation"
data-cites="khatami_fluctuation-dissipation_2013 dalessio_quantum_2016"><sup><a
href="#ref-khatami_fluctuation-dissipation_2013"
role="doc-biblioref"><strong>khatami_fluctuation-dissipation_2013?</strong></a>,<a
href="#ref-dalessio_quantum_2016"
role="doc-biblioref"><strong>dalessio_quantum_2016?</strong></a></sup></span>.</p>
<p>An alternate view on ETH is the statement that in thermalising
systems individual eigenstates look thermal when viewed locally. Take a
eigenstate <span class="math inline">\(|\alpha\rangle\)</span> with
energy <span class="math inline">\(E_\alpha\)</span> and as before
define an effective temperature with <span
class="math inline">\(E_\alpha = Tr[H e^{-\beta
\hat{H}}/\mathcal{Z}]\)</span>. This statement of the ETH says that if
we partition the system into subsystems A and B with a limit taken as B
becomes very large, B will act as a heat bath for A. Specifically the
reduced density matrix <span class="math inline">\(\rho_A = Tr_B
\op{\alpha}{\alpha}\)</span> is equal to the thermal density matrix:</p>
<p><span class="math display">\[\rho_A = Tr_B |\alpha\rangle \langle
\alpha| = \mathcal{Z}^{-1} Tr_B [e^{-\beta \hat{H}}] \]</span></p>
<p>Intuitively, for thermalisation to happen, the degrees of freedom
must be sufficiently well coupled that energy transport occurs. This
condition is broken by systems with localised states so a lack of
thermalisation is often used as a diagnostic tool for localisation.</p>
<h3 id="anderson-localisation">Anderson Localisation</h3>
<p>Localisation was first studied by Anderson in 1958<span
class="citation" data-cites="anderson_absence_1958-1"><sup><a
href="#ref-anderson_absence_1958-1"
role="doc-biblioref"><strong>anderson_absence_1958-1?</strong></a></sup></span>
in the context of non-interacting fermions subject to a static or
quenched disorder potential <span class="math inline">\(V_j\)</span>
drawn uniformly from the interval <span
class="math inline">\([-W,W]\)</span>:</p>
<p><span class="math display">\[
H = -t\sum_{\expval{jk}} c^\dagger_j c_k + \sum_j V_j c_j^\dagger c_j
\]</span></p>
<p>At sufficiently strong disorder the Anderson model exhibits
exponentially localised eigenfunctions <span
class="math inline">\(\psi(x) = f(x) e^{-x/\lambda}\)</span> which
cannot contribute to diffusive transport processes. Except in 1D where
any disorder strength is sufficient. Intuitively this happens because
hopping processes between nearby sites become off-resonant, hindering
the hybridisation that would normally lead to extended Bloch states<span
class="citation" data-cites="kramer_localization_1993"><sup><a
href="#ref-kramer_localization_1993"
role="doc-biblioref"><strong>kramer_localization_1993?</strong></a></sup></span>.</p>
<p>In one and two dimensions, all the states in the Anderson model are
localised. In three dimensions there are mobility edges. Mobility edges
are critical energies in the spectrum which separate delocalised states
in a band from localised states which form a band tail<span
class="citation" data-cites="abanin_recent_2017"><sup><a
href="#ref-abanin_recent_2017"
role="doc-biblioref"><strong>abanin_recent_2017?</strong></a></sup></span>}.
An argument due to Lifshitz shows that the density of state of the band
tail should decay exponentially and localised and extended stats cannot
co-exist at the same energy as they would hybridise into extended
states<span class="citation"
data-cites="kramer_localization_1993"><sup><a
href="#ref-kramer_localization_1993"
role="doc-biblioref"><strong>kramer_localization_1993?</strong></a></sup></span>}.</p>
<p>It was thought that mobility edges could not exist in 1D because all
the states localised in the presence of any amount of disorder. This is
true for uncorrelated potentials<span class="citation"
data-cites="goldshtein_pure_1977"><sup><a
href="#ref-goldshtein_pure_1977"
role="doc-biblioref"><strong>goldshtein_pure_1977?</strong></a></sup></span>}.
However, it was shown that if the disorder potential <span
class="math inline">\(V_j\)</span> contains spatial correlations
mobility edges do exist in 1D<span class="citation"
data-cites="izrailev_localization_1999"><sup><a
href="#ref-izrailev_localization_1999"
role="doc-biblioref"><strong>izrailev_localization_1999?</strong></a></sup></span>,
izrailev_anomalous_2012}. Ref.<span class="citation"
data-cites="croy_anderson_2011"><sup><a href="#ref-croy_anderson_2011"
role="doc-biblioref"><strong>croy_anderson_2011?</strong></a></sup></span>}
extends this work to look at power law decay of the correlations: <span
class="math display">\[ C(l) = \expval{V_i V_{i+l}} \propto l^{-\alpha}
\]</span> % Figure <span
class="math inline">\(\ref{fig:anderson_dos}\)</span> shows numerical
calculations of the Localisation length (see later) and density of
states for the power law correlated Anderson model. At the unperturbed
band edges <span class="math inline">\(\abs{E} = 2\)</span>, the states
transition from extended to localised. The behaviour close to the edge
takes a universal scaling form with exponents dependant on <span
class="math inline">\(\alpha\)</span>.</p>
<h3 id="many-body-localisation">Many Body Localisation</h3>
<p>A related phenomena known as many body localisation (MBL) was found
in interacting systems with quenched disorder. A simple example comes
from adding density-density interactions to the Anderson model:</p>
<p><span class="math display">\[
H = -t\sum_{\expval{jk}} c^\dagger_j c_k + \sum_j V_j c_j^\dagger c_j +
U\sum_{jk} n_j n_k
\]</span> % with <span class="math inline">\(n_j = c^\dagger_j
c_j\)</span> Here, in contrast to the Anderson model, localisation
phenomena can be proven robust to weak perturbations of the
Hamiltonian<span class="citation"
data-cites="imbrie_many-body_2016"><sup><a
href="#ref-imbrie_many-body_2016"
role="doc-biblioref"><strong>imbrie_many-body_2016?</strong></a></sup></span>}.</p>
<p>MBL is defined by the emergence of an extensive number of quasi-local
operators called local integrals of motions (LIOMs) or l-bits. Following
Ref.<span class="citation" data-cites="abanin_recent_2017"><sup><a
href="#ref-abanin_recent_2017"
role="doc-biblioref"><strong>abanin_recent_2017?</strong></a></sup></span>},
using a spin system with variables <span
class="math inline">\(\sigma^z_i\)</span>, any operator can be written
in the general form:</p>
<p><span class="math display">\[ \tau^z_i = \sigma^z_i +
\sum_{\alpha\beta kl} f_{kl}^{\alpha\beta} \sigma^\alpha_{i+k}
\sigma_z\beta_{i+k} + ...\]</span> % what defines a MBL system is that
there exist extensively many <span
class="math inline">\(\tau^z_i\)</span> for which the coefficients decay
exponentially with distance <span
class="math inline">\(f_{kl}^{\alpha\beta} \propto
e^{-\max(\abs{l},\abs{k}) / \xi}\)</span>. These LIOMs commute with the
Hamiltonian and each other <span class="math inline">\([\hat{H},
\tau^z_i] = [\tau^z_i, \tau^z_j] = 0\)</span>. It is this extensive
number of conserved local charges that leads to the localisation
properties of MBL. It also has implications for the way entanglement
grows over time in MBL systems.</p>
<p>Since the Hamiltonian commutes with all the LIOMs and they are a
complete operator basis, the Hamiltonian can be written as:</p>
<p><span class="math display">\[\hat{H} = \sum_{i} h_i \tau^z_i +
\sum_{ij} J_{ij} \tau^z_i \tau^z_j + \sum_{ijk} J_{ij} \tau^z_i \tau^z_j
\tau^z_k+ ...\]</span> % Where again the couplings decay exponentially,
albeit with a different length scale <span
class="math inline">\(\Bar{\xi}\)</span>. From this form we see that
distant l-bits can only become entangled on a timescale of:</p>
<p><span class="math display">\[ t_{\mathrm{ent}}(r) \propto
\frac{\hbar}{J_0} e^{r/\Bar{\xi}} \]</span> % and hence quantum
correlations and entanglement propagates logarithmically in MBL
systems<span class="citation"
data-cites="imbrie_diagonalization_2016"><sup><a
href="#ref-imbrie_diagonalization_2016"
role="doc-biblioref"><strong>imbrie_diagonalization_2016?</strong></a></sup></span>}.</p>
<h3 id="disorder-free-localisation">Disorder Free localisation</h3>
<p>Both Anderson localisation and MBL depend on the presence of quenched
disorder. Recently the idea of disorder-free localisation has gained
traction, asking whether the disorder necessary to generate localisation
can be generated entirely from the dynamics of a model itself.</p>
<p>The idea was first proposed in the context of Helium mixtures<span
class="citation" data-cites="kagan1984localization"><sup><a
href="#ref-kagan1984localization"
role="doc-biblioref">1</a></sup></span>} and then extended to
heavy-light mixtures in which multiple species with large mass ratios
interact, the idea being that the heavier particles act as an effective
disorder potential for the lighter ones, inducing localisation. Two such
models<span class="citation"
data-cites="yao_quasi-many-body_2016"><sup><a
href="#ref-yao_quasi-many-body_2016"
role="doc-biblioref"><strong>yao_quasi-many-body_2016?</strong></a></sup></span>,schiulaz_dynamics_2015}
instead find that the models thermalise exponentially slowly in system
size, which Ref.<span class="citation"
data-cites="yao_quasi-many-body_2016"><sup><a
href="#ref-yao_quasi-many-body_2016"
role="doc-biblioref"><strong>yao_quasi-many-body_2016?</strong></a></sup></span>}
dubs Quasi-MBL. A. Smith, J. Knolle et al instead looked at models
containing an extensive number of conserved quantities and demonstrated
true disorder free localisation<span class="citation"
data-cites="smith_disorder-free_2017"><sup><a
href="#ref-smith_disorder-free_2017"
role="doc-biblioref"><strong>smith_disorder-free_2017?</strong></a></sup></span>}.</p>
<h3 id="diagnostics-of-localisation">Diagnostics of Localisation</h3>
<h4 id="inverse-participation-ratio">Inverse Participation Ratio</h4>
<p>The inverse participation ratio is defined for a normalised wave
function <span class="math inline">\(\psi_i = \psi(x_i), \sum_i
\abs{\psi_i}^2 = 1\)</span> as its fourth moment<span class="citation"
data-cites="kramer_localization_1993"><sup><a
href="#ref-kramer_localization_1993"
role="doc-biblioref"><strong>kramer_localization_1993?</strong></a></sup></span>}:
<span class="math display">\[
P^{-1} = \sum_i \abs{\psi_i}^4
\]</span> % It acts as a measure of the portion of space occupied by the
wave function. For localised states it will be independent of system
size while for plane wave states in d dimensions $P = L^d $. States may
also be intermediate between localised and extended, described by their
fractal dimensionality <span class="math inline">\(d &gt; d* &gt;
0\)</span>: <span class="math display">\[
P(L) \goeslike L^{d*}
\]</span> % For extended states <span class="math inline">\(d* =
0\)</span> while for localised ones <span class="math inline">\(d* =
0\)</span>. In this work we take use an energy resolved IPR<span
class="citation" data-cites="antipov_interaction-tuned_2016-1"><sup><a
href="#ref-antipov_interaction-tuned_2016-1"
role="doc-biblioref"><strong>antipov_interaction-tuned_2016-1?</strong></a></sup></span>:
<span class="math display">\[
DOS(\omega) = \sum_n \delta(\omega - \epsilon_n)
IPR(\omega) = DOS(\omega)^{-1} \sum_{n,i} \delta(\omega - \epsilon_n)
\abs{\psi_{n,i}}^4
\]</span> Where <span class="math inline">\(\psi_{n,i}\)</span> is the
wavefunction corresponding to the energy <span
class="math inline">\(\epsilon_n\)</span> at the ith site. In practice
we bin the energies and IPRs into a fine energy grid and use Lorentzian
smoothing if necessary.</p>
<h4 id="transfer-matrix-approach">Transfer Matrix Approach</h4>
<p>The transfer matrix method (TMM) can be used to calculate the
localisation length of the eigenstates of a system. Following Refs.<span
class="citation" data-cites="kramer_localization_1993"><sup><a
href="#ref-kramer_localization_1993"
role="doc-biblioref"><strong>kramer_localization_1993?</strong></a></sup></span>,
smith_dynamical_2018}, for bi-linear, 1D Hamiltonians the method
represents the action of <span class="math inline">\(H\)</span> on a
state <span class="math inline">\(\ket{\psi} = \sum_i \psi_i
\ket{i}\)</span> with energy E by a matrix equation: <span
class="math display">\[
H &amp;= - \sum_i (c^\dagger_i c_{i+1} + c^\dagger_{i+1} c_{i}) - \sum_i
h_i c^\dagger_i c_i \\
E\ket{\psi} &amp;= H \ket{\psi} \\
\label{eq:tmm_difference} E \psi_i &amp;= -(\psi_{i+1} + \psi_{i-1}) -
h_i \psi_i
\]</span> % Where Eq. <span
class="math inline">\(\ref{eq:tmm_difference}\)</span> can be
represented by a matrix equation: <span class="math display">\[
\begin{pmatrix}
\psi_{i+1}\\
\psi_{i}
\end{pmatrix}
=
\begin{pmatrix}
-(E + h_i) &amp; -1\\
1 &amp; 0
\end{pmatrix}
\begin{pmatrix}
\psi_{i}\\
\psi_{i-1}
\end{pmatrix}
= T_i
\begin{pmatrix}
\psi_{i}\\
\psi_{i-1}
\end{pmatrix}
\]</span> % Defining a product along the chain: <span
class="math display">\[Q_L = \prod_{i=0}^L T_i\]</span> % Oseledec’s
theorem proves that there exists a limiting matrix <span
class="math inline">\(\Gamma\)</span>: <span class="math display">\[
\Gamma = \lim_{L \to \infty} (Q_L Q_L^\dagger)^{\frac{1}{2L}}
\]</span> % with eigenvalues <span
class="math inline">\(\exp(\gamma_i)\)</span> where <span
class="math inline">\(\gamma_i\)</span> are the Lyapunov exponents of
<span class="math inline">\(Q_L\)</span>. The smallest Lyapunov exponent
is the inverse of the localisation length of the state. In practice one
takes <span class="math inline">\(Q_L\)</span> with L equal to the
system size, finds the smallest eigenvalue q and estimates the
localisation length by: <span class="math display">\[
\lambda = \frac{L}{\ln{q}}
\]</span> % As noted by<span class="citation"
data-cites="smith_dynamical_2018"><sup><a
href="#ref-smith_dynamical_2018"
role="doc-biblioref"><strong>smith_dynamical_2018?</strong></a></sup></span>
this method can be numerically unstable because the matrix elements
diverge or vanish exponentially. To get around this, the authors break
the matrix multiplication into chunks and the logarithms of the
eigenvalues of each are stored separately.</p>
<h2 id="numerical-methods">Numerical Methods</h2>
<p>In this section we will define the Markov Chain Monte Carlo (MCMC)
method in general then detail its application to the FK model. We will
then cover methods applicable to the measurements obtained from MCMC.
These include calculation of the density of states and energy resolved
inverse participation ratio as well as phase transition diagnostics such
as the Binder cumulant.</p>
<h3 id="markov-chain-monte-carlo-1">Markov Chain Monte Carlo}</h3>
<p>Markov Chain Monte Carlo (MCMC) is a technique for evaluating thermal
expectation values <span class="math inline">\(\expval{O}\)</span> with
respect to some physical system defined by a set of states <span
class="math inline">\(\{x: x \in S\}\)</span> and a free energy <span
class="math inline">\(F(x)\)</span><span class="citation"
data-cites="krauth_introduction_1996"><sup><a
href="#ref-krauth_introduction_1996"
role="doc-biblioref"><strong>krauth_introduction_1996?</strong></a></sup></span>}.
The thermal expectation value is defined via a Boltzmann weighted sum
over the entire states: <span class="math display">\[
\tex{O} &amp;= \frac{1}{\Z} \sum_{x \in S} O(x) P(x) \\
P(x) &amp;= \frac{1}{\Z} e^{-\beta F(x)} \\
\Z &amp;= \sum_{x \in S} e^{-\beta F(x)}
\]</span></p>
<p>When the state space is too large to evaluate this sum directly, MCMC
defines a stochastic algorithm which generates a random walk <span
class="math inline">\(\{x_0\ldots x_i\ldots x_N\}\)</span> whose
distribution <span class="math inline">\(p(x_i)\)</span> approaches a
target distribution <span class="math inline">\(P(x)\)</span> in the
large N limit.</p>
<p><span class="math display">\[\lim_{i\to\infty} p(x_i) =
P(x)\]</span></p>
<p>In this case the target distribution will be the thermal one <span
class="math inline">\(P(x) \rightarrow \Z^{-1} e^{-\beta F(x)}\)</span>.
The major benefit of the method being that as long as one can express
the desired <span class="math inline">\(P(x)\)</span> up to a
multiplicative constant, MCMC can be applied. Since <span
class="math inline">\(e^{-\beta F(x)}\)</span> is relatively easy to
evaluate, MCMC provides a useful method for finite temperature
physics.</p>
<p>Once the random walk has been carried out for many steps, the
expectation values of <span class="math inline">\(O\)</span> can be
estimated from the MCMC samples: <span class="math display">\[
\tex{O} = \sum_{i = 0}^{N} O(x_i) + \mathcal{O}(\frac{1}{\sqrt{N}})
\]</span> The the samples in the random walk are correlated so the
samples effectively contain less information than <span
class="math inline">\(N\)</span> independent samples would. As a
consequence the variance is larger than the <span
class="math inline">\(\tex{O^2} - \tex{O}^2\)</span> form it would have
if the estimates were uncorrelated. Methods of estimating the true
variance of <span class="math inline">\(\tex{O}\)</span> and decided how
many steps are needed will be considered later.</p>
<p>Markov chains are defined by a transition function $(x_{i+1} x_i) $
giving the probability that a chain in state <span
class="math inline">\(x_i\)</span> at time <span
class="math inline">\(i\)</span> will transition to a state <span
class="math inline">\(x_{i+1}\)</span>. Since we must transition
somewhere at each step, this comes with the normalisation condition that
<span class="math inline">\(\sum\limits_x \T(x&#39; \rightarrow x) =
1\)</span>.</p>
<p>If we define an ensemble of Markov chains and consider the
distributions we get a sequence of distributions <span
class="math inline">\(\{p_0(x), p_1(x), p_2(x)\ldots\}\)</span> with
<span class="math display">\[p_{i+1}(x) &amp;= \sum_{x&#39; \in S}
p_i(x&#39;) \T(x&#39; \rightarrow x)\]</span> <span
class="math inline">\(p_o(x)\)</span> might be a delta function on one
particular starting state which would then be smoothed out by the
transition function repeatedly.</p>
<p>As we’d like to draw samples from the target distribution <span
class="math inline">\(P(x)\)</span> the trick is to choose $(x_{i+1}
x_i) $ such that :</p>
<p><span class="math display">\[
P(x) &amp;= \sum_{x&#39;} P(x&#39;) \T(x&#39; \rightarrow x)
\]</span> In other words the MCMC dynamics defined by <span
class="math inline">\(\T\)</span> must be constructed to have the target
distribution as their only fixed point. This condition is called the
global balance condition. Along with some more technical considerations
such as ergodcity which won’t be considered here, global balance
suffices to ensure that a MCMC method is correct.</p>
<p>A sufficient but not necessary condition for global balance to hold
is detailed balance:</p>
<p><span class="math display">\[
P(x) \T(x \rightarrow x&#39;) = P(x&#39;) \T(x&#39; \rightarrow x)
\]</span> % In practice most algorithms are constructed to satisfy
detailed balance though there are arguments that relaxing the condition
can lead to faster algorithms<span class="citation"
data-cites="kapfer_sampling_2013"><sup><a
href="#ref-kapfer_sampling_2013"
role="doc-biblioref"><strong>kapfer_sampling_2013?</strong></a></sup></span>}.</p>
<p>The goal of MCMC is then to choose <span
class="math inline">\(\T\)</span> so that it has the desired thermal
distribution <span class="math inline">\(P(x)\)</span> as its fixed
point and that it converges quickly onto it. This boils down to
requiring that the matrix representation of <span
class="math inline">\(T_{ij} = \T(x_i \to x_j)\)</span> has an
eigenvector equal to <span class="math inline">\(P_i = P(x_i)\)</span>
with eigenvalue 1 and all other eigenvalues with magnitude less than
one. The convergence time depends on the magnitude of the second largest
eigenvalue.</p>
<p>In order to actually choose new states according to <span
class="math inline">\(\T\)</span> one chooses states from a proposal
distribution <span class="math inline">\(q(x_i \to x&#39;)\)</span> that
can be directly sampled from. For instance, this might mean flipping a
single random spin in a spin chain, in which case <span
class="math inline">\(q(x&#39;\to x_i)\)</span> is the uniform
distribution on states reachable by one spin flip from <span
class="math inline">\(x_i\)</span>. The proposal <span
class="math inline">\(x&#39;\)</span> is then accepted or rejected with
an acceptance probability <span class="math inline">\(\A(x&#39;\to
x_{i+1})\)</span>, if the proposal is rejected then <span
class="math inline">\(x_{i+1} = x_{i}\)</span>. Now <span
class="math inline">\(\T(x\to x&#39;) = q(x\to x&#39;)\A(x \to
x&#39;)\)</span>.</p>
<p>The Metropolis-Hasting algorithm is a slight extension of the
original Metropolis algorithm that allows for non-symmetric proposal
distributions $q(xx’) q(x’x) $. It can be derived starting from detailed
balance<span class="citation"
data-cites="krauth_introduction_1996"><sup><a
href="#ref-krauth_introduction_1996"
role="doc-biblioref"><strong>krauth_introduction_1996?</strong></a></sup></span>}:
<span class="math display">\[
P(x)\T(x \to x&#39;) &amp;= P(x&#39;)\T(x&#39; \to x) \\
P(x)q(x \to x&#39;)\A(x \to x&#39;) &amp;= P(x&#39;)q(x&#39; \to
x)\A(x&#39; \to x) \\
\label{eq:db2} \frac{\A(x \to x&#39;)}{\A(x&#39; \to x)} &amp;=
\frac{P(x&#39;)q(x&#39; \to x)}{P(x)q(x \to x&#39;)} = f(x, x&#39;)\\
\]</span> % The Metropolis-Hastings algorithm is the choice: <span
class="math display">\[\label{eq:mh} \A(x \to x&#39;) = \min\left(1,
f(x,x&#39;)\right)\]</span> % Noting that <span
class="math inline">\(f(x,x&#39;) = 1/f(x&#39;,x)\)</span>, Eq. <span
class="math inline">\(\ref{eq:mh}\)</span> can be seen to satisfy Eq.
<span class="math inline">\(\ref{eq:db2}\)</span> by considering the two
cases <span class="math inline">\(f(x,x&#39;) &gt; 1\)</span> and <span
class="math inline">\(f(x,x&#39;) &lt; 1\)</span>.</p>
<p>By choosing the proposal distribution such that <span
class="math inline">\(f(x,x&#39;)\)</span> is as close as possible to
one, the rate of rejections can be reduced and the algorithm sped
up.</p>
%
<p>%Thinning, burn in, multiple runs %Simulated annealing and Parallel
Tempering</p>
<h3 id="applying-mcmc-to-the-fk-model">Applying MCMC to the FK
model}</h3>
<p>MCMC can be applied to sample over the classical degrees of freedom
of the model. We take the full Hamiltonian and split it into a classical
and a quantum part: <span class="math display">\[
H_{\mathrm{FK}} &amp;= -\sum_{&lt;ij&gt;} c^\dagger_{i}c_{j} + U
\sum_{i} (c^\dagger_{i}c_{i} - 1/2)( n_i - 1/2) \\
&amp;+ \sum_{ij} J_{ij} (n_i - 1/2) (n_j - 1/2) - \mu \sum_i
(c^\dagger_{i}c_{i} + n_i)\\
H_q &amp;= -\sum_{&lt;ij&gt;} c^\dagger_{i}c_{j} + \sum_{i}
\left(U(n_i - 1/2) - \mu\right) c^\dagger_{i}c_{i}\\
H_c &amp;= \sum_i \mu n_i - \frac{U}{2}(n_i - 1/2) +
\sum_{ij}J_{ij}(n_i - 1/2)(n_j - 1/2)
\]</span> % There are <span class="math inline">\(2^N\)</span> possible
ion configurations <span class="math inline">\(\{ n_i \}\)</span>, we
define <span class="math inline">\(n^k_i\)</span> to be the occupation
of the ith site of the kth configuration. The quantum part of the free
energy can then be defined through the quantum partition function <span
class="math inline">\(\Z^k\)</span> associated with each ionic state
<span class="math inline">\(n^k_i\)</span>: <span
class="math display">\[
F^k &amp;= -1/\beta \ln{\Z^k} \\
\]</span> % Such that the overall partition function is: <span
class="math display">\[
\Z &amp;= \sum_k e^{- \beta H^k} Z^k \\
&amp;= \sum_k e^{-\beta (H^k + F^k)} \\
\]</span> % Because fermions are limited to occupation numbers of 0 or 1
<span class="math inline">\(Z^k\)</span> simplifies nicely. If <span
class="math inline">\(m^j_i = \{0,1\}\)</span> is defined as the
occupation of the level with energy <span
class="math inline">\(\epsilon^k_i\)</span> then the partition function
is a sum over all the occupation states labelled by j: <span
class="math display">\[
Z^k &amp;= \Tr e^{-\beta F^k} = \sum_j e^{-\beta \sum_i m^j_i
\epsilon^k_i}\\
&amp;= \sum_j \prod_i e^{- \beta m^j_i \epsilon^k_i}= \prod_i
\sum_j e^{- \beta m^j_i \epsilon^k_i}\\
&amp;= \prod_i (1 + e^{- \beta \epsilon^k_i})\\
F^k &amp;= -1/\beta \sum_k \ln{(1 + e^{- \beta \epsilon^k_i})}
\]</span> % Observables can then be calculated from the partition
function, for examples the occupation numbers:</p>
<p><span class="math display">\[
\tex{N} &amp;= \frac{1}{\beta} \frac{1}{Z} \frac{\partial Z}{\partial
\mu} = - \frac{\partial F}{\partial \mu}\\
&amp;= \frac{1}{\beta} \frac{1}{Z} \frac{\partial}{\partial \mu}
\sum_k e^{-\beta (H^k + F^k)}\\
&amp;= 1/Z \sum_k (N^k_{\mathrm{ion}} + N^k_{\mathrm{electron}})
e^{-\beta (H^k + F^k)}\\
\]</span> % with the definitions:</p>
<p><span class="math display">\[
N^k_{\mathrm{ion}} &amp;= - \frac{\partial H^k}{\partial \mu} = \sum_i
n^k_i\\
N^k_{\mathrm{electron}} &amp;= - \frac{\partial F^k}{\partial \mu} =
\sum_i \left(1 + e^{\beta \epsilon^k_i}\right)^{-1}\\
\]</span> % The MCMC algorithm consists of performing a random walk over
the states <span class="math inline">\(\{ n^k_i \}\)</span>. In the
simplest case the proposal distribution corresponds to flipping a random
site from occupied to unoccupied or vice versa, since this proposal is
symmetric the acceptance function becomes: <span class="math display">\[
P(k) &amp;= \Z^{-1} e^{-\beta(H^k + F^k)} \\
\A(k \to k&#39;) &amp;= \min\left(1, \frac{P(k&#39;)}{P(k)}\right) =
\min\left(1, e^{\beta(H^{k&#39;} + F^{k&#39;})-\beta(H^k + F^k)}\right)
\]</span> % At each step <span class="math inline">\(F^k\)</span> is
calculated by diagonalising the tri-diagonal matrix representation of
<span class="math inline">\(H_q\)</span> with open boundary conditions.
Observables are simply averages over the their value at each step of the
random walk. The full spectrum and eigenbasis is too large to save to
disk so usually running averages of key observables are taken as the
walk progresses.</p>
<p>In a MCMC method a key property is the proportion of the time that
proposals are accepted, the acceptance rate. If this rate is too low the
random walk is trying to take overly large steps in energy space which
problematic because it means very few new samples will be generated. If
it is too high it implies the steps are too small, a problem because
then the walk will take longer to explore the state space and the
samples will be highly correlated. Ideal values for the acceptance rate
can be calculated under certain assumptions<span class="citation"
data-cites="roberts_weak_1997"><sup><a href="#ref-roberts_weak_1997"
role="doc-biblioref"><strong>roberts_weak_1997?</strong></a></sup></span>}.
Here we monitor the acceptance rate and if it is too high we re-run the
MCMC with a modified proposal distribution that has a chance to propose
moves that flip multiple sites at a time.</p>
<p>In addition we exploit the particle-hole symmetry of the problem by
occasionally proposing a flip of the entire state. This works because
near half-filling, flipping the occupations of all the sites will
produce a state at or near the energy of the current one.</p>
<p>The matrix diagonalisation is the most computationally expensive step
of the process, a speed up can be obtained by modifying the proposal
distribution to depend on the classical part of the energy, a trick
gleaned from Ref.<span class="citation"
data-cites="krauth_introduction_1996"><sup><a
href="#ref-krauth_introduction_1996"
role="doc-biblioref"><strong>krauth_introduction_1996?</strong></a></sup></span>}:
<span class="math display">\[
q(k \to k&#39;) &amp;= \min\left(1, e^{\beta (H^{k&#39;} - H^k)}\right)
\\
\A(k \to k&#39;) &amp;= \min\left(1, e^{\beta(F^{k&#39;}- F^k)}\right)
\]</span> % This allows the method to reject some states without
performing the diagonalisation at no cost to the accuracy of the MCMC
method.</p>
<p>An extension of this idea is to try to define a classical model with
a similar free energy dependence on the classical state as the full
quantum, Ref.<span class="citation"
data-cites="huang_accelerated_2017"><sup><a
href="#ref-huang_accelerated_2017"
role="doc-biblioref"><strong>huang_accelerated_2017?</strong></a></sup></span>}
does this with restricted Boltzmann machines whose form is very similar
to a classical spin model.</p>
<p>In order to reduce the effects of the boundary conditions and the
finite size of the system we redefine and normalise the coupling matrix
to have 0 derivative at its furthest extent rather than cutting off
abruptly.</p>
<p><span class="math display">\[
J&#39;(x) &amp;= \abs{\frac{L}{\pi}\sin \frac{\pi x}{L}}^{-\alpha} \\
J(x) &amp;= \frac{J_0 J&#39;(x)}{\sum_y J&#39;(y)}
\]</span> % The scaling ensures that, in the ordered phase, the overall
potential felt by each site due to the rest of the system is independent
of system size.</p>
<p>The Binder cumulant is defined as: <span class="math display">\[U_B =
1 - \frac{\tex{\mu_4}}{3\tex{\mu_2}^2}\]</span> % where <span
class="math display">\[\mu_n = \tex{(m - \tex{m})^n}\]</span> % are the
central moments of the order parameter m: <span class="math display">\[m
= \sum_i (-1)^i (2n_i - 1) / N\]</span> % The Binder cumulant evaluated
against temperature can be used as a diagnostic for the existence of a
phase transition. If multiple such curves are plotted for different
system sizes, a crossing indicates the location of a critical point<span
class="citation" data-cites="binder_finite_1981"><sup><a
href="#ref-binder_finite_1981"
role="doc-biblioref"><strong>binder_finite_1981?</strong></a></sup></span>,
musial_monte_2002}.</p>
<h2 id="markov-chain-monte-carlo-in-practice">Markov Chain Monte-Carlo
in Practice}</h2>
<h3 id="quick-intro-to-mcmc">Quick Intro to MCMC}</h3>
<p>The main paper relies on extensively to evaluate thermal expectation
values within the model by walking over states of the classical spin
system <span class="math inline">\(S_i\)</span>. For a classical system,
the thermal expectation value of some operator <span
class="math inline">\(O\)</span> is defined by a Boltzmann weighted sum
over the classical state space: <span class="math display">\[
\tex{O} &amp;= \frac{1}{\Z} \sum_{\s \in S} O(x) P(x) \\
P(x) &amp;= \frac{1}{\Z} e^{-\beta F(x)} \\
\Z &amp;= \sum_{\s \in S} e^{-\beta F(x)}
\]</span> While for a quantum system these sums are replaced by
equivalent traces. The obvious approach to evaluate these sums
numerically would be to directly loop over all the classical states in
the system and perform the sum. But we all know know why this isn’t
feasible: the state space is too large! Indeed even if we could do it,
it would still be computationally wasteful since at low temperatures the
sums are dominated by low energy excitations about the ground states of
the system. Even worse, in our case we must fully solve the fermionic
system via exact diagonalisation for each classical state in the sum, a
very expensive operation!~.}</p>
<p> sidesteps these issues by defining a random walk that focuses on the
states with the greatest Boltzmann weight. At low temperatures this
means we need only visit a few low energy states to make good estimates
while at high temperatures the weights become uniform so a small number
of samples distributed across the state space suffice. However we will
see that the method is not without difficulties of its own.</p>
<p>%MCMC from an ensemble point of view In implementation can be boiled
down to choosing a transition function $(_{t} _t+1) $ where <span
class="math inline">\(\s\)</span> are vectors representing classical
spin configurations. We start in some initial state <span
class="math inline">\(\s_0\)</span> and then repeatedly jump to new
states according to the probabilities given by <span
class="math inline">\(\T\)</span>. This defines a set of random walks
<span class="math inline">\(\{\s_0\ldots \s_i\ldots \s_N\}\)</span>.
Fig.~<span class="math inline">\(\ref{fig:single}\)</span> shows this in
practice: we have a (rather small) ensemble of <span
class="math inline">\(M = 2\)</span> walkers starting at the same point
in state space and then spreading outwards by flipping spins along the
way.</p>
<p>In pseudo-code one could write the MCMC simulation for a single
walker as:</p>
<p>Where the function here produces a state with probability determined
by the and the transition function <span
class="math inline">\(\T\)</span>.</p>
<p>If we ran many such walkers in parallel we could then approximate the
distribution <span class="math inline">\(p_t(\s; \s_0)\)</span> which
tells us where the walkers are likely to be after they’ve evolved for
<span class="math inline">\(t\)</span> steps from an initial state <span
class="math inline">\(\s_0\)</span>. We need to carefully choose <span
class="math inline">\(\T\)</span> such that after a large number of
steps <span class="math inline">\(k\)</span> (the convergence time) the
probability <span class="math inline">\(p_t(\s;\s_0)\)</span> approaches
the thermal distribution <span class="math inline">\(P(\s; \beta) =
\Z^{-1} e^{-\beta F(\s)}\)</span>. This turns out to be quite easy to
achieve using the Metropolis-Hasting algorithm.</p>
<h3 id="convergence-time">Convergence Time}</h3>
<p>Considering <span class="math inline">\(p(\s)\)</span> as a vector
<span class="math inline">\(\vec{p}\)</span> whose jth entry is the
probability of the jth state <span class="math inline">\(p_j =
p(\s_j)\)</span>, and writing <span class="math inline">\(\T\)</span> as
the matrix with entries <span class="math inline">\(T_{ij} = \T(\s_j
\rightarrow \s_i)\)</span> we can write the update rule for the ensemble
probability as: <span class="math display">\[\vec{p}_{t+1} = \T
\vec{p}_t \implies \vec{p}_{t} = \T^t \vec{p}_0\]</span> where <span
class="math inline">\(\vec{p}_0\)</span> is vector which is one on the
starting state and zero everywhere else. Since all states must
transition to somewhere with probability one: <span
class="math inline">\(\sum_i T_{ij} = 1\)</span>.</p>
<p>Matrices that satisfy this are called stochastic matrices exactly
because they model these kinds of Markov processes. It can be shown that
they have real eigenvalues, and ordering them by magnitude, that <span
class="math inline">\(\lambda_0 = 1\)</span> and <span
class="math inline">\(0 &lt; \lambda_{i\neq0} &lt; 1\)</span>.
%https://en.wikipedia.org/wiki/Stochastic_matrix Assuming <span
class="math inline">\(\T\)</span> has been chosen correctly, its single
eigenvector with eigenvalue 1 will be the thermal distribution so
repeated application of the transition function eventually leads there,
while memory of the initial conditions decays exponentially with a
convergence time <span class="math inline">\(k\)</span> determined by
<span class="math inline">\(\lambda_1\)</span>. In practice this means
that one throws away the data from the beginning of the random walk in
order reduce the dependence on the initial conditions and be close
enough to the target distribution.</p>
<h3 id="auto-correlation-time">Auto-correlation Time}</h3>
<p>At this stage one might think we’re done. We can indeed draw
independent samples from <span class="math inline">\(P(\s;
\beta)\)</span> by starting from some arbitrary initial state and doing
<span class="math inline">\(k\)</span> steps to arrive at a sample.
However a key insight is that after the convergence time, every state
generated is a sample from <span class="math inline">\(P(\s;
\beta)\)</span>! They are not, however, independent samples. In
Fig.~<span class="math inline">\(\ref{fig:raw}\)</span> it is already
clear that the samples of the order parameter m have some
auto-correlation because only a few spins are flipped each step but even
when the number of spins flipped per step is increased, Fig.~<span
class="math inline">\(\ref{fig:m_autocorr}\)</span> shows that it can be
an important effect near the phase transition. Let’s define the
auto-correlation time <span class="math inline">\(\tau(O)\)</span>
informally as the number of MCMC samples of some observable O that are
statistically equal to one independent sample.~ for a more rigorous
definition involving a sum over the auto-correlation function.} The
auto-correlation time is generally shorter than the convergence time so
it therefore makes sense from an efficiency standpoint to run a single
walker for many MCMC steps rather than to run a huge ensemble for <span
class="math inline">\(k\)</span> steps each.</p>
<p>Once the random walk has been carried out for many steps, the
expectation values of <span class="math inline">\(O\)</span> can be
estimated from the MCMC samples <span
class="math inline">\(\s_i\)</span>: <span class="math display">\[
\tex{O} = \sum_{i = 0}^{N} O(\s_i) + \mathcal{O}(\frac{1}{\sqrt{N}})
\]</span> The the samples are correlated so the N of them effectively
contains less information than <span class="math inline">\(N\)</span>
independent samples would, in fact roughly <span
class="math inline">\(N/\tau\)</span> effective samples. As a
consequence the variance is larger than the <span
class="math inline">\(\qex{O^2} - \qex{O}^2\)</span> form it would have
if the estimates were uncorrelated. There are many methods in the
literature for estimating the true variance of <span
class="math inline">\(\qex{O}\)</span> and deciding how many steps are
needed but my approach has been to run a small number of parallel
chains, which are independent, in order to estimate the statistical
error produced. This is a slightly less computationally efficient
because it requires throwing away those <span
class="math inline">\(k\)</span> steps generated before convergence
multiple times but it is a conceptually simple workaround.</p>
<p>In summary, to do efficient simulations we want to reduce both the
convergence time and the auto-correlation time as much as possible. In
order to explain how, we need to introduce the Metropolis-Hasting (MH)
algorithm and how it gives an explicit form for the transition
function.</p>
<h3 id="the-metropolis-hastings-algorithm">The Metropolis-Hastings
Algorithm}</h3>
<p>MH breaks up the transition function into a proposal distribution
<span class="math inline">\(q(\s \to \s&#39;)\)</span> and an acceptance
function <span class="math inline">\(\A(\s \to \s&#39;)\)</span>. <span
class="math inline">\(q\)</span> needs to be something that we can
directly sample from, and in our case generally takes the form of
flipping some number of spins in <span
class="math inline">\(\s\)</span>, i.e if we’re flipping a single random
spin in the spin chain, <span class="math inline">\(q(\s \to
\s&#39;)\)</span> is the uniform distribution on states reachable by one
spin flip from <span class="math inline">\(\s\)</span>. This also gives
the nice symmetry property that <span class="math inline">\(q(\s \to
\s&#39;) = q(\s&#39; \to \s)\)</span>.</p>
<p>The proposal <span class="math inline">\(\s&#39;\)</span> is then
accepted or rejected with an acceptance probability <span
class="math inline">\(\A(\s \to \s&#39;)\)</span>, if the proposal is
rejected then <span class="math inline">\(\s_{i+1} = \s_{i}\)</span>.
Hence:</p>
<p><span class="math display">\[\T(x\to x&#39;) = q(x\to x&#39;)\A(x \to
x&#39;)\]</span></p>
<p>When the proposal distribution is symmetric as ours is, it cancels
out in the expression for the acceptance function and the
Metropolis-Hastings algorithm is simply the choice: <span
class="math display">\[ \A(x \to x&#39;) = \min\left(1,
e^{-\beta\;\Delta F}\right)\]</span> Where <span
class="math inline">\(F\)</span> is the overall free energy of the
system, including both the quantum and classical sector.</p>
<p>To implement the acceptance function in practice we pick a random
number in the unit interval and accept if it is less than <span
class="math inline">\(e^{-\beta\;\Delta F}\)</span>:</p>
<p>This has the effect of always accepting proposed states that are
lower in energy and sometimes accepting those that are higher in energy
than the current state.</p>
<h3 id="choosing-the-proposal-distribution">Choosing the proposal
distribution}</h3>
<p>Now we can discuss how to minimise the auto-correlations. The general
principle is that one must balance the proposal distribution between two
extremes. Choose overlay small steps, like flipping only a single spin
and the acceptance rate will be high because <span
class="math inline">\(\Delta F\)</span> will usually be small, but each
state will be very similar to the previous and the auto-correlations
will be high too, making sampling inefficient. On the other hand,
overlay large steps, like randomising a large portion of the spins each
step, will result in very frequent rejections, especially at low
temperatures.</p>
<p>I evaluated a few different proposal distributions for use with the
FK model.</p>
<p>Fro Figure~<span class="math inline">\(\ref{fig:comparison}\)</span>
we see that even at moderately high temperatures <span
class="math inline">\(T &gt; T_c\)</span> flipping one or two sites is
the best choice. However for some simulations at very high temperature
flipping more spins is warranted. Tuning the proposal distribution
automatically seems like something that would not yield enough benefit
for the additional complexity it would require.</p>
<h3 id="two-step-trick">Two Step Trick</h3>
<p>Our method already relies heavily on the split between the classical
and quantum sector to derive a sign problem free MCMC algorithm but it
turns out that there is a further trick we can play with it. The free
energy term is the sum of an easy to compute classical energy and a more
expensive quantum free energy, we can split the acceptance function into
two in such as way as to avoid having to compute the full exact
diagonalisation some of the time:</p>
<div class="sourceCode" id="cb2"><pre
class="sourceCode python"><code class="sourceCode python"><span id="cb2-1"><a href="#cb2-1" aria-hidden="true" tabindex="-1"></a></span>
<span id="cb2-2"><a href="#cb2-2" aria-hidden="true" tabindex="-1"></a>current_state <span class="op">=</span> initial_state</span>
<span id="cb2-3"><a href="#cb2-3" aria-hidden="true" tabindex="-1"></a></span>
<span id="cb2-4"><a href="#cb2-4" aria-hidden="true" tabindex="-1"></a><span class="cf">for</span> i <span class="kw">in</span> <span class="bu">range</span>(N_steps):</span>
<span id="cb2-5"><a href="#cb2-5" aria-hidden="true" tabindex="-1"></a> new_state <span class="op">=</span> proposal(current_state)</span>
<span id="cb2-6"><a href="#cb2-6" aria-hidden="true" tabindex="-1"></a></span>
<span id="cb2-7"><a href="#cb2-7" aria-hidden="true" tabindex="-1"></a> df_classical <span class="op">=</span> classical_free_energy_change(current_state, new_state, parameters)</span>
<span id="cb2-8"><a href="#cb2-8" aria-hidden="true" tabindex="-1"></a> <span class="cf">if</span> exp(<span class="op">-</span>beta <span class="op">*</span> df_classical) <span class="op">&lt;</span> uniform(<span class="dv">0</span>,<span class="dv">1</span>):</span>
<span id="cb2-9"><a href="#cb2-9" aria-hidden="true" tabindex="-1"></a> f_quantum <span class="op">=</span> quantum_free_energy(current_state, new_state, parameters)</span>
<span id="cb2-10"><a href="#cb2-10" aria-hidden="true" tabindex="-1"></a> </span>
<span id="cb2-11"><a href="#cb2-11" aria-hidden="true" tabindex="-1"></a> <span class="cf">if</span> exp(<span class="op">-</span> beta <span class="op">*</span> df_quantum) <span class="op">&lt;</span> uniform(<span class="dv">0</span>,<span class="dv">1</span>):</span>
<span id="cb2-12"><a href="#cb2-12" aria-hidden="true" tabindex="-1"></a> current_state <span class="op">=</span> new_state</span>
<span id="cb2-13"><a href="#cb2-13" aria-hidden="true" tabindex="-1"></a> </span>
<span id="cb2-14"><a href="#cb2-14" aria-hidden="true" tabindex="-1"></a> states[i] <span class="op">=</span> current_state</span>
<span id="cb2-15"><a href="#cb2-15" aria-hidden="true" tabindex="-1"></a> </span></code></pre></div>
<p>lets cite Figure<a href="#fig:binder">1</a></p>
<p>lets cite to person<span class="citation"
data-cites="trebstKitaevMaterials2022"><sup><a
href="#ref-trebstKitaevMaterials2022"
role="doc-biblioref">2</a></sup></span>. and then multple<span
class="citation"
data-cites="banerjeeProximateKitaevQuantum2016 trebstKitaevMaterials2022"><sup><a
href="#ref-trebstKitaevMaterials2022" role="doc-biblioref">2</a>,<a
href="#ref-banerjeeProximateKitaevQuantum2016"
role="doc-biblioref">3</a></sup></span>. what is we surround it by
spaces?<span class="citation"
data-cites="trebstKitaevMaterials2022"><sup><a
href="#ref-trebstKitaevMaterials2022"
role="doc-biblioref">2</a></sup></span></p>
<div class="sourceCode" id="cb3"><pre
class="sourceCode python"><code class="sourceCode python"></code></pre></div>
<div class="sourceCode" id="cb4"><pre
class="sourceCode python"><code class="sourceCode python"></code></pre></div>
<p></ij></ij></ij></ij></ij></ij></ij></p>
<div id="refs" class="references csl-bib-body" data-line-spacing="2"
role="doc-bibliography">
<div id="ref-kagan1984localization" class="csl-entry"
role="doc-biblioentry">
<div class="csl-left-margin">1. </div><div
class="csl-right-inline">Kagan, Y. &amp; Maksimov, L. Localization in a
system of interacting particles diffusing in a regular crystal.
<em>Zhurnal Eksperimental’noi i Teoreticheskoi Fiziki</em>
<strong>87</strong>, 348–365 (1984).</div>
</div>
<div id="ref-trebstKitaevMaterials2022" class="csl-entry"
role="doc-biblioentry">
<div class="csl-left-margin">2. </div><div
class="csl-right-inline">Trebst, S. &amp; Hickey, C. <a
href="https://doi.org/10.1016/j.physrep.2021.11.003">Kitaev
materials</a>. <em>Physics Reports</em> <strong>950</strong>, 1–37
(2022).</div>
</div>
<div id="ref-banerjeeProximateKitaevQuantum2016" class="csl-entry"
role="doc-biblioentry">
<div class="csl-left-margin">3. </div><div
class="csl-right-inline">Banerjee, A. <em>et al.</em> <a
href="https://doi.org/10.1038/nmat4604">Proximate <span>Kitaev Quantum
Spin Liquid Behaviour</span> in {\alpha}-<span>RuCl</span>$_3$</a>.
<em>Nature Mater</em> <strong>15</strong>, 733–740 (2016).</div>
</div>
</div>
</main>
</body>
</html>
Reference in New Issue View Git Blame Copy Permalink