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

.. only:: html

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

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

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

.. _sphx_glr_examples_01_triangle-basic.py:


Basic Triangle Example
======================

In this example we'll create some basic geometries and turn them into meshes.
to illustrate some of the mesh generation features that Triangle provides in
combination with polygon, point, and linestring geometries represented by
geopandas.

.. GENERATED FROM PYTHON SOURCE LINES 11-18

.. code-block:: Python

    import geopandas as gpd
    import matplotlib.pyplot as plt
    import numpy as np
    import shapely.geometry as sg

    import pandamesh as pm








.. GENERATED FROM PYTHON SOURCE LINES 19-24

A simple rectangular mesh
-------------------------

The most simple example is perhaps a rectangle. We'll create a vector
geometry, store this in a geodataframe, and associate a cell size.

.. GENERATED FROM PYTHON SOURCE LINES 24-36

.. code-block:: Python


    polygon = sg.Polygon(
        [
            [0.0, 0.0],
            [10.0, 0.0],
            [10.0, 10.0],
            [0.0, 10.0],
        ]
    )
    gdf = gpd.GeoDataFrame(geometry=[polygon])
    gdf["cellsize"] = 2.0








.. GENERATED FROM PYTHON SOURCE LINES 37-41

We'll use this polygon to generate a mesh. We start by initializing a
TriangleMesher, which is a simple wrapper around the Python bindings to the
Triangle C-library. This wrapper extracts the coordinates and presents them
in the appropriate manner for triangle.

.. GENERATED FROM PYTHON SOURCE LINES 41-46

.. code-block:: Python


    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles)




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_001.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_001.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 47-51

Defaults
--------

The TriangleMesher class is initialized with a number of default parameters:

.. GENERATED FROM PYTHON SOURCE LINES 51-54

.. code-block:: Python


    print(mesher)





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

 .. code-block:: none

    TriangleMesher
        xoff = 5.0
        yoff = 5.0
        vertices = np.ndarray with shape (4, 2)
        segments = np.ndarray with shape (4, 2)
        regions = np.ndarray with shape (1, 4)
        holes = None
        minimum_angle = 20.0
        conforming_delaunay = True
        suppress_exact_arithmetic = False
        maximum_steiner_points = None
        delaunay_algorithm = DelaunayAlgorithm.DIVIDE_AND_CONQUER
        consistency_check = False




.. GENERATED FROM PYTHON SOURCE LINES 55-56

We can change a parameter, and see what effects this has on the mesh:

.. GENERATED FROM PYTHON SOURCE LINES 56-61

.. code-block:: Python


    mesher.conforming_delaunay = False
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles)




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_002.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_002.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 62-64

To generate a mesh with smaller cell sizes, we adjust the geodataframe, and
recreate the mesher.

.. GENERATED FROM PYTHON SOURCE LINES 64-69

.. code-block:: Python


    gdf["cellsize"] = 1.0
    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles)



.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_003.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_003.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 70-75

Multiple cell size zones
------------------------

Multiple zones of cell sizes are supported, as every polygon can be associated
with a cell size in the geodataframe.

.. GENERATED FROM PYTHON SOURCE LINES 75-90

.. code-block:: Python


    polygon2 = sg.Polygon(
        [
            [10.0, 0.0],
            [20.0, 0.0],
            [20.0, 10.0],
            [10.0, 10.0],
        ]
    )
    gdf = gpd.GeoDataFrame(geometry=[polygon, polygon2])
    gdf["cellsize"] = [2.0, 1.0]

    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles)



.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_004.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_004.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 91-95

Polygons with holes ("donut" geometries)
----------------------------------------

Holes in polygons work as expected:

.. GENERATED FROM PYTHON SOURCE LINES 95-107

.. code-block:: Python


    outer = [(0.0, 0.0), (10.0, 0.0), (10.0, 10.0), (0.0, 10.0)]
    inner = [(3.0, 3.0), (7.0, 3.0), (7.0, 7.0), (3.0, 7.0)]

    donut = sg.Polygon(shell=outer, holes=[inner])
    gdf = gpd.GeoDataFrame(geometry=[donut])
    gdf["cellsize"] = [2.0]

    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles)




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_005.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_005.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 108-114

Local refinement
----------------

To do local refinement, we need to ensure there is no overlap between the
polygons. The coordinates of the hole of the outer polygon should match
exactly with the coordinates of the exterior boundary of the inner polygon.

