.. 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 :ref:`Go to the end ` 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 6-8 .. code-block:: Python import numpy as np .. GENERATED FROM PYTHON SOURCE LINES 9-11 Importing Necessary Libraries """""""""""""""""""""""""""""" .. GENERATED FROM PYTHON SOURCE LINES 11-14 .. code-block:: Python import gempy as gp import gempy_viewer as gpv .. GENERATED FROM PYTHON SOURCE LINES 15-48 Importing and Defining Input Data """"""""""""""""""""""""""""""""" :obj:`gempy.core.data.GeoModel` GemPy uses Python objects to store the data that builds the geological model. The main data classes include: - :obj:`gempy.core.data.GeoModel` - :obj:`gempy.core.data.StructuralFrame` - :obj:`gempy.core.data.StructuralGroup` - :obj:`gempy.core.data.StructuralElement` - :obj:`gempy.core.data.SurfacePointsTable` - :obj:`gempy.core.data.OrientationsTable` - :obj:`gempy.core.data.Grid` You can also create data from raw CSV files (comma-separated values). This could be useful if you are exporting model data from a different program or creating it in a spreadsheet software like Microsoft Excel or LibreOffice Calc. In this tutorial, we'll use CSV files to generate input data. You can find these example files in the `gempy data` repository on GitHub. The data consists of x, y, and z positional values for all surface points and orientation measurements. Additional data includes poles, azimuth and polarity (or the gradient components). Surface points are assigned a formation, which can be a lithological unit (like "Sandstone") or a structural feature (like "Main Fault"). It's important to note that, in GemPy, interface position points mark the **bottom** of a layer. If you need points to represent the top of a formation (for example, when modeling an intrusion), you can define an inverted orientation measurement. While generating data from CSV files, we also need to define the model's real extent in x, y, and z. This extent defines the area used for interpolation and many of the plotting functions. We also set a resolution to establish a regular grid right away. This resolution will dictate the number of voxels used during modeling. We're using a medium resolution of 50x50x50 here, which results in 125,000 voxels. The model extent should enclose all relevant data in a representative space. As our model voxels are prisms rather than cubes, the resolution can differ from the extent. However, it is recommended to avoid going beyond 100 cells in each direction (1,000,000 voxels) to prevent excessive computational costs. .. GENERATED FROM PYTHON SOURCE LINES 50-64 .. code-block:: Python data_path = 'https://raw.githubusercontent.com/cgre-aachen/gempy_data/master/' geo_model: gp.data.GeoModel = gp.create_geomodel( project_name='Tutorial_ch1_1_Basics', extent=[0, 2000, 0, 2000, 0, 750], refinement=6, # * Here we define the number of octree levels. If octree levels are defined, the resolution is ignored. importer_helper=gp.data.ImporterHelper( path_to_orientations=data_path + "/data/input_data/getting_started/simple_fault_model_orientations.csv", path_to_surface_points=data_path + "/data/input_data/getting_started/simple_fault_model_points.csv", hash_surface_points="4cdd54cd510cf345a583610585f2206a2936a05faaae05595b61febfc0191563", hash_orientations="7ba1de060fc8df668d411d0207a326bc94a6cdca9f5fe2ed511fd4db6b3f3526" ) ) .. rst-class:: sphx-glr-script-out .. code-block:: none Surface points hash: 4cdd54cd510cf345a583610585f2206a2936a05faaae05595b61febfc0191563 Orientations hash: 7ba1de060fc8df668d411d0207a326bc94a6cdca9f5fe2ed511fd4db6b3f3526 .. GENERATED FROM PYTHON SOURCE LINES 65-78 .. admonition:: New in GemPy 3! GemPy 3 has introduced the ``ImporterHelper`` class to streamline importing data from various sources. This class simplifies the process of passing multiple arguments needed for importing data and will likely see further extensions in the future. Currently, one of its uses is to handle `pooch` arguments for downloading data from the internet. Reviewing the Imported Data """"""""""""""""""""""""""" Now that the `geo_model` object is set up and the data imported from the CSV files, we review the data imported using the properties `surface_points_copy` and `orientations_copy`. Using `structural_frame.element_id_name_map`, we can see which ID corresponds to which structural element name in the data. .. GENERATED FROM PYTHON SOURCE LINES 81-83 .. code-block:: Python geo_model.surface_points_copy .. raw:: html
XYZidnugget
700.001000.00300.00497450340.00
600.001000.00200.00497450340.00
500.001000.00100.00497450340.00
1000.001000.00600.00497450340.00
1100.001000.00700.00497450340.00
1000.0050.00350.001829037370.00
1000.00150.00333.331829037370.00
1000.00300.00333.331829037370.00
1000.00500.00366.671829037370.00
1000.001000.00433.331829037370.00
...............
1100.001700.00473.334262594380.00
1100.001950.00490.004262594380.00
0.001000.00433.334262594380.00
300.001000.00400.004262594380.00
600.001000.00366.674262594380.00
1300.001000.00566.674262594380.00
1600.001000.00550.004262594380.00
1900.001000.00566.674262594380.00
1700.00500.00533.334262594380.00
1700.001500.00516.674262594380.00


