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

.. only:: html

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

        Click :ref:`here <sphx_glr_download_integrations_gempy_export_MOOSE.py>`
        to download the full example code

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

.. _sphx_glr_integrations_gempy_export_MOOSE.py:


Export a geological model from GemPy to use in MOOSE
====================================================

.. GENERATED FROM PYTHON SOURCE LINES 8-10

.. code-block:: python3

    import gempy as gp








.. GENERATED FROM PYTHON SOURCE LINES 11-19

Creating a geological model
---------------------------

The procedure of generating a geological model is presented in detail in
`Chapter
1-1 <https://nbviewer.jupyter.org/github/cgre-aachen/gempy/blob/master/notebooks/tutorials/ch1-1_Basics.ipynb>`__
of the GemPy tutorials, so it will only be briefly presented here


.. GENERATED FROM PYTHON SOURCE LINES 21-22

Initiate a model

.. GENERATED FROM PYTHON SOURCE LINES 22-25

.. code-block:: python3

    geo_model = gp.create_model('tutorial_moose_exp')
    data_path = 'https://raw.githubusercontent.com/cgre-aachen/gempy_data/master/'








.. GENERATED FROM PYTHON SOURCE LINES 26-27

Import data from CSV-files with setting the resolution and model extent

.. GENERATED FROM PYTHON SOURCE LINES 27-32

.. code-block:: python3

    gp.init_data(geo_model, [0, 2000., 0, 2000., 0, 2000.], [50, 50, 80],
                 path_o=data_path + "/data/input_data/tut_chapter1/simple_fault_model_orientations.csv",
                 path_i=data_path + "/data/input_data/tut_chapter1/simple_fault_model_points.csv",
                 default_values=True)





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

 Out:

 .. code-block:: none

    Active grids: ['regular']

    tutorial_moose_exp  2021-04-18 11:39



.. GENERATED FROM PYTHON SOURCE LINES 33-34

present the units and series

.. GENERATED FROM PYTHON SOURCE LINES 34-36

.. code-block:: python3

    geo_model.surfaces






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <style  type="text/css" >
    #T_7b74d_row0_col3{
                background-color:  #015482;
            }#T_7b74d_row1_col3{
                background-color:  #9f0052;
            }#T_7b74d_row2_col3{
                background-color:  #ffbe00;
            }#T_7b74d_row3_col3{
                background-color:  #728f02;
            }#T_7b74d_row4_col3{
                background-color:  #443988;
            }#T_7b74d_row5_col3{
                background-color:  #ff3f20;
            }</style><table id="T_7b74d_" ><thead>    <tr>        <th class="blank level0" ></th>        <th class="col_heading level0 col0" >surface</th>        <th class="col_heading level0 col1" >series</th>        <th class="col_heading level0 col2" >order_surfaces</th>        <th class="col_heading level0 col3" >color</th>        <th class="col_heading level0 col4" >id</th>    </tr></thead><tbody>
                    <tr>
                            <th id="T_7b74d_level0_row0" class="row_heading level0 row0" >0</th>
                            <td id="T_7b74d_row0_col0" class="data row0 col0" >Shale</td>
                            <td id="T_7b74d_row0_col1" class="data row0 col1" >Default series</td>
                            <td id="T_7b74d_row0_col2" class="data row0 col2" >1</td>
                            <td id="T_7b74d_row0_col3" class="data row0 col3" >#015482</td>
                            <td id="T_7b74d_row0_col4" class="data row0 col4" >1</td>
                </tr>
                <tr>
                            <th id="T_7b74d_level0_row1" class="row_heading level0 row1" >1</th>
                            <td id="T_7b74d_row1_col0" class="data row1 col0" >Sandstone_1</td>
                            <td id="T_7b74d_row1_col1" class="data row1 col1" >Default series</td>
                            <td id="T_7b74d_row1_col2" class="data row1 col2" >2</td>
                            <td id="T_7b74d_row1_col3" class="data row1 col3" >#9f0052</td>
                            <td id="T_7b74d_row1_col4" class="data row1 col4" >2</td>
                </tr>
                <tr>
                            <th id="T_7b74d_level0_row2" class="row_heading level0 row2" >2</th>
                            <td id="T_7b74d_row2_col0" class="data row2 col0" >Siltstone</td>
                            <td id="T_7b74d_row2_col1" class="data row2 col1" >Default series</td>
                            <td id="T_7b74d_row2_col2" class="data row2 col2" >3</td>
                            <td id="T_7b74d_row2_col3" class="data row2 col3" >#ffbe00</td>
                            <td id="T_7b74d_row2_col4" class="data row2 col4" >3</td>
                </tr>
                <tr>
                            <th id="T_7b74d_level0_row3" class="row_heading level0 row3" >3</th>
                            <td id="T_7b74d_row3_col0" class="data row3 col0" >Sandstone_2</td>
                            <td id="T_7b74d_row3_col1" class="data row3 col1" >Default series</td>
                            <td id="T_7b74d_row3_col2" class="data row3 col2" >4</td>
                            <td id="T_7b74d_row3_col3" class="data row3 col3" >#728f02</td>
                            <td id="T_7b74d_row3_col4" class="data row3 col4" >4</td>
                </tr>
                <tr>
                            <th id="T_7b74d_level0_row4" class="row_heading level0 row4" >4</th>
                            <td id="T_7b74d_row4_col0" class="data row4 col0" >Main_Fault</td>
                            <td id="T_7b74d_row4_col1" class="data row4 col1" >Default series</td>
                            <td id="T_7b74d_row4_col2" class="data row4 col2" >5</td>
                            <td id="T_7b74d_row4_col3" class="data row4 col3" >#443988</td>
                            <td id="T_7b74d_row4_col4" class="data row4 col4" >5</td>
                </tr>
                <tr>
                            <th id="T_7b74d_level0_row5" class="row_heading level0 row5" >5</th>
                            <td id="T_7b74d_row5_col0" class="data row5 col0" >basement</td>
                            <td id="T_7b74d_row5_col1" class="data row5 col1" >Basement</td>
                            <td id="T_7b74d_row5_col2" class="data row5 col2" >1</td>
                            <td id="T_7b74d_row5_col3" class="data row5 col3" >#ff3f20</td>
                            <td id="T_7b74d_row5_col4" class="data row5 col4" >6</td>
                </tr>
        </tbody></table>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 37-38

