
.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "examples\06-PyPoscar\plot_rdf_cutoff_pyposcar.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

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

        :ref:`Go to the end <sphx_glr_download_examples_06-PyPoscar_plot_rdf_cutoff_pyposcar.py>`
        to download the full example code.

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

.. _sphx_glr_examples_06-PyPoscar_plot_rdf_cutoff_pyposcar.py:


.. _ref_example_rdf:

Analyzing Radial Distribution Functions (RDF)
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

In this example, we'll explore the radial distribution functions (RDF) of an atomic structure. The RDF provides insights into the probability of finding an atom at a certain distance from another. This can be useful for understanding the local environment of atoms.

We'll use the `pyprocar` package to:

1. Parse the POSCAR file and obtain atomic positions.
2. Compute the general KDE curve for all distances.
3. Compute KDE curves for specific atomic species interactions.
4. Visualize the results using `matplotlib`.

Let's dive in!

.. GENERATED FROM PYTHON SOURCE LINES 18-33

.. code-block:: Python


    import itertools
    import os

    import matplotlib.pyplot as plt
    import numpy as np
    from matplotlib.gridspec import GridSpec

    import pyprocar.pyposcar as p
    from pyprocar.utils import DATA_DIR

    colors = ["b", "r", "y", "g", "m", "k", "c"]

    data_dir = os.path.join(DATA_DIR, "examples", "PyPoscar", "04-rdf")








.. GENERATED FROM PYTHON SOURCE LINES 34-36

Parsing the POSCAR file
+++++++++++++++++++++++

.. GENERATED FROM PYTHON SOURCE LINES 36-40

.. code-block:: Python


    a = p.poscar.Poscar(os.path.join(data_dir, "POSCAR-AGNR-defect.vasp"))
    a.parse()








.. GENERATED FROM PYTHON SOURCE LINES 41-43

Computing RDF and KDE Curves
++++++++++++++++++++++++++++++

.. GENERATED FROM PYTHON SOURCE LINES 43-70

.. code-block:: Python


    my_rdf = p.rdf.RDF(poscar=a)
    thresholdSp = my_rdf.neighbor_thresholdSp

    # Calculates KDE curves for each species interaction
    spCurves = my_rdf.KDE_CurveSp()
    # All calculations share the same kde space
    domain = my_rdf.KDE_space
    print("CutOff values for each interaction")
    for inter, cutoff in zip(my_rdf.interactions, thresholdSp):
        print(inter, " ", cutoff)


    print(
        "As you can see most interactions have meaaningfull values,"
        " but there is no way to defene nearest neighbors for others (e.g. H-B)."
        " A large cutoff is provided in these cases."
    )

    print(
        "Note: The class for finding nearest neighbors uses a database to ignore these large distances as `nearest neighbors`"
    )

    NeighborsClass = p.latticeUtils.Neighbors(poscar=a)
    NeighborsClass.estimateMaxBondDist()
    Neighbor_threshold = NeighborsClass.d_Max





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

 .. code-block:: none

    CutOff values for each interaction
    ['C' 'C']   1.9500000000000002
    ['C' 'H']   1.6500000000000001
    ['C' 'B']   1.9500000000000002
    ['C' 'N']   1.9500000000000002
    ['H' 'H']   2.2
    ['H' 'B']   5.95
    ['H' 'N']   5.95
    ['B' 'N']   1.4772259999999975
    As you can see most interactions have meaaningfull values, but there is no way to defene nearest neighbors for others (e.g. H-B). A large cutoff is provided in these cases.
    Note: The class for finding nearest neighbors uses a database to ignore these large distances as `nearest neighbors`




.. GENERATED FROM PYTHON SOURCE LINES 71-73

Visualizing the Results
++++++++++++++++++++++++

.. GENERATED FROM PYTHON SOURCE LINES 73-155

