{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"nbsphinx": "hidden"
},
"outputs": [],
"source": [
"%matplotlib inline\n",
"%config InlineBackend.figure_format = 'svg'\n",
"\n",
"import warnings \n",
"warnings.filterwarnings(\"ignore\") #ignore some matplotlib warnings\n",
"\n",
"import numpy as np\n",
"\n",
"from triqs.plot.mpl_interface import plt\n",
"plt.rcParams[\"figure.figsize\"] = (6,5) # set default size for all figures\n",
"\n",
"from triqs.utility.redirect import start_redirect\n",
"start_redirect()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Perturbation order histograms\n",
"\n",
"To analyze the behaviour of the Markov chain in CTHYB the perturbation order histograms are a central tool. This is a simple example showing how to sample and plot these histograms.\n",
"\n",
"As an example we solve the one-orbital Anderson impurity at inverse temperature $\\beta$, embedded in a flat (Wilson) bath with non-interacting Green's function $G^{-1}_{0,\\sigma} (i\\omega_n) = i \\omega_n - \\epsilon_f - V^2 \\Gamma_\\sigma(i \\omega_n)$, and local interaction $H_\\mathrm{int} = U n_\\uparrow n_\\downarrow$, where $U$ is the Hubbard interaction, $\\epsilon_f$ is the local energy. The bath $\\Gamma$, and the bath has bandwidth $D$ and hybridization $V$. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Starting on 1 Nodes at : 2019-06-05 17:20:04.584678\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
"╔╦╗╦═╗╦╔═╗ ╔═╗ ┌─┐┌┬┐┬ ┬┬ ┬┌┐ \n",
" ║ ╠╦╝║║═╬╗╚═╗ │ │ ├─┤└┬┘├┴┐\n",
" ╩ ╩╚═╩╚═╝╚╚═╝ └─┘ ┴ ┴ ┴ ┴ └─┘\n",
"\n",
"The local Hamiltonian of the problem:\n",
"(-1,0)*c_dag('down',0)*c('down',0) + (-1,0)*c_dag('up',0)*c('up',0) + (2,0)*c_dag('down',0)*c_dag('up',0)*c('up',0)*c('down',0)\n",
"Using autopartition algorithm to partition the local Hilbert space\n",
"Found 4 subspaces.\n",
"\n",
"Warming up ...\n",
"\n",
"Accumulating ...\n",
"17:20:04 0% ETA 00:00:14 cycle 667 of 100000\n",
"17:20:06 14% ETA 00:00:12 cycle 14816 of 100000\n",
"17:20:09 31% ETA 00:00:10 cycle 31517 of 100000\n",
"17:20:12 51% ETA 00:00:07 cycle 51752 of 100000\n",
"17:20:16 78% ETA 00:00:03 cycle 78933 of 100000\n",
"\n",
"\n",
"[Rank 0] Collect results: Waiting for all mpi-threads to finish accumulating...\n",
"[Rank 0] Timings for all measures:\n",
"Measure | seconds \n",
"Average sign | 0.0121935 \n",
"Perturbation order | 0.017634 \n",
"Perturbation order (down) | 0.0144553 \n",
"Perturbation order (up) | 0.011125 \n",
"Total measure time | 0.0554078 \n",
"[Rank 0] Acceptance rate for all moves:\n",
"Move set Insert two operators: 0.056355\n",
" Move Insert Delta_up: 0.0565787\n",
" Move Insert Delta_down: 0.0561309\n",
"Move set Remove two operators: 0.0564614\n",
" Move Remove Delta_up: 0.0566373\n",
" Move Remove Delta_down: 0.0562848\n",
"Move Shift one operator: 0.758509\n",
"[Rank 0] Warmup lasted: 0.0125826 seconds [00:00:00]\n",
"[Rank 0] Simulation lasted: 14.9179 seconds [00:00:14]\n",
"[Rank 0] Number of measures: 100000\n",
"Total number of measures: 100000\n",
"Average sign: (1,0)\n"
]
}
],
"source": [
"from triqs.operators import n\n",
"from h5 import HDFArchive\n",
"from triqs.gf import inverse, iOmega_n, Wilson\n",
"\n",
"from triqs_cthyb import Solver\n",
"\n",
"U, e_f, D, V, beta = 2., -1., 1., 1., 20.\n",
"Sz = 0.5 * ( n('up', 0) - n('down', 0) )\n",
"\n",
"S = Solver(\n",
" beta=beta, gf_struct=[('up',[0]), ('down',[0])],\n",
" n_tau=400, n_iw=50,)\n",
"\n",
"S.G0_iw << inverse(iOmega_n - e_f - V**2 * Wilson(D))\n",
"\n",
"S.solve(\n",
" h_int=U*n('up',0)*n('down',0),\n",
" n_cycles=100000,\n",
" length_cycle=20,\n",
" n_warmup_cycles=100,\n",
" measure_G_tau=False,\n",
" measure_pert_order=True,\n",
" )"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Plotting perturbation order histograms\n",
"\n",
"The measured histograms are available as the member properties ``S.perturbation_order`` and ``S.perturbation_order_total`` of the cthyb solver.\n",
"\n",
"For an ergodic Markov chain with sufficient number of warmup sweeps, the perturbation order distribution should be approximately Gaussian (or Possonian for low average orders)."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n",
"\n",
"\n"
],
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"from triqs.plot.mpl_interface import oplot, oplotr, plt\n",
"\n",
"oplot(S.perturbation_order_total, label='total')\n",
"for b in S.perturbation_order:\n",
" oplot(S.perturbation_order[b], label='block {:s}'.format(b))\n",
"\n",
"plt.semilogy([], [])\n",
"plt.xlim([0, 40]);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The total average perturbation order is important information to determine the ``length_cycle`` parameter that controls the number of attempted moves between measurements. To ensure low correlation between measurements the ``length_cycle`` should be high enough that the attempted moves can significanty change the configuration."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Autor: H. U.R. Strand (2019)"
]
}
],
"metadata": {
"celltoolbar": "Edit Metadata",
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3"
}
},
"nbformat": 4,
"nbformat_minor": 2
}