combine units in series and make two series, as the fault needs its own

.. GENERATED FROM PYTHON SOURCE LINES 38-46

.. 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)

    # set the fault series to be fault object
    geo_model.set_is_fault(['Fault_Series'], change_color=False)






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <table border="1" class="dataframe">
      <thead>
        <tr style="text-align: right;">
          <th></th>
          <th>order_series</th>
          <th>BottomRelation</th>
          <th>isActive</th>
          <th>isFault</th>
          <th>isFinite</th>
        </tr>
      </thead>
      <tbody>
        <tr>
          <th>Fault_Series</th>
          <td>1</td>
          <td>Fault</td>
          <td>True</td>
          <td>True</td>
          <td>False</td>
        </tr>
        <tr>
          <th>Strat_Series</th>
          <td>2</td>
          <td>Erosion</td>
          <td>True</td>
          <td>False</td>
          <td>False</td>
        </tr>
      </tbody>
    </table>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 47-48

check whether series were assigned correctly

.. GENERATED FROM PYTHON SOURCE LINES 48-50

.. code-block:: python3

    geo_model.surfaces






.. raw:: html

    <div class="output_subarea output_html rendered_html output_result">
    <style  type="text/css" >
    #T_9fd1c_row0_col3{
                background-color:  #443988;
            }#T_9fd1c_row1_col3{
                background-color:  #015482;
            }#T_9fd1c_row2_col3{
                background-color:  #9f0052;
            }#T_9fd1c_row3_col3{
                background-color:  #ffbe00;
            }#T_9fd1c_row4_col3{
                background-color:  #728f02;
            }#T_9fd1c_row5_col3{
                background-color:  #ff3f20;
            }</style><table id="T_9fd1c_" ><thead>    <tr>        <th class="blank level0" ></th>        <th class="col_heading level0 col0" >surface</th>        <th class="col_heading level0 col1" >series</th>        <th class="col_heading level0 col2" >order_surfaces</th>        <th class="col_heading level0 col3" >color</th>        <th class="col_heading level0 col4" >id</th>    </tr></thead><tbody>
                    <tr>
                            <th id="T_9fd1c_level0_row0" class="row_heading level0 row0" >4</th>
                            <td id="T_9fd1c_row0_col0" class="data row0 col0" >Main_Fault</td>
                            <td id="T_9fd1c_row0_col1" class="data row0 col1" >Fault_Series</td>
                            <td id="T_9fd1c_row0_col2" class="data row0 col2" >1</td>
                            <td id="T_9fd1c_row0_col3" class="data row0 col3" >#443988</td>
                            <td id="T_9fd1c_row0_col4" class="data row0 col4" >1</td>
                </tr>
                <tr>
                            <th id="T_9fd1c_level0_row1" class="row_heading level0 row1" >0</th>
                            <td id="T_9fd1c_row1_col0" class="data row1 col0" >Shale</td>
                            <td id="T_9fd1c_row1_col1" class="data row1 col1" >Strat_Series</td>
                            <td id="T_9fd1c_row1_col2" class="data row1 col2" >1</td>
                            <td id="T_9fd1c_row1_col3" class="data row1 col3" >#015482</td>
                            <td id="T_9fd1c_row1_col4" class="data row1 col4" >2</td>
                </tr>
                <tr>
                            <th id="T_9fd1c_level0_row2" class="row_heading level0 row2" >1</th>
                            <td id="T_9fd1c_row2_col0" class="data row2 col0" >Sandstone_1</td>
                            <td id="T_9fd1c_row2_col1" class="data row2 col1" >Strat_Series</td>
                            <td id="T_9fd1c_row2_col2" class="data row2 col2" >2</td>
                            <td id="T_9fd1c_row2_col3" class="data row2 col3" >#9f0052</td>
                            <td id="T_9fd1c_row2_col4" class="data row2 col4" >3</td>
                </tr>
                <tr>
                            <th id="T_9fd1c_level0_row3" class="row_heading level0 row3" >2</th>
                            <td id="T_9fd1c_row3_col0" class="data row3 col0" >Siltstone</td>
                            <td id="T_9fd1c_row3_col1" class="data row3 col1" >Strat_Series</td>
                            <td id="T_9fd1c_row3_col2" class="data row3 col2" >3</td>
                            <td id="T_9fd1c_row3_col3" class="data row3 col3" >#ffbe00</td>
                            <td id="T_9fd1c_row3_col4" class="data row3 col4" >4</td>
                </tr>
                <tr>
                            <th id="T_9fd1c_level0_row4" class="row_heading level0 row4" >3</th>
                            <td id="T_9fd1c_row4_col0" class="data row4 col0" >Sandstone_2</td>
                            <td id="T_9fd1c_row4_col1" class="data row4 col1" >Strat_Series</td>
                            <td id="T_9fd1c_row4_col2" class="data row4 col2" >4</td>
                            <td id="T_9fd1c_row4_col3" class="data row4 col3" >#728f02</td>
                            <td id="T_9fd1c_row4_col4" class="data row4 col4" >5</td>
                </tr>
                <tr>
                            <th id="T_9fd1c_level0_row5" class="row_heading level0 row5" >5</th>
                            <td id="T_9fd1c_row5_col0" class="data row5 col0" >basement</td>
                            <td id="T_9fd1c_row5_col1" class="data row5 col1" >Strat_Series</td>
                            <td id="T_9fd1c_row5_col2" class="data row5 col2" >5</td>
                            <td id="T_9fd1c_row5_col3" class="data row5 col3" >#ff3f20</td>
                            <td id="T_9fd1c_row5_col4" class="data row5 col4" >6</td>
                </tr>
        </tbody></table>
    </div>
    <br />
    <br />