.. code-block:: Python


    gs = GridSpec(2, 2, height_ratios=[1, 2])

    ax1 = plt.subplot(gs[0, 0])
    ax2 = plt.subplot(gs[1, 0])
    ax3 = plt.subplot(gs[:, 1])

    # fig,axes = plt.subplots(2,2)


    # This is the kde curve for all the distances in POSCAR
    ax1.plot(domain, my_rdf.KDE_Curve())
    ax1.set_title("General KDE curve")

    cc = itertools.cycle(colors)
    for i, curve in enumerate(spCurves):
        label = my_rdf.interactions[i]
        c = next(cc)
        ax2.plot(domain, curve, label=label, color=c)
        ax2.set_title("Species KDE curves")
        ax2.legend()


    Neighbor_threshold = np.ndarray.flatten(Neighbor_threshold)
    Neighbor_threshold = set(Neighbor_threshold)
    print("The DataBase Thresholds are:")
    print(Neighbor_threshold)

    # Distances for Histogram

    distance_matrix = my_rdf.distances
    distance_matrix = np.array(distance_matrix)
    distances = np.ndarray.flatten(distance_matrix)
    bins = 100
    n, bins, patches = ax3.hist(distances, bins=bins, edgecolor="black", alpha=0.7)

    ax3.set_xlabel("Distance Between Atoms")

    # Plotting CutOff point vertical lines

    labels = my_rdf.interactions
    labels = labels.tolist()

    Plotlabels = []
    for interaction in labels:
        Plotlabels.append(interaction[0] + "-" + interaction[1])

    cc = itertools.cycle(colors)
    var_cutoffs = []
    for i, cutOff in enumerate(thresholdSp):
        c = next(cc)
        var_cutoffs.append(
            ax3.plot([cutOff, cutOff], [0, 1000], color=c, label=my_rdf.interactions[i])
        )
    leg1 = ax3.legend(
        title="CutOff points",
        handles=[x[0] for x in var_cutoffs],
        labels=Plotlabels,
        loc="upper left",
    )

    cc = itertools.cycle(colors)
    data_cutoffs = []


    for i, cutOff in enumerate(Neighbor_threshold):
        c = next(cc)
        data_cutoffs.append(
            ax3.plot([cutOff, cutOff], [0, 1000], color="grey", linestyle="dotted")
        )

    # leg2 = ax3.legend(title = 'DataBased CutOff', handles = [x[0] for x in data_cutoffs], labels = color_label,loc = 'upper right')

    ax3.add_artist(leg1)
    # ax3.add_artist(leg2)

    ax1.grid(True)
    ax2.grid(True)
    ax3.grid(True)

    plt.tight_layout()
    plt.show()



.. image-sg:: /examples/06-PyPoscar/images/sphx_glr_plot_rdf_cutoff_pyposcar_001.png
   :alt: General KDE curve, Species KDE curves
   :srcset: /examples/06-PyPoscar/images/sphx_glr_plot_rdf_cutoff_pyposcar_001.png
   :class: sphx-glr-single-img


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

 .. code-block:: none

    C:\Users\lllang\Desktop\Current_Projects\pyprocar\examples\06-PyPoscar\plot_rdf_cutoff_pyposcar.py:91: MatplotlibDeprecationWarning: Passing label as a length 2 sequence when plotting a single dataset is deprecated in Matplotlib 3.9 and will error in 3.11.  To keep the current behavior, cast the sequence to string before passing.
      ax2.plot(domain, curve, label=label, color=c)
    The DataBase Thresholds are:
    {0.7725483399593904, 1.8348023074035522, 1.9434419177103412, 1.7744469683442248, 1.3036753236814713, 1.4123149339882604, 1.2433199846221439, 2.0520815280171307, 1.883086578651014, 1.7140916292848973}
    C:\Users\lllang\Desktop\Current_Projects\pyprocar\examples\06-PyPoscar\plot_rdf_cutoff_pyposcar.py:125: MatplotlibDeprecationWarning: Passing label as a length 2 sequence when plotting a single dataset is deprecated in Matplotlib 3.9 and will error in 3.11.  To keep the current behavior, cast the sequence to string before passing.
      ax3.plot([cutOff, cutOff], [0, 1000], color=c, label=my_rdf.interactions[i])





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

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


.. _sphx_glr_download_examples_06-PyPoscar_plot_rdf_cutoff_pyposcar.py:

.. only:: html

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

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

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

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

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

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

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


.. only:: html

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

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