.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "tutorials/ch1_fundamentals/ch1_1_basics.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_tutorials_ch1_fundamentals_ch1_1_basics.py: 1.1 -Basics of geological modeling with GemPy ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ .. GENERATED FROM PYTHON SOURCE LINES 7-8 Importing GemPy .. GENERATED FROM PYTHON SOURCE LINES 8-20 .. code-block:: python3 import gempy as gp # Importing auxiliary libraries import numpy as np import pandas as pd import matplotlib.pyplot as plt import os # Setting options np.random.seed(1515) pd.set_option('precision', 2) .. GENERATED FROM PYTHON SOURCE LINES 21-69 Importing and creating a set of input data """""""""""""""""""""""""""""""""""""""""" The data used for the construction of a model in GemPy is stored in Python objects. The main data classes are: :: - Surface_points - Orientations - Grid - Surfaces - Series - Additional data - Faults We will see each of this class in further detail in the future. Most of data can also be generated from raw data that comes in the form of CSV-files (CSV = comma-separated values). Such files might be attained by exporting model data from a different program such as GeoModeller or by simply creating it in spreadsheet software such as Microsoft Excel or LibreOffice Calc. In this tutorial, all input data is created by importing such CSV-files. These exemplary files can be found in the ``input_data`` folder in the root folder of GemPy. The data comprises :math:`x`-, :math:`y`- and :math:`z`-positional values for all surface points and orientation measurements. For the latter, poles, azimuth and polarity are additionally included. Surface points are furthermore assigned a formation. This might be a lithological unit such as "Sandstone" or a structural feature such as "Main Fault". It is decisive to remember that, in GemPy, interface position points mark the **bottom** of a layer. If such points are needed to resemble a top of a formation (e.g. when modeling an intrusion), this can be achieved by defining a respectively inverted orientation measurement. As we generate our ``Data`` from CSV-files, we also have to define our model's real extent in :math:`x`, :math:`y` and :math:`z`, as well as declare a desired resolution for each axis. This resolution will in turn determine the number of voxels used during modeling. Here, we rely on a medium resolution of 50x50x50, amounting to 125,000 voxels. The model extent should be chosen in a way that it contains all relevant data in a representative space. As our model voxels are not cubes, but prisms, the resolution can take a different shape than the extent. We don't recommend going much higher than 100 cells in every direction (1,000,000 voxels), as higher resolutions will become increasingly expensive to compute. .. GENERATED FROM PYTHON SOURCE LINES 71-73 .. code-block:: python3 geo_model = gp.create_model('Tutorial_ch1_1_Basics') .. GENERATED FROM PYTHON SOURCE LINES 74-83 .. code-block:: python3 data_path = 'https://raw.githubusercontent.com/cgre-aachen/gempy_data/master/' # Importing the data from CSV-files and setting extent and resolution gp.init_data(geo_model, [0, 2000., 0, 2000., 0, 750.], [50, 50, 50], path_o=data_path + "/data/input_data/getting_started/" "simple_fault_model_orientations.csv", path_i=data_path + "/data/input_data/getting_started/" "simple_fault_model_points.csv", default_values=True) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Active grids: ['regular'] Tutorial_ch1_1_Basics 2021-04-18 11:28 .. GENERATED FROM PYTHON SOURCE LINES 84-86 .. code-block:: python3 geo_model.surfaces .. raw:: html
surface series order_surfaces color id
0 Shale Default series 1 #015482 1
1 Sandstone_1 Default series 2 #9f0052 2
2 Siltstone Default series 3 #ffbe00 3
3 Sandstone_2 Default series 4 #728f02 4
4 Main_Fault Default series 5 #443988 5
5 basement Basement 1 #ff3f20 6


.. GENERATED FROM PYTHON SOURCE LINES 87-91 The input data can then be listed using the command ``get_data``. Note that the order of formations and respective allocation to series is still completely arbitrary. We will fix this in the following. .. GENERATED FROM PYTHON SOURCE LINES 94-96 .. code-block:: python3 gp.get_data(geo_model, 'surface_points').head() .. raw:: html
X Y Z smooth surface
0 1000 50 450.00 2.00e-06 Shale
1 1000 150 433.33 2.00e-06 Shale
2 1000 300 433.33 2.00e-06 Shale
3 1000 500 466.67 2.00e-06 Shale
4 1000 1000 533.33 2.00e-06 Shale


.. GENERATED FROM PYTHON SOURCE LINES 97-99 .. code-block:: python3 gp.get_data(geo_model, 'orientations').head() .. raw:: html
X Y Z G_x G_y G_z smooth surface
0 1000 1000 300 0.32 1.00e-12 0.95 0.01 Shale
1 400 1000 420 0.32 1.00e-12 0.95 0.01 Sandstone_2
2 500 1000 300 -0.95 1.00e-12 0.32 0.01 Main_Fault


.. GENERATED FROM PYTHON SOURCE LINES 100-143 Declaring the sequential order of geological formations ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - TODO update this section We want our geological units to appear in the correct order relative to age. Such order might for example be given by a depositional sequence of stratigraphy, unconformities due to erosion or other lithological genesis events such as igneous intrusions. A similar age-related order is to be declared for the faults in our model. In GemPy, the function ``set_series`` is used to assign formations to different sequential series via declaration in a Python dictionary. Defining the correct order of series is vital to the construction of the model! If you are using Python >3.6, the age-related order will already be defined by the order of key entries, i.e. the first entry is the youngest series, the last one the oldest. For older versions of Python, you will have to specify the correct order as a separate list attribute "``order_series``\ " (see cell below). You can assign several surfaces to one series. The order of the units within such as series is only relevant for the color code, thus we recommend to be consistent. You can define this order via another attribute "``order_formations``/ " or by using the specific command ``set_order_formations``. (If the order of the pile differs from the final result the color of the interfaces and input data will be different.) Every fault is treated as an independent series and have to be at set at the **top of the pile**. The relative order between the distinct faults defines the tectonic relation between them (first entry is the youngest). In a model with simple sequential stratigraphy, all layer formations can be assigned to one single series without a problem. All unit boundaries and their order would then be given by interface points. However, to model more complex lithostratigraphical relations and interactions, the definition of separate series becomes important. For example, you would need to declare a "newer" series to model an unconformity or an intrusion that disturbs older stratigraphy. By default we create a simple sequence inferred by the data: .. GENERATED FROM PYTHON SOURCE LINES 146-156 Our example model comprises four main layers (plus an underlying basement that is automatically generated by GemPy) and one main normal fault displacing those layers. Assuming a simple stratigraphy where each younger unit was deposited onto the underlying older one, we can assign these layer formations to one series called "Strat\_Series". For the fault, we declare a respective "Fault\_Series" as the first key entry in the ``set_series`` dictionary. We could give any other names to these series, the formations however have to be referred to as named in the input data. .. GENERATED FROM PYTHON SOURCE LINES 158-160 .. code-block:: python3 geo_model.surfaces .. raw:: html
surface series order_surfaces color id
0 Shale Default series 1 #015482 1
1 Sandstone_1 Default series 2 #9f0052 2
2 Siltstone Default series 3 #ffbe00 3
3 Sandstone_2 Default series 4 #728f02 4
4 Main_Fault Default series 5 #443988 5
5 basement Basement 1 #ff3f20 6


.. GENERATED FROM PYTHON SOURCE LINES 161-167 .. code-block:: python3 gp.map_stack_to_surfaces(geo_model, {"Fault_Series": 'Main_Fault', "Strat_Series": ('Sandstone_2', 'Siltstone', 'Shale', 'Sandstone_1', 'basement')}, remove_unused_series=True) .. raw:: html
surface series order_surfaces color id
4 Main_Fault Fault_Series 1 #443988 1
0 Shale Strat_Series 1 #015482 2
1 Sandstone_1 Strat_Series 2 #9f0052 3
2 Siltstone Strat_Series 3 #ffbe00 4
3 Sandstone_2 Strat_Series 4 #728f02 5
5 basement Strat_Series 5 #ff3f20 6


.. GENERATED FROM PYTHON SOURCE LINES 168-170 .. code-block:: python3 geo_model.surfaces .. raw:: html
surface series order_surfaces color id
4 Main_Fault Fault_Series 1 #443988 1
0 Shale Strat_Series 1 #015482 2
1 Sandstone_1 Strat_Series 2 #9f0052 3
2 Siltstone Strat_Series 3 #ffbe00 4
3 Sandstone_2 Strat_Series 4 #728f02 5
5 basement Strat_Series 5 #ff3f20 6


.. GENERATED FROM PYTHON SOURCE LINES 171-173 .. code-block:: python3 geo_model.stack .. raw:: html
order_series BottomRelation isActive isFault isFinite
Fault_Series 1 Erosion True False False
Strat_Series 2 Erosion True False False


.. GENERATED FROM PYTHON SOURCE LINES 174-176 .. code-block:: python3 geo_model.set_is_fault(['Fault_Series']) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Fault colors changed. If you do not like this behavior, set change_color to False. .. raw:: html
order_series BottomRelation isActive isFault isFinite
Fault_Series 1 Fault True True False
Strat_Series 2 Erosion True False False


.. GENERATED FROM PYTHON SOURCE LINES 177-179 .. code-block:: python3 geo_model.faults.faults_relations_df .. raw:: html
Fault_Series Strat_Series
Fault_Series False True
Strat_Series False False


.. GENERATED FROM PYTHON SOURCE LINES 180-182 .. code-block:: python3 geo_model.faults .. raw:: html
order_series BottomRelation isActive isFault isFinite
Fault_Series 1 Fault True True False
Strat_Series 2 Erosion True False False


.. GENERATED FROM PYTHON SOURCE LINES 183-185 .. code-block:: python3 geo_model.faults.faults_relations_df .. raw:: html
Fault_Series Strat_Series
Fault_Series False True
Strat_Series False False


.. GENERATED FROM PYTHON SOURCE LINES 186-196 Returning information from our input data ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Our model input data, here named "*geo\_model*", contains all the information that is essential for the construction of our model. You can access different types of information by using ``gp.get_data`` or simply by accessiong the atrribues. We can, for example, return the coordinates of our modeling grid via: .. GENERATED FROM PYTHON SOURCE LINES 198-200 .. code-block:: python3 geo_model.grid .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Grid Object. Values: array([[ 20. , 20. , 7.5], [ 20. , 20. , 22.5], [ 20. , 20. , 37.5], ..., [1980. , 1980. , 712.5], [1980. , 1980. , 727.5], [1980. , 1980. , 742.5]]) .. GENERATED FROM PYTHON SOURCE LINES 201-214 As mentioned before, GemPy's core algorithm is based on interpolation of two types of data: - surface\_points and - orientation measurements (if you want to know more on how this this interpolation algorithm works, checkout our paper: https://www.geosci-model-dev.net/12/1/2019/gmd-12-1-2019.pdf). We introduced the function ``get\_data`` above. You can also specify which kind of data you want to call, by declaring the string attribute "*dtype*" to be either ``'surface_points'`` (interfaces) or ``'orientations'``\ . Interfaces Dataframe: ~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 216-218 .. code-block:: python3 gp.get_data(geo_model, 'surface_points').head() .. raw:: html
X Y Z smooth surface
52 700 1000 300.0 2.00e-06 Main_Fault
53 600 1000 200.0 2.00e-06 Main_Fault
54 500 1000 100.0 2.00e-06 Main_Fault
55 1000 1000 600.0 2.00e-06 Main_Fault
56 1100 1000 700.0 2.00e-06 Main_Fault


.. GENERATED FROM PYTHON SOURCE LINES 219-222 Orientations Dataframe: ~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 224-226 .. code-block:: python3 gp.get_data(geo_model, 'orientations') .. raw:: html
X Y Z G_x G_y G_z smooth surface
2 500 1000 300 -0.95 1.00e-12 0.32 0.01 Main_Fault
0 1000 1000 300 0.32 1.00e-12 0.95 0.01 Shale
1 400 1000 420 0.32 1.00e-12 0.95 0.01 Sandstone_2


.. GENERATED FROM PYTHON SOURCE LINES 227-239 Notice that now all **surfaces** have been assigned to a **series** and are displayed in the correct order (from young to old). Visualizing input data ~~~~~~~~~~~~~~~~~~~~~~ We can also visualize our input data. This might for example be useful to check if all points and measurements are defined the way we want them to. Using the function ``plot_data``\ , we attain a 2D projection of our data points onto a plane of chosen *direction* (we can choose this attribute to be either :math:`x`, :math:`y` or :math:`z`). .. GENERATED FROM PYTHON SOURCE LINES 241-244 .. code-block:: python3 plot = gp.plot_2d(geo_model, show_lith=False, show_boundaries=False) plt.show() .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_001.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 245-258 Using ``plot_data_3D``\ , we can also visualize this data in 3D. Note that direct 3D visualization in GemPy requires `the Visualization Toolkit `__ (VTK) to be installed. All 3D plots in GemPy are interactive. This means that we can drag and drop any data point and measurement. The perpendicular axis views in VTK are particularly useful to move points solely on a desired 2D plane. Any changes will then be stored permanently in the "InputData" dataframe. If we want to reset our data points, we will then need to reload our original input data. Executing the cell below will open a new window with a 3D interactive plot of our data. .. GENERATED FROM PYTHON SOURCE LINES 261-263 .. code-block:: python3 gpv = gp.plot_3d(geo_model, image=False, plotter_type='basic') .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_002.png :alt: ch1 1 basics :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 264-277 Model generation ~~~~~~~~~~~~~~~~ Once we have made sure that we have defined all our primary information as desired in our object ``DataManagement.InputData`` (named ``geo_data`` in these tutorials), we can continue with the next step towards creating our geological model: preparing the input data for interpolation. This is done by generating an ``InterpolatorData`` object (named ``interp_data`` in these tutorials) from our ``InputData`` object via the following function: .. GENERATED FROM PYTHON SOURCE LINES 279-284 .. code-block:: python3 gp.set_interpolator(geo_model, compile_theano=True, theano_optimizer='fast_compile', ) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Setting kriging parameters to their default values. Compiling theano function... Level of Optimization: fast_compile Device: cpu Precision: float64 Number of faults: 1 Compilation Done! Kriging values: values range 2926.17 $C_o$ 203869.05 drift equations [3, 3] .. GENERATED FROM PYTHON SOURCE LINES 285-300 This function rescales the extent and coordinates of the original data (and store it in the attribute ``geo_data_res`` which behaves as a usual ``InputData`` object) and adds mathematical parameters that are needed for conducting the interpolation. The computation of this step may take a while, as it also compiles a theano function which is required for the model computation. However, should this not be needed, we can skip it by declaring ``compile_theano = False`` in the function. Furthermore, this preparation process includes an assignment of numbers to each formation. Note that GemPy's always creates a default basement formation as the last formation number. Afterwards, numbers are allocated from youngest to oldest as defined by the sequence of series and formations. On the property ``formations`` on our interpolation data, we can find out which number has been assigned to which formation: .. GENERATED FROM PYTHON SOURCE LINES 303-309 The parameters used for the interpolation can be returned using the function ``get_kriging_parameters``. These are generated automatically from the original data, but can be changed if needed. However, users should be careful doing so, if they do not fully understand their significance. .. GENERATED FROM PYTHON SOURCE LINES 311-313 .. code-block:: python3 gp.get_data(geo_model, 'kriging') .. raw:: html
values
range 2926.17
$C_o$ 203869.05
drift equations [3, 3]


.. GENERATED FROM PYTHON SOURCE LINES 314-338 At this point, we have all we need to compute our full model via ``compute_model``. By default, this will return two separate solutions in the form of arrays. The first gives information on the lithological formations, the second on the fault network in the model. These arrays consist of two subarrays as entries each: 1. Lithology block model solution: - Entry [0]: This array shows what kind of lithological formation is found in each voxel, as indicated by a respective formation\_number. - Entry [1]: Potential field array that represents the orientation of lithological units and layers in the block model. 2. Fault network block model solution: - Entry [0]: Array in which all fault-separated areas of the model are represented by a distinct number contained in each voxel. - Entry [1]: Potential field array related to the fault network in the block model. Below, we illustrate these different model solutions and how they can be used. .. GENERATED FROM PYTHON SOURCE LINES 341-343 .. code-block:: python3 sol = gp.compute_model(geo_model) .. GENERATED FROM PYTHON SOURCE LINES 344-346 .. code-block:: python3 sol .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Lithology ids [6. 6. 6. ... 2. 2. 2.] .. GENERATED FROM PYTHON SOURCE LINES 347-349 .. code-block:: python3 geo_model.solutions .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Lithology ids [6. 6. 6. ... 2. 2. 2.] .. GENERATED FROM PYTHON SOURCE LINES 350-356 Direct model visualization in GemPy ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Model solutions can be easily visualized in 2D sections in GemPy directly. Let's take a look at our lithology block: .. GENERATED FROM PYTHON SOURCE LINES 358-360 .. code-block:: python3 gp.plot_2d(geo_model, show_data=True) plt.show() .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_003.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 364-375 With ``cell_number=25`` and remembering that we defined our resolution to be 50 cells in each direction, we have chosen a section going through the middle of our block. We have moved 25 cells in ``direction='y'``, the plot thus depicts a plane parallel to the :math:`x`- and :math:`y`-axes. Setting ``plot_data=True``, we could plot original data together with the results. Changing the values for ``cell_number`` and ``direction``, we can move through our 3D block model and explore it by looking at different 2D planes. We can do the same with out lithological scalar-field solution: .. GENERATED FROM PYTHON SOURCE LINES 377-380 .. code-block:: python3 gp.plot_2d(geo_model, show_data=False, show_scalar=True, show_lith=False) plt.show() .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_004.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 381-384 .. code-block:: python3 gp.plot_2d(geo_model, series_n=1, show_data=False, show_scalar=True, show_lith=False) plt.show() .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_005.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 385-390 This illustrates well the fold-related deformation of the stratigraphy, as well as the way the layers are influenced by the fault. The fault network modeling solutions can be visualized in the same way: .. GENERATED FROM PYTHON SOURCE LINES 392-394 .. code-block:: python3 geo_model.solutions.scalar_field_at_surface_points .. rst-class:: sphx-glr-script-out Out: .. code-block:: none array([[0.03075848, 0. , 0. , 0. , 0. ], [0. , 0.77174354, 0.72471042, 0.80357372, 0.83598092]]) .. GENERATED FROM PYTHON SOURCE LINES 395-398 .. code-block:: python3 gp.plot_2d(geo_model, show_block=True, show_lith=False) plt.show() .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_006.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 399-402 .. code-block:: python3 gp.plot_2d(geo_model, series_n=1, show_block=True, show_lith=False) plt.show() .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_007.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 403-413 Marching cubes and vtk visualization ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ In addition to 2D sections we can extract surfaces to visualize in 3D renderers. Surfaces can be visualized as 3D triangle complexes in VTK (see function plot\_surfaces\_3D below). To create these triangles, we need to extract respective vertices and simplices from the potential fields of lithologies and faults. This process is automatized in GemPy with the function ``get_surface``\ . .. GENERATED FROM PYTHON SOURCE LINES 415-418 .. code-block:: python3 ver, sim = gp.get_surfaces(geo_model) gpv = gp.plot_3d(geo_model, image=False, plotter_type='basic') .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_008.png :alt: ch1 1 basics :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 419-426 Using the rescaled interpolation data, we can also run our 3D VTK visualization in an interactive mode which allows us to alter and update our model in real time. Similarly to the interactive 3D visualization of our input data, the changes are permanently saved (in the ``InterpolationInput.dataframe`` object). Additionally, the resulting changes in the geological models are re-computed in real time. .. GENERATED FROM PYTHON SOURCE LINES 429-431 Adding topography ~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 431-433 .. code-block:: python3 geo_model.set_topography(d_z=(350, 750)) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Active grids: ['regular' 'topography'] Grid Object. Values: array([[ 20. , 20. , 7.5 ], [ 20. , 20. , 22.5 ], [ 20. , 20. , 37.5 ], ..., [2000. , 1918.36734694, 423.48951452], [2000. , 1959.18367347, 430.25455308], [2000. , 2000. , 431.07163663]]) .. GENERATED FROM PYTHON SOURCE LINES 434-444 .. code-block:: python3 gp.compute_model(geo_model) gp.plot_2d(geo_model, show_topography=True) plt.show() # sphinx_gallery_thumbnail_number = 9 gpv = gp.plot_3d(geo_model, plotter_type='basic', show_topography=True, show_surfaces=True, show_lith=True, image=False) .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_009.png :alt: ch1 1 basics :class: sphx-glr-single-img .. image:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_010.png :alt: Cell Number: mid Direction: y :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out Out: .. code-block:: none /WorkSSD/PythonProjects/gempy/gempy/core/solution.py:173: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray. self.geological_map = np.array( .. GENERATED FROM PYTHON SOURCE LINES 445-452 Compute at a given location ~~~~~~~~~~~~~~~~~~~~~~~~~~~ This is done by modifying the grid to a custom grid and recomputing. Notice that the results are given as *grid + surfaces\_points\_ref + surface\_points\_rest locations* .. GENERATED FROM PYTHON SOURCE LINES 454-457 .. code-block:: python3 x_i = np.array([[3, 5, 6]]) sol = gp.compute_model(geo_model, at=x_i) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Active grids: ['custom'] /WorkSSD/PythonProjects/gempy/gempy/core/solution.py:168: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray. self.custom = np.array( .. GENERATED FROM PYTHON SOURCE LINES 458-459 Therefore if we just want the value at **x\_i**: .. GENERATED FROM PYTHON SOURCE LINES 461-463 .. code-block:: python3 sol.custom .. rst-class:: sphx-glr-script-out Out: .. code-block:: none array([array([[6.]]), array([[0.18630133], [0.63163565]])], dtype=object) .. GENERATED FROM PYTHON SOURCE LINES 464-465 This return the id, and the scalar field values for each series .. GENERATED FROM PYTHON SOURCE LINES 467-470 Save the model ~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 472-482 GemPy uses Python [pickle] for fast storing temporary objects (https://docs.python.org/3/library/pickle.html). However, module version consistency is required. For loading a pickle into GemPy, you have to make sure that you are using the same version of pickle and dependent modules (e.g.: ``Pandas``, ``NumPy``) as were used when the data was originally stored. For long term-safer storage we can export the ``pandas.DataFrames`` to csv by using: .. GENERATED FROM PYTHON SOURCE LINES 484-485 .. code-block:: python3 gp.save_model(geo_model) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none True .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 0 minutes 11.939 seconds) .. _sphx_glr_download_tutorials_ch1_fundamentals_ch1_1_basics.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: ch1_1_basics.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: ch1_1_basics.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_