.. GENERATED FROM PYTHON SOURCE LINES 51-57

Model generation
----------------

After loading in the data, we set it up for interpolation and compute
the model.


.. GENERATED FROM PYTHON SOURCE LINES 59-60

set up interpolator

.. GENERATED FROM PYTHON SOURCE LINES 60-65

.. code-block:: python3

    gp.set_interpolator(geo_model,
                        compile_theano=True,
                        theano_optimizer='fast_compile',
                        verbose=[])





.. 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               3464.1
    $C_o$            285714.29
    drift equations     [3, 3]

    <gempy.core.interpolator.InterpolatorModel object at 0x7fcb8b722e50>



.. GENERATED FROM PYTHON SOURCE LINES 66-67

compute the model

.. GENERATED FROM PYTHON SOURCE LINES 67-69

.. code-block:: python3

    gp.compute_model(geo_model, compute_mesh=False);





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

 Out:

 .. code-block:: none



    Lithology ids 
      [6. 6. 6. ... 2. 2. 2.] 




.. GENERATED FROM PYTHON SOURCE LINES 70-71

have a look at the data and computed model

.. GENERATED FROM PYTHON SOURCE LINES 71-73

.. code-block:: python3

    gp.plot_3d(geo_model)




.. image:: /integrations/images/sphx_glr_gempy_export_MOOSE_001.png
    :alt: gempy export MOOSE
    :class: sphx-glr-single-img


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

 Out:

 .. code-block:: none


    <gempy.plot.vista.GemPyToVista object at 0x7fcb89256100>



.. GENERATED FROM PYTHON SOURCE LINES 74-83

Exporting the Model to MOOSE
----------------------------

The voxel-model above already is the same as a model discretized in a
hexahedral grid, so my immediately be used as input in a simulation
tool, e.g. `MOOSE <https://mooseframework.org/>`__. For this, we need to
access to the unit IDs assigned to each voxel in GemPy. The array
containing these IDs is called ``lith_block``.


.. GENERATED FROM PYTHON SOURCE LINES 85-88