.. GENERATED FROM PYTHON SOURCE LINES 114-124

.. code-block:: Python


    refined = sg.Polygon(inner)

    gdf = gpd.GeoDataFrame(geometry=[donut, refined])
    gdf["cellsize"] = [2.0, 0.5]

    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles)




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_006.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_006.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 125-131

Force points into the triangulation
-----------------------------------

We may also force points into the triangulation, by adding points to the
geodataframe. Let's assume we'd like to a series of points at x=1.0, at a
distance of 0.5.

.. GENERATED FROM PYTHON SOURCE LINES 131-140

.. code-block:: Python


    y = np.arange(0.5, 10.0, 0.5)
    x = np.full(y.size, 1.0)
    points = gpd.points_from_xy(x, y)

    gdf = gpd.GeoDataFrame(geometry=[donut, refined, *points])
    gdf["cellsize"] = [2.0, 0.5] + (len(points) * [np.nan])
    gdf.plot(facecolor="none")




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_007.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_007.png
   :class: sphx-glr-single-img


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

 .. code-block:: none


    <Axes: >



.. GENERATED FROM PYTHON SOURCE LINES 141-143

We can now see the points forced in the triangulation, by plotting the
contents of the geodataframe on top of the generated mesh:

.. GENERATED FROM PYTHON SOURCE LINES 143-150

.. code-block:: Python


    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()

    fig, ax = plt.subplots()
    pm.plot(vertices, triangles, ax=ax)
    gdf.plot(facecolor="none", edgecolor="red", ax=ax)



.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_008.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_008.png
   :class: sphx-glr-single-img


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

 .. code-block:: none


    <Axes: >



.. GENERATED FROM PYTHON SOURCE LINES 151-156

Force linestrings into the triangulation
----------------------------------------

We may do the same with linestrings. Here, we will add a vertical line at
x = 9.0.

.. GENERATED FROM PYTHON SOURCE LINES 156-173

.. code-block:: Python


    line = sg.LineString(
        [
            [9.0, 2.0],
            [9.0, 8.0],
        ]
    )
    gdf = gpd.GeoDataFrame(geometry=[donut, refined, line, *points])
    gdf["cellsize"] = [2.0, 0.5, np.nan] + (len(points) * [np.nan])

    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()

    fig, ax = plt.subplots()
    pm.plot(vertices, triangles, ax=ax)
    gdf.plot(facecolor="none", edgecolor="red", ax=ax)




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_009.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_009.png
   :class: sphx-glr-single-img


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

 .. code-block:: none


    <Axes: >



.. GENERATED FROM PYTHON SOURCE LINES 174-178

Specify cell size along line string
-----------------------------------

Finally, we may also specify the cell size along the line.

.. GENERATED FROM PYTHON SOURCE LINES 178-190

.. code-block:: Python


    line = sg.LineString([(2.0, 8.0), (8.0, 2.0)])
    gdf = gpd.GeoDataFrame(geometry=[polygon, line])
    gdf["cellsize"] = [2.0, 0.5]

    fig, ax = plt.subplots()

    mesher = pm.TriangleMesher(gdf)
    vertices, triangles = mesher.generate()
    pm.plot(vertices, triangles, ax=ax)
    gdf.plot(facecolor="none", edgecolor="red", ax=ax)




.. image-sg:: /examples/images/sphx_glr_01_triangle-basic_010.png
   :alt: 01 triangle basic
   :srcset: /examples/images/sphx_glr_01_triangle-basic_010.png
   :class: sphx-glr-single-img


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

 .. code-block:: none


    <Axes: >



.. GENERATED FROM PYTHON SOURCE LINES 191-198

Conclusion
----------

In real use, the vector geometries will be more complex, and not based on
just a few coordinate pairs. Such cases are presented in the other examples,
but the same principles apply: we may use polygons with associated cell
sizes, and linestrings and points to steer the triangulation.


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

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


.. _sphx_glr_download_examples_01_triangle-basic.py:

.. only:: html

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

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

      :download:`Download Jupyter notebook: 01_triangle-basic.ipynb <01_triangle-basic.ipynb>`

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

      :download:`Download Python source code: 01_triangle-basic.py <01_triangle-basic.py>`

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

      :download:`Download zipped: 01_triangle-basic.zip <01_triangle-basic.zip>`


.. only:: html

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

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