.. GENERATED FROM PYTHON SOURCE LINES 84-87 .. code-block:: Python geo_model.orientations_copy .. raw:: html
XYZG_xG_yG_zidnugget
500.001000.00300.00-0.95-0.000.32497450340.01
400.001000.00420.000.320.000.952470168680.01
1000.001000.00300.000.320.000.953323137030.01


.. GENERATED FROM PYTHON SOURCE LINES 88-117 Declaring the Sequential Order of Structural Elements (Geological Formations) """"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" In our model, we want the geological units to appear in the correct chronological order. This order could be determined by a sequence of stratigraphic deposition, unconformities due to erosion, or other lithological genesis events like igneous intrusions. A similar age-related order is declared for faults in our model. In GemPy, we use the function `gempy.map_stack_to_surfaces` to assign formations or faults to different sequential series by declaring them in a Python dictionary. The correct ordering of series is crucial for model construction! It's possible to assign several surfaces to one series. The order of units within a series only affects the color code, so we recommend maintaining consistency. The order can be defined by simply changing the order of the lists within the `gempy.core.data.StructuralFrame.structural_groups` and `gempy.core.data.StructuralGroups.elements` attributes. Faults are treated as independent groups and must be younger than the groups they affect. The relative order between different faults defines their tectonic relationship (the first entry is the youngest). For a model with simple sequential stratigraphy, all layer formations can be assigned to one series without an issue. All unit boundaries and their order would then be determined by interface points. However, to model more complex lithostratigraphical relations and interactions, separate series definitions become important. For example, modeling an unconformity or an intrusion that disrupts older stratigraphy would require declaring a younger series. By default, a simple sequence/group is created inferred from the data as shown above. 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 structural group called "Strat_Series". For the fault, we declare a respective "Fault_Series" as the first key entry in the mapping dictionary. We could give any other names to these series; the formations, however, have to be referred to as named in the input data. In the following, we map the "Main Fault" to the "Fault Series" and the individual formations to the "Strat Series". .. GENERATED FROM PYTHON SOURCE LINES 119-128 .. code-block:: Python gp.map_stack_to_surfaces( gempy_model=geo_model, mapping_object= # TODO: This mapping I do not like it too much. We should be able to do it passing the data objects directly { "Fault_Series": 'Main_Fault', "Strat_Series": ('Sandstone_2', 'Siltstone', 'Shale', 'Sandstone_1') } ) .. raw:: html
Structural Groups: StructuralGroup:
Name:Fault_Series
Structural Relation:StackRelationType.ERODE
Elements:
StructuralElement:
Name:Main_Fault

StructuralGroup:
Name:Strat_Series
Structural Relation:StackRelationType.ERODE
Elements:
StructuralElement:
Name:Sandstone_2

StructuralElement:
Name:Siltstone

StructuralElement:
Name:Shale

StructuralElement:
Name:Sandstone_1
Fault Relations:
Fault_Seri...Strat_Seri...
Fault_Series
Strat_Series
True
False


.. GENERATED FROM PYTHON SOURCE LINES 129-131 Note how the structural frame still indicates the "Fault Series" group to have a relation type "erode". We still need to tell GemPy that we want this group to be a fault. We do this using the function `set_is_fault`. .. GENERATED FROM PYTHON SOURCE LINES 133-138 .. code-block:: Python gp.set_is_fault( frame=geo_model.structural_frame, fault_groups=['Fault_Series'] ) .. raw:: html
Structural Groups: StructuralGroup:
Name:Fault_Series
Structural Relation:StackRelationType.FAULT
Elements:
StructuralElement:
Name:Main_Fault

StructuralGroup:
Name:Strat_Series
Structural Relation:StackRelationType.ERODE
Elements:
StructuralElement:
Name:Sandstone_2

StructuralElement:
Name:Siltstone

StructuralElement:
Name:Shale

StructuralElement:
Name:Sandstone_1
Fault Relations:
Fault_Seri...Strat_Seri...
Fault_Series
Strat_Series
True
False


.. GENERATED FROM PYTHON SOURCE LINES 139-147 Now, all surfaces have been assigned to a series and are displayed in the correct order (from young to old). Returning Information from Our Input Data ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Our model input data, named "geo_model", contains all essential information for constructing our model. You can access different types of information by accessing the attributes. For instance, you can retrieve the coordinates of our modeling grid as follows: .. GENERATED FROM PYTHON SOURCE LINES 150-152 .. code-block:: Python geo_model.grid .. rst-class:: sphx-glr-script-out .. code-block:: none Grid(values=array([[ 15.625 , 15.625 , 5.859375], [ 15.625 , 15.625 , 17.578125], [ 15.625 , 15.625 , 29.296875], ..., [1984.375 , 1984.375 , 720.703125], [1984.375 , 1984.375 , 732.421875], [1984.375 , 1984.375 , 744.140625]], shape=(262144, 3)), length=array([], dtype=float64), _octree_grid=RegularGrid(resolution=array([64, 64, 64]), extent=array([ 0., 2000., 0., 2000., 0., 750.]), values=array([[ 15.625 , 15.625 , 5.859375], [ 15.625 , 15.625 , 17.578125], [ 15.625 , 15.625 , 29.296875], ..., [1984.375 , 1984.375 , 720.703125], [1984.375 , 1984.375 , 732.421875], [1984.375 , 1984.375 , 744.140625]], shape=(262144, 3)), mask_topo=array([], shape=(0, 3), dtype=bool), _transform=None), _dense_grid=None, _custom_grid=None, _topography=None, _sections=None, _centered_grid=None, _transform=None, _octree_levels=-1) .. GENERATED FROM PYTHON SOURCE LINES 153-160 Visualizing input data ~~~~~~~~~~~~~~~~~~~~~~ Since we have already imported our input data, we can go ahead and visualize it in 2D and 3D. This can be useful to check if all points and orientations are defined the way we want them to be. Using the function `gempy_viewer.plot2d`, we obtain a 2D projection of our data points onto a plane of the chosen *direction* (we can choose this attribute to be either x, y or z, with the default being y). Below, we can freely switch the direction and check out the projection from different sides to get a good idea. .. GENERATED FROM PYTHON SOURCE LINES 164-166 .. code-block:: Python plot = gpv.plot_2d(geo_model, show_lith=False, show_boundaries=False) .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_001.png :alt: Cell Number: mid Direction: y :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 167-169 Beyond 2D, however, we can also visualize our input data in full 3D using `gempy_viewer.plot_3d`. Note that direct 3D visualization in GemPy requires [the Visualization Toolkit](https://www.vtk.org/) (VTK) to be installed. .. GENERATED FROM PYTHON SOURCE LINES 172-174 .. code-block:: Python gpv.plot_3d(geo_model, image=False, plotter_type='basic') .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_002.png :alt: ch1 1 basics :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 175-188 Model Generation ~~~~~~~~~~~~~~~~ Once we've correctly defined all our primary information in our `gempy.core.data.GeoModel` object (referred to as `geo_model` in these tutorials), we can proceed with the model computation step. We can go ahead and save the solution of a specific computation as we do below, but solutions are also stored within the `gempy.core.data.GeoModel` object for future reference. .. admonition:: New in GemPy 3! Numpy and TensorFlow backend Unlike previous versions, GemPy 3 doesn't rely on `theano` or `asera`. Instead, it utilizes `numpy` or `tensorflow`. Consequently, we no longer need to recompile all theano functions (TensorFlow uses eager execution; we found no notable speed difference after profiling the XLA compiler). .. GENERATED FROM PYTHON SOURCE LINES 190-194 The parameters used for the interpolation are stored in `gempy.core.data.GeoModel.interpolation_options`. These parameters have sensible default values that you can modify if necessary. However, we advise caution when changing these parameters unless you fully understand their implications. .. GENERATED FROM PYTHON SOURCE LINES 196-197 Display the current interpolation options .. GENERATED FROM PYTHON SOURCE LINES 197-199 .. code-block:: Python geo_model.interpolation_options .. rst-class:: sphx-glr-script-out .. code-block:: none InterpolationOptions(kernel_options=KernelOptions(range=1.7, c_o=10.0, uni_degree=1, i_res=4.0, gi_res=2.0, number_dimensions=3, kernel_function=AvailableKernelFunctions.cubic, kernel_solver=Solvers.DEFAULT, compute_condition_number=False, optimizing_condition_number=False, condition_number=None), evaluation_options=EvaluationOptions(_number_octree_levels=6, _number_octree_levels_surface=4, octree_curvature_threshold=-1.0, octree_error_threshold=1.0, octree_min_level=2, mesh_extraction=True, mesh_extraction_masking_options=, mesh_extraction_fancy=True, evaluation_chunk_size=500000, compute_scalar_gradient=False, verbose=False), debug=True, cache_mode=, cache_model_name='Tutorial_ch1_1_Basics', block_solutions_type=, sigmoid_slope=5000000, debug_water_tight=False, temp_interpolation_values=TempInterpolationValues(current_octree_level=0)) .. GENERATED FROM PYTHON SOURCE LINES 200-204 With all our prerequisites in place, we can now compute our complete geological model using :func:`gempy.compute_model`. This function returns a :obj:`gempy.core.data.Solutions` object. The following sections illustrate these different model solutions and how to utilize them. .. GENERATED FROM PYTHON SOURCE LINES 206-207 Compute the geological model and get the solutions .. GENERATED FROM PYTHON SOURCE LINES 207-210 .. code-block:: Python sol = gp.compute_model(geo_model) sol .. rst-class:: sphx-glr-script-out .. code-block:: none Setting Backend To: AvailableBackends.numpy /home/leguark/TeamCity/work/3a8738c25f60c3c9/venv/lib/python3.10/site-packages/gempy_engine/modules/activator/_soft_segment.py:95: RuntimeWarning: overflow encountered in exp return 1.0 / (1.0 + bt.t.exp(x)) Chunking done: 20 chunks Chunking done: 20 chunks /home/leguark/TeamCity/work/3a8738c25f60c3c9/venv/lib/python3.10/site-packages/gempy_engine/modules/activator/_soft_segment.py:95: RuntimeWarning: overflow encountered in exp return 1.0 / (1.0 + bt.t.exp(x)) .. raw:: html
Solutions: 6 Octree Levels, 5 DualContouringMeshes


.. GENERATED FROM PYTHON SOURCE LINES 211-212 Solutions are also stored within the :obj:`gempy.core.data.GeoModel` object, for future reference. .. GENERATED FROM PYTHON SOURCE LINES 214-216 .. code-block:: Python geo_model.solutions .. raw:: html
Solutions: 6 Octree Levels, 5 DualContouringMeshes


.. GENERATED FROM PYTHON SOURCE LINES 217-223 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 225-227 .. code-block:: Python gpv.plot_2d(geo_model, show_data=True, cell_number="mid", direction='y') .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_003.png :alt: Cell Number: mid Direction: y :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 228-238 With ``cell_number=mid``, we have chosen a section going through the middle of our block. We have moved in ``direction='y'``, the plot thus depicts a plane parallel to the :math:`x`- and :math:`y`-axes. Setting ``show_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 the underlying scalar-field solution: .. GENERATED FROM PYTHON SOURCE LINES 241-249 .. code-block:: Python gpv.plot_2d( model=geo_model, series_n=0, # This will plot the scalar field used for the fault show_data=False, show_scalar=True, show_lith=False ) .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_004.png :alt: Cell Number: mid Direction: y :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 250-258 .. code-block:: Python gpv.plot_2d( model=geo_model, series_n=1, # This will plot the scalar field used for the stratigraphy show_data=False, show_scalar=True, show_lith=False ) .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_005.png :alt: Cell Number: mid Direction: y :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_005.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 259-274 Dual Contouring 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 using dual contouring in the :obj:`gempy_engine`. .. admonition:: New in GemPy 3! Dual Contouring GemPy 3 uses dual contouring to extract surfaces from the scalar fields. The method is completely coded in :obj:`gempy_engine` what also enables further improvements in the midterm. This method is more efficient to use together with octrees and suited better the new capabilities of gempy3. .. GENERATED FROM PYTHON SOURCE LINES 276-278 .. code-block:: Python gpv.plot_3d(geo_model, show_data=False, image=False, plotter_type='basic') .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_006.png :alt: ch1 1 basics :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 279-284 Adding topography ~~~~~~~~~~~~~~~~~~~~~~~~~~~ In gempy we can add more grid types for different purposes. We will explore this concept in more detail in the next tutorials (:doc:`ch1_3a_grids`). For now, we will just add a topography grid to our model. This grid allows us to intersect the surfaces as well as compute a high resolution geological mal. .. GENERATED FROM PYTHON SOURCE LINES 287-306 .. code-block:: Python gp.set_topography_from_random( grid=geo_model.grid, fractal_dimension=1.2, d_z=np.array([350, 750]), topography_resolution=np.array([50, 50]), ) gp.compute_model(geo_model) gpv.plot_2d(geo_model, show_topography=True) gpv.plot_3d( model=geo_model, plotter_type='basic', show_topography=True, show_surfaces=True, show_lith=True, image=False ) .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_007.png :alt: ch1 1 basics :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_007.png :class: sphx-glr-single-img .. image-sg:: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_008.png :alt: Cell Number: mid Direction: y :srcset: /tutorials/ch1_fundamentals/images/sphx_glr_ch1_1_basics_008.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Active grids: GridTypes.NONE|TOPOGRAPHY|OCTREE Setting Backend To: AvailableBackends.numpy /home/leguark/TeamCity/work/3a8738c25f60c3c9/venv/lib/python3.10/site-packages/gempy_engine/modules/activator/_soft_segment.py:95: RuntimeWarning: overflow encountered in exp return 1.0 / (1.0 + bt.t.exp(x)) Chunking done: 20 chunks Chunking done: 20 chunks /home/leguark/TeamCity/work/3a8738c25f60c3c9/venv/lib/python3.10/site-packages/gempy_engine/modules/activator/_soft_segment.py:95: RuntimeWarning: overflow encountered in exp return 1.0 / (1.0 + bt.t.exp(x)) .. GENERATED FROM PYTHON SOURCE LINES 307-312 Compute at a given location ~~~~~~~~~~~~~~~~~~~~~~~~~~~ This is done by modifying the grid to a custom grid and recomputing. .. GENERATED FROM PYTHON SOURCE LINES 314-321 .. code-block:: Python x_i = np.array([[1000, 1000, 1000]]) lith_values_at_coords: np.ndarray = gp.compute_model_at( gempy_model=geo_model, at=x_i ) lith_values_at_coords .. rst-class:: sphx-glr-script-out .. code-block:: none Active grids: GridTypes.NONE|TOPOGRAPHY|CUSTOM|OCTREE Setting Backend To: AvailableBackends.numpy /home/leguark/TeamCity/work/3a8738c25f60c3c9/venv/lib/python3.10/site-packages/gempy_engine/modules/activator/_soft_segment.py:95: RuntimeWarning: overflow encountered in exp return 1.0 / (1.0 + bt.t.exp(x)) Chunking done: 20 chunks Chunking done: 20 chunks /home/leguark/TeamCity/work/3a8738c25f60c3c9/venv/lib/python3.10/site-packages/gempy_engine/modules/activator/_soft_segment.py:95: RuntimeWarning: overflow encountered in exp return 1.0 / (1.0 + bt.t.exp(x)) array([2.]) .. GENERATED FROM PYTHON SOURCE LINES 322-323 Therefore if we just want the value at **x\_i**: .. GENERATED FROM PYTHON SOURCE LINES 325-327 .. code-block:: Python geo_model.solutions.raw_arrays.custom .. rst-class:: sphx-glr-script-out .. code-block:: none array([2.]) .. GENERATED FROM PYTHON SOURCE LINES 328-335 .. admonition:: Work in progress GemPy3 model serialization is currently being redisigned. Therefore, at the current version, there is not a build in method to save the model. However, since now the data model should be completely robust, you should be able to save the :obj:`gempy.core.data.GeoModel` and all its attributes using the standard python library [pickle](https://docs.python.org/3/library/pickle.html) sphinx_gallery_thumbnail_number = -2 .. rst-class:: sphx-glr-timing **Total running time of the script:** (1 minutes 25.098 seconds) .. _sphx_glr_download_tutorials_ch1_fundamentals_ch1_1_basics.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: ch1_1_basics.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: ch1_1_basics.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: ch1_1_basics.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_