.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "tutorials/c_advanced/d04_topology.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_tutorials_c_advanced_d04_topology.py: Analyzing Geomodel Topology ============================== Extracting adjacency graphs and topology relationships from a computed model This tutorial uses the ``gempy_plugins`` topology analysis module to derive an adjacency graph between the unique geobodies of a faulted model, then visualizes and queries that graph. .. note:: This tutorial relies on ``gempy_plugins``, a separate package maintained in its own repository rather than by the core GemPy developers. .. GENERATED FROM PYTHON SOURCE LINES 14-24 .. code-block:: Python import gempy as gp import gempy_viewer as gpv from gempy_plugins.topology_analysis import topology as tp import os import warnings warnings.filterwarnings("ignore") .. GENERATED FROM PYTHON SOURCE LINES 25-33 Load example Model ^^^^^^^^^^^^^^^^^^ First let's set up a very simple example model. For that we initialize the geo_model object with the correct model extent and the resolution we like. Then we load our data points from csv files and set the series and order the formations (stratigraphic pile). .. GENERATED FROM PYTHON SOURCE LINES 35-62 .. code-block:: Python data_path = os.path.abspath('../../') geo_model = gp.create_geomodel( project_name='Model_Tutorial6', extent=[0, 3000, 0, 20, 0, 2000], resolution=[50, 10, 67], refinement=1, # * For this model is better not to use octrees because we want to see what is happening in the scalar fields importer_helper=gp.data.ImporterHelper( path_to_orientations=data_path + "/data/input_data/tut_chapter6/ch6_data_fol.csv", path_to_surface_points=data_path + "/data/input_data/tut_chapter6/ch6_data_interf.csv", ) ) gp.map_stack_to_surfaces( gempy_model=geo_model, mapping_object= { "fault": "Fault", "Rest": ('Layer 2', 'Layer 3', 'Layer 4', 'Layer 5') } ) gp.set_is_fault(geo_model, ['fault']) geo_model.interpolation_options.mesh_extraction = False gp.compute_model(geo_model) .. rst-class:: sphx-glr-script-out .. code-block:: none Setting Backend To: AvailableBackends.PYTORCH GPU enabled. Using device: cuda GPU device count: 1 Current GPU device: 0 .. GENERATED FROM PYTHON SOURCE LINES 63-66 .. code-block:: Python gpv.plot_2d(geo_model, cell_number=[5]) .. image-sg:: /tutorials/c_advanced/images/sphx_glr_d04_topology_001.png :alt: Cell Number: 5 Direction: y :srcset: /tutorials/c_advanced/images/sphx_glr_d04_topology_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 67-82 Analyzing Topology ^^^^^^^^^^^^^^^^^^ The gempy_plugins topology module lets us analyze the topology of a model. All we need for this is our geo_model object, the lithology block, and the fault block. We pass those into ``tp.compute_topology``, which is the starting point for several useful things: - an adjacency graph **G**, representing the topological relationships of the model - the **centroids** of all the unique topological regions in the model (x,y,z coordinates of their center) - from these, look-up tables between lithology id's and node labels (and vice versa), and adjacency queries between specific geobodies .. GENERATED FROM PYTHON SOURCE LINES 85-88 .. code-block:: Python edges, centroids = tp.compute_topology(geo_model) .. GENERATED FROM PYTHON SOURCE LINES 89-94 The first output of the topology function is the ``set`` of edges representing topology relationships between unique geobodies of the block model. An edge is represented by a ``tuple`` of two ``int`` geobody (or node) labels: .. GENERATED FROM PYTHON SOURCE LINES 96-99 .. code-block:: Python edges .. rst-class:: sphx-glr-script-out .. code-block:: none {(9, 10), (4, 10), (1, 2), (3, 4), (1, 8), (3, 10), (2, 3), (2, 9), (1, 7), (4, 5), (3, 9), (5, 10), (6, 7), (8, 9), (1, 6), (7, 8), (2, 8)} .. GENERATED FROM PYTHON SOURCE LINES 100-104 The second output is the centroids ``dict``, mapping the unique geobody id's (graph node id's) to the geobody centroid position in grid coordinates: .. GENERATED FROM PYTHON SOURCE LINES 106-109 .. code-block:: Python centroids .. rst-class:: sphx-glr-script-out .. code-block:: none {np.int64(1): array([35.27893175, 4.5 , 50.19485658]), np.int64(2): array([36.46666667, 4.5 , 29.14444444]), np.int64(3): array([37.59756098, 4.5 , 21.62195122]), np.int64(4): array([38.84563758, 4.5 , 14.00671141]), np.int64(5): array([39.09550562, 4.5 , 5.37640449]), np.int64(6): array([ 9.79081633, 4.5 , 60.10204082]), np.int64(7): array([10.17687075, 4.5 , 51.02721088]), np.int64(8): array([11.37804878, 4.5 , 43.47560976]), np.int64(9): array([12.51098901, 4.5 , 35.90659341]), np.int64(10): array([13.659857 , 4.5 , 15.34320735])} .. GENERATED FROM PYTHON SOURCE LINES 110-113 After computing the model topology, we can overlay the topology graph over a model section: .. GENERATED FROM PYTHON SOURCE LINES 116-122 Visualizing topology ~~~~~~~~~~~~~~~~~~~~ 2-D Visualization of the Topology Graph ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 125-131 .. code-block:: Python gpv.plot_topology( regular_grid=geo_model.grid.regular_grid, edges=edges, centroids=centroids ) .. image-sg:: /tutorials/c_advanced/images/sphx_glr_d04_topology_002.png :alt: d04 topology :srcset: /tutorials/c_advanced/images/sphx_glr_d04_topology_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 132-140 .. code-block:: Python plot_2d = gpv.plot_2d(geo_model, cell_number=[5], show=False) gpv.plot_topology( regular_grid=geo_model.grid.regular_grid, edges=edges, centroids=centroids, ax=plot_2d.axes[0] ) .. image-sg:: /tutorials/c_advanced/images/sphx_glr_d04_topology_003.png :alt: Cell Number: 5 Direction: y :srcset: /tutorials/c_advanced/images/sphx_glr_d04_topology_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 141-147 Adjacency Matrix ~~~~~~~~~~~~~~~~ Another way to encode and visualize the geomodel topology is using an adjacency graph: .. GENERATED FROM PYTHON SOURCE LINES 149-152 .. code-block:: Python M = tp.get_adjacency_matrix(geo_model, edges, centroids) print(M) .. rst-class:: sphx-glr-script-out .. code-block:: none [[False True False False False True True True False False] [ True False True False False False False True True False] [False True False True False False False False True True] [False False True False True False False False False True] [False False False True False False False False False True] [ True False False False False False True False False False] [ True False False False False True False True False False] [ True True False False False False True False True False] [False True True False False False False True False True] [False False True True True False False False True False]] .. GENERATED FROM PYTHON SOURCE LINES 153-156 .. code-block:: Python tp.plot_adjacency_matrix(geo_model, M) .. image-sg:: /tutorials/c_advanced/images/sphx_glr_d04_topology_004.png :alt: Topology Adjacency Matrix :srcset: /tutorials/c_advanced/images/sphx_glr_d04_topology_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 157-160 Look-up tables ~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 163-168 The ``topology`` asset provides several look-up tables to work with the unique geobody topology id's. Mapping node id's back to lithology / surface id's: .. GENERATED FROM PYTHON SOURCE LINES 170-174 .. code-block:: Python lith_lot = tp.get_lot_node_to_lith_id(geo_model, centroids) lith_lot .. rst-class:: sphx-glr-script-out .. code-block:: none {np.int64(1): np.int64(2), np.int64(2): np.int64(3), np.int64(3): np.int64(4), np.int64(4): np.int64(5), np.int64(5): np.int64(6), np.int64(6): np.int64(2), np.int64(7): np.int64(3), np.int64(8): np.int64(4), np.int64(9): np.int64(5), np.int64(10): np.int64(6)} .. GENERATED FROM PYTHON SOURCE LINES 175-177 Figuring out which nodes are in which fault block: .. GENERATED FROM PYTHON SOURCE LINES 179-183 .. code-block:: Python fault_lot = tp.get_lot_node_to_fault_block(geo_model, centroids) fault_lot .. rst-class:: sphx-glr-script-out .. code-block:: none {np.int64(1): np.int64(0), np.int64(2): np.int64(0), np.int64(3): np.int64(0), np.int64(4): np.int64(0), np.int64(5): np.int64(0), np.int64(6): np.int64(1), np.int64(7): np.int64(1), np.int64(8): np.int64(1), np.int64(9): np.int64(1), np.int64(10): np.int64(1)} .. GENERATED FROM PYTHON SOURCE LINES 184-187 We can also easily map the lithology id to the corresponding topology id's: .. GENERATED FROM PYTHON SOURCE LINES 189-192 .. code-block:: Python tp.get_lot_lith_to_node_id(lith_lot) .. rst-class:: sphx-glr-script-out .. code-block:: none {np.int64(2): [np.int64(1), np.int64(6)], np.int64(3): [np.int64(2), np.int64(7)], np.int64(4): [np.int64(3), np.int64(8)], np.int64(5): [np.int64(4), np.int64(9)], np.int64(6): [np.int64(5), np.int64(10)]} .. GENERATED FROM PYTHON SOURCE LINES 193-196 Detailed node labeling ~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 198-199 sphinx_gallery_thumbnail_number = 4 .. GENERATED FROM PYTHON SOURCE LINES 199-200 .. code-block:: Python dedges, dcentroids = tp.get_detailed_labels(geo_model, edges, centroids) .. GENERATED FROM PYTHON SOURCE LINES 201-209 .. code-block:: Python plot_2d = gpv.plot_2d(geo_model, cell_number=[5], show=False) gpv.plot_topology( regular_grid=geo_model.grid.regular_grid, edges=dedges, centroids=dcentroids, ax=plot_2d.axes[0] ) .. image-sg:: /tutorials/c_advanced/images/sphx_glr_d04_topology_005.png :alt: Cell Number: 5 Direction: y :srcset: /tutorials/c_advanced/images/sphx_glr_d04_topology_005.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 210-212 .. code-block:: Python dedges .. rst-class:: sphx-glr-script-out .. code-block:: none {('6_0', '6_1'), ('4_0', '5_0'), ('2_0', '3_0'), ('4_0', '5_1'), ('4_0', '6_1'), ('2_1', '3_1'), ('3_0', '4_1'), ('5_0', '6_0'), ('3_0', '5_1'), ('3_0', '4_0'), ('2_0', '3_1'), ('4_1', '5_1'), ('3_1', '4_1'), ('2_0', '2_1'), ('5_1', '6_1'), ('5_0', '6_1'), ('2_0', '4_1')} .. GENERATED FROM PYTHON SOURCE LINES 213-216 .. code-block:: Python dcentroids .. rst-class:: sphx-glr-script-out .. code-block:: none {'2_0': array([35.27893175, 4.5 , 50.19485658]), '3_0': array([36.46666667, 4.5 , 29.14444444]), '4_0': array([37.59756098, 4.5 , 21.62195122]), '5_0': array([38.84563758, 4.5 , 14.00671141]), '6_0': array([39.09550562, 4.5 , 5.37640449]), '2_1': array([ 9.79081633, 4.5 , 60.10204082]), '3_1': array([10.17687075, 4.5 , 51.02721088]), '4_1': array([11.37804878, 4.5 , 43.47560976]), '5_1': array([12.51098901, 4.5 , 35.90659341]), '6_1': array([13.659857 , 4.5 , 15.34320735])} .. GENERATED FROM PYTHON SOURCE LINES 217-220 Checking adjacency ~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 223-227 So let's say we want to check if the purple layer (id 5) is connected across the fault to the yellow layer (id 3). For this we can make easy use of the detailed labeling and the ``check_adjacency`` function: .. GENERATED FROM PYTHON SOURCE LINES 229-232 .. code-block:: Python tp.check_adjacency(dedges, "5_1", "3_0") .. rst-class:: sphx-glr-script-out .. code-block:: none True .. GENERATED FROM PYTHON SOURCE LINES 233-236 We can also check all geobodies that are adjacent to the purple layer (id 5) on the left side of the fault (fault id 1): .. GENERATED FROM PYTHON SOURCE LINES 238-240 .. code-block:: Python tp.get_adjacencies(dedges, "5_1") .. rst-class:: sphx-glr-script-out .. code-block:: none {'6_1', '4_1', '4_0', '3_0'} .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.768 seconds) .. _sphx_glr_download_tutorials_c_advanced_d04_topology.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: d04_topology.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: d04_topology.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: d04_topology.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_