.. _sphx_glr_auto_examples_plot_atom.py:


==================================
Interacting atom spectral function
==================================

For the case of contact interaction in the single orbital case
the atomic Green function as given by the Lehmann Representation.




.. rst-class:: sphx-glr-horizontal


    *

      .. image:: /auto_examples/images/sphx_glr_plot_atom_001.png
            :scale: 47

    *

      .. image:: /auto_examples/images/sphx_glr_plot_atom_002.png
            :scale: 47

    *

      .. image:: /auto_examples/images/sphx_glr_plot_atom_003.png
            :scale: 47





.. code-block:: python

    # author: Óscar Nájera

    from __future__ import division, absolute_import, print_function
    import numpy as np
    import matplotlib.pyplot as plt
    from dmft.common import matsubara_freq, gw_invfouriertrans
    from slaveparticles.quantum import fermion
    from slaveparticles.quantum.operators import gf_lehmann, diagonalize


    def hamiltonian(U, mu):
        r"""Generate a single orbital isolated atom Hamiltonian in particle-hole
        symmetry. Include chemical potential for grand Canonical calculations

        .. math::
            \mathcal{H} - \mu N = -\frac{U}{2}(n_\uparrow - n_\downarrow)^2
            - \mu(n_\uparrow + n_\downarrow)

        """
        d_up, d_dw = [fermion.destruct(2, sigma) for sigma in range(2)]
        sigma_z = d_up.T * d_up - d_dw.T * d_dw
        H = - U / 2 * sigma_z * sigma_z - mu * (d_up.T * d_up + d_dw.T * d_dw)
        return H, d_up, d_dw


    def gf(w, U, mu, beta):
        """Calculate by Lehmann representation the green function"""
        H, d_up, d_dw = hamiltonian(U, mu)
        e, v = diagonalize(H.todense())
        g_up = gf_lehmann(e, v, d_up.T, beta, w)
        return g_up


    beta = 50
    U = 1.
    mu_v = np.array([0, 0.2, 0.45, 0.5, 0.65])
    c_v = ['b', 'g', 'r', 'k', 'm']

    f, axw = plt.subplots(2, sharex=True)
    f.subplots_adjust(hspace=0)
    w = np.linspace(-1.5, 1.5, 500) + 1j * 1e-2
    for mu, c in zip(mu_v, c_v):
        gw = gf(w, U, mu, beta)
        axw[0].plot(w.real, gw.real, c, label=r'$\mu={}$'.format(mu))
        axw[1].plot(w.real, -1 * gw.imag / np.pi, c)
    axw[0].legend()
    axw[0].set_title(r'Real Frequencies Green functions, $\beta={}$'.format(beta))
    axw[0].set_ylabel(r'$\Re e G(\omega)$')
    axw[1].set_ylabel(r'$A(\omega)$')
    axw[1].set_xlabel(r'$\omega$')


    gwp, axwn = plt.subplots(2, sharex=True)
    gwp.subplots_adjust(hspace=0)
    gtp, axt = plt.subplots()
    wn = matsubara_freq(beta, 64)
    tau = np.linspace(0, beta, 2**10)
    for mu, c in zip(mu_v, c_v):
        giw = gf(1j * wn, U, mu, beta)
        axwn[0].plot(wn, giw.real, c + 's-', label=r'$\mu={}$'.format(mu))
        axwn[1].plot(wn, giw.imag, c + 'o-')

        gt = gw_invfouriertrans(giw, tau, wn)
        axt.plot(tau, gt, label=r'$\mu={}$'.format(mu))

    axwn[0].legend()
    axwn[0].set_title(r'Matsubara Green functions, $\beta={}$'.format(beta))
    axwn[1].set_xlabel(r'$\omega_n$')
    axwn[0].set_ylabel(r'$\Re e G(i\omega_n)$')
    axwn[1].set_ylabel(r'$\Im m G(i\omega_n)$')

    axt.set_ylim(top=0.05)
    axt.legend(loc=0)
    axt.set_title(r'Imaginary time Green functions, $\beta={}$'.format(beta))
    axt.set_xlabel(r'$\tau$')
    axt.set_ylabel(r'$G(\tau)$')

**Total running time of the script:** ( 0 minutes  0.531 seconds)



.. container:: sphx-glr-footer


  .. container:: sphx-glr-download

     :download:`Download Python source code: plot_atom.py <plot_atom.py>`



  .. container:: sphx-glr-download

     :download:`Download Jupyter notebook: plot_atom.ipynb <plot_atom.ipynb>`

.. rst-class:: sphx-glr-signature

    `Generated by Sphinx-Gallery <https://sphinx-gallery.readthedocs.io>`_