.. code-block:: python3

    ids = geo_model.solutions.lith_block
    print(ids)





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

 Out:

 .. code-block:: none

    [6. 6. 6. ... 2. 2. 2.]




.. GENERATED FROM PYTHON SOURCE LINES 89-95

This array has the shape of ``(x,)`` and would be immediately useful, if
GemPy and the chosen simulation code would *populate* a grid in the same
way. Of course, however, that is not the case. This is why we have to
restructure the ``lith_block`` array, so it can be read correctly by
MOOSE.


.. GENERATED FROM PYTHON SOURCE LINES 97-98

model resolution

.. GENERATED FROM PYTHON SOURCE LINES 98-103

.. code-block:: python3

    nx, ny, nz = geo_model.grid.regular_grid.resolution

    # model extent
    xmin, xmax, ymin, ymax, zmin, zmax = geo_model.grid.regular_grid.extent








.. GENERATED FROM PYTHON SOURCE LINES 104-109

These two parameters are important to, a) restructure ``lith_block``,
and b) write the input file for MOOSE correctly. For a), we need to
reshape ``lith_block`` again to its three dimensions and *re-flatten* it
in a *MOOSE-conform* way.


.. GENERATED FROM PYTHON SOURCE LINES 111-112

reshape to 3D array

.. GENERATED FROM PYTHON SOURCE LINES 112-116

.. code-block:: python3

    units = ids.reshape((nx, ny, nz))
    # flatten MOOSE conform
    units = units.flatten('F')








.. GENERATED FROM PYTHON SOURCE LINES 117-157

| The importance of ``nx, ny, nz`` is apparent from the cell above. But
  what about ``xmin``, …, ``zmax``?
| A MOOSE input-file for mesh generation has the following syntax:

.. code:: python

   [MeshGenerators]
     [./gmg]
       type = GeneratedMeshGenerator
       dim = 3
       nx = 50
       ny = 50
       nz = 80
       xmin = 0.0
       xmax = 2000.0
       yim = 0.0
       ymax = 2000.0
       zmin = 0.0
       zmax = 2000.0
       block_id = '1 2 3 4 5 6'
       block_name = 'Main_Fault Sandstone_2 Siltstone Shale Sandstone_1 basement'
     [../]

     [./subdomains]
       type = ElementSubdomainIDGenerator
       input = gmg
       subdomain_ids = ' ' # here you paste the transformed lith_block vector
     [../]
   []

   [Mesh]
     type = MeshGeneratorMesh
   []

So these parameters are required inputs in the ``[MeshGenerators]``
object in the MOOSE input file. ``GemPy`` has a method to directly
create such an input file, stored in ``gempy.utils.export.py``.

The following cell shows how to call the method:


.. GENERATED FROM PYTHON SOURCE LINES 159-164

.. code-block:: python3


    # sphinx_gallery_thumbnail_path = '_static/GemPy_model_combined.png'
    import gempy.utils.export as export
    export.export_moose_input(geo_model, path='')





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

 Out:

 .. code-block:: none

    Successfully exported geological model as moose input to ./




.. GENERATED FROM PYTHON SOURCE LINES 165-207

This method automatically stores a file
``geo_model_units_moose_input.i`` at the specified path. Either this
input file could be extended with parameters to directly run a
simulation, or it is used just for creating a mesh. In the latter case,
the next step would be, to run the compiled MOOSE executable witch the
optional flag ``--mesh-only``.

E.g. with using the `PorousFlow
module <https://mooseframework.inl.gov/modules/porous_flow/>`__:

.. code:: bash

   $path_to_moose/moose/modules/porous_flow/porous_flow-opt -i pct_voxel_mesh.i --mesh-only

How to compile MOOSE is described in their
`documentation <https://mooseframework.inl.gov/getting_started/index.html>`__.

The now generated mesh with the name
``geo_model_units_moose_input_in.e`` can be used as input for another
MOOSE input file, which contains the main simulation parameters. To call
the file with the grid, the following part has to be added in the MOOSE
simulation input file:

.. code:: python

   [Mesh]
     file = geo_model_units_moose_input_in.e
   []

.. raw:: html

   <hr>

The final output of the simulation may also be such an ``.e``, which
can, for instance, be opened with
`paraview <https://www.paraview.org/>`__. A simulated temperature field
(purely conductive) of the created model would look like this:

.. figure:: https://raw.githubusercontent.com/Japhiolite/a-Moose-and-you/master/imgs/GemPy_model_combined.png
   :alt: gempy_temperature




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

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


.. _sphx_glr_download_integrations_gempy_export_MOOSE.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: gempy_export_MOOSE.py <gempy_export_MOOSE.py>`



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

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


.. only:: html

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

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