.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "auto_examples/giese/giese_comparison.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_auto_examples_giese_giese_comparison.py>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_auto_examples_giese_giese_comparison.py:


Comparison of Giese and PTtools solvers
=======================================

.. GENERATED FROM PYTHON SOURCE LINES 5-88



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


    *

      .. image-sg:: /auto_examples/giese/images/sphx_glr_giese_comparison_001.png
         :alt: PTtools, Giese et al.
         :srcset: /auto_examples/giese/images/sphx_glr_giese_comparison_001.png
         :class: sphx-glr-multi-img

    *

      .. image-sg:: /auto_examples/giese/images/sphx_glr_giese_comparison_002.png
         :alt: giese comparison
         :srcset: /auto_examples/giese/images/sphx_glr_giese_comparison_002.png
         :class: sphx-glr-multi-img


.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    /home/docs/checkouts/readthedocs.org/user_builds/pttools/checkouts/latest/pttools/models/const_cs.py:681: RuntimeWarning:

    invalid value encountered in scalar multiply







|

.. code-block:: Python


    import logging

    import matplotlib.pyplot as plt
    import numpy as np

    from examples.utils import save
    from pttools.analysis.parallel import create_bubbles
    from pttools.bubble.bubble_quantities import get_kappa
    from pttools.bubble.giese import kappaNuMuModel
    from pttools.models import ConstCSModel
    from pttools.speedup import GITHUB_ACTIONS

    logger = logging.getLogger(__name__)


    def main():
        if kappaNuMuModel is None:
            logger.warning("Giese et al. code is not available. Skipping comparison.")
        alpha_ns = np.array([0.01, 0.03, 0.1, 0.3, 1, 3])
        colors = ["b", "y", "r", "g", "purple", "grey"]
        n_v_walls = 20 if GITHUB_ACTIONS else 50
        v_walls = np.linspace(0.2, 0.95, n_v_walls)
        a_s = 5
        a_b = 1
        V_s = 1
        models = [
            ConstCSModel(css2=1 / 3, csb2=1 / 3, a_s=a_s, a_b=a_b, V_s=V_s, alpha_n_min=alpha_ns[0]),
            ConstCSModel(css2=1 / 3, csb2=1 / 4, a_s=a_s, a_b=a_b, V_s=V_s, alpha_n_min=alpha_ns[0]),
            ConstCSModel(css2=1 / 4, csb2=1 / 3, a_s=a_s, a_b=a_b, V_s=V_s, alpha_n_min=alpha_ns[0]),
            ConstCSModel(css2=1 / 4, csb2=1 / 4, a_s=a_s, a_b=a_b, V_s=V_s, alpha_n_min=alpha_ns[0])
        ]
        logger.info(f"Minimum alpha_ns: %s", [model.alpha_n_min for model in models])
        for model in models:
            logger.info("Model parameters: %s", model.params_str())

        figsize_x = 10
        fig1: plt.Figure = plt.figure(figsize=(figsize_x, 6))
        fig2: plt.Figure = plt.figure(figsize=(figsize_x, 6))
        axs1 = fig1.subplots(1, 2)
        ax1 = axs1[0]
        ax2 = axs1[1]
        ax3 = fig2.add_subplot()

        kappas_pttools = np.empty((len(models), alpha_ns.size, v_walls.size))
        kappas_giese = np.empty((len(models), alpha_ns.size, v_walls.size))
        for i_model, model in enumerate(models):
            ls = "--" if i_model in [2, 3] else "-"
            bubbles_pttools, kappas_pttools[i_model, :, :] = create_bubbles(
                model=model, v_walls=v_walls, alpha_ns=alpha_ns, func=get_kappa,
                bubble_kwargs={"allow_invalid": False}, allow_bubble_failure=True
            )
            if kappaNuMuModel is not None:
                bubbles_giese, kappas_giese[i_model, :, :] = create_bubbles(
                    model=model, v_walls=v_walls, alpha_ns=alpha_ns, func=get_kappa,
                    bubble_kwargs={"allow_invalid": False, "use_giese_solver": True}, allow_bubble_failure=True
                )
            for i_alpha_n, (alpha_n, color) in enumerate(zip(alpha_ns, colors)):
                kpt = kappas_pttools[i_model, i_alpha_n, :]
                ax1.plot(v_walls, kpt, ls=ls, color=color, alpha=0.5)
                if kappaNuMuModel is not None:
                    kg = kappas_giese[i_model, i_alpha_n, :]
                    ax2.plot(v_walls, kg, ls=ls, color=color, alpha=0.5)
                    ax3.plot(v_walls, (kpt - kg)/kg, color=color, alpha=0.5)

        for ax in axs1.flat:
            ax.set_xlabel(r"$v_\text{wall}$")
            ax.set_ylabel(r"$\kappa$")
            ax.set_xlim(v_walls.min(), v_walls.max())
            ax.set_yscale("log")
        ax1.set_title("PTtools")
        ax2.set_title("Giese et al.")
        ax3.set_xlim(v_walls.min(), v_walls.max())
        ax3.set_yscale("log")

        return fig1, fig2


    if __name__ == "__main__":
        fig, fig2 = main()
        save(fig, "giese_comparison")
        save(fig2, "giese_comparison_diff")
        plt.show()


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (1 minutes 33.157 seconds)

**Estimated memory usage:**  316 MB


.. _sphx_glr_download_auto_examples_giese_giese_comparison.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-jupyter

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

    .. container:: sphx-glr-download sphx-glr-download-python

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

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: giese_comparison.zip <giese_comparison.zip>`


.. only:: html

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

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_