.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/howto_conversions/plot_E_feflowlib_H_simulation.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_auto_examples_howto_conversions_plot_E_feflowlib_H_simulation.py: Workflow: Hydraulic model - conversion, simulation and post-processing ====================================================================== .. sectionauthor:: Julian Heinze (Helmholtz Centre for Environmental Research GmbH - UFZ) In this example we show how a simple flow/hydraulic FEFLOW model can be converted to a pyvista.UnstructuredGrid and then be simulated in ot. .. GENERATED FROM PYTHON SOURCE LINES 12-13 0. Necessary imports .. GENERATED FROM PYTHON SOURCE LINES 13-25 .. code-block:: Python import tempfile from pathlib import Path import numpy as np import pyvista as pv import ogstools as ot from ogstools.examples import feflow_model_box_Neumann from ogstools.feflowlib import ( FeflowModel, ) .. GENERATED FROM PYTHON SOURCE LINES 26-28 1. Load a FEFLOW model (.fem) as a FeflowModel object to further work it. During the initialisation, the FEFLOW file is converted. .. GENERATED FROM PYTHON SOURCE LINES 28-39 .. code-block:: Python temp_dir = Path(tempfile.mkdtemp("feflow_test_simulation")) feflow_model = FeflowModel(feflow_model_box_Neumann, temp_dir / "boxNeumann") pv.global_theme.colorbar_orientation = "vertical" feflow_model.mesh.plot( show_edges=True, off_screen=True, scalars="P_HEAD", cpos=[0, 1, 0.5], scalar_bar_args={"position_x": 0.1, "position_y": 0.25}, ) print(feflow_model.mesh) .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_001.png :alt: plot E feflowlib H simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_001.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /builds/ogs/tools/ogstools/docs/auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_001.vtksz .. rst-class:: sphx-glr-script-out .. code-block:: none UnstructuredGrid (0x780b3253e080) N Cells: 11462 N Points: 6768 X Bounds: 0.000e+00, 1.000e+02 Y Bounds: 0.000e+00, 1.000e+02 Z Bounds: -1.000e+02, 0.000e+00 N Arrays: 23 .. GENERATED FROM PYTHON SOURCE LINES 40-41 2. Extract the subdomains conditions. .. GENERATED FROM PYTHON SOURCE LINES 41-51 .. code-block:: Python subdomains = feflow_model.subdomains # Since there can be multiple boundary conditions in the subdomains, # they are plotted iteratively. plotter = pv.Plotter(shape=(len(subdomains), 1)) for i, (name, boundary_condition) in enumerate(subdomains.items()): # topsurface refers to a cell based boundary condition. if name != "topsurface": plotter.subplot(i, 0) plotter.add_mesh(boundary_condition, scalars=name) plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_002.png :alt: plot E feflowlib H simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_002.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /builds/ogs/tools/ogstools/docs/auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_002.vtksz .. GENERATED FROM PYTHON SOURCE LINES 52-53 3. Define endtime and time stepping in the project-file. .. GENERATED FROM PYTHON SOURCE LINES 53-57 .. code-block:: Python feflow_model.setup_prj( end_time=int(1e8), time_stepping=[(1, 10), (5, 100), (10, 1000), (10, 1e6), (1, 1e7)], ) .. GENERATED FROM PYTHON SOURCE LINES 58-59 4. Run the model .. GENERATED FROM PYTHON SOURCE LINES 59-60 .. code-block:: Python feflow_model.run() .. rst-class:: sphx-glr-script-out .. code-block:: none OGS finished with project file /tmp/tmpf371zdkcfeflow_test_simulation/boxNeumann.prj. Execution took 3.4284188747406006 s Project file written to output. .. GENERATED FROM PYTHON SOURCE LINES 61-62 5. Read the results and plot them. .. GENERATED FROM PYTHON SOURCE LINES 62-76 .. code-block:: Python ms = ot.MeshSeries(temp_dir / "boxNeumann.pvd") # Read the last timestep: ogs_sim_res = ms.mesh(ms.timesteps[-1]) """ It is also possible to read the file directly with pyvista: ogs_sim_res = pv.read(temp_dir / "boxNeumann_ts_1_t_1.000000.vtu") """ ogs_sim_res.plot( show_edges=True, off_screen=True, scalars="HEAD_OGS", cpos=[0, 1, 0.5], scalar_bar_args={"position_x": 0.1, "position_y": 0.25}, ) .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_003.png :alt: plot E feflowlib H simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_003.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /builds/ogs/tools/ogstools/docs/auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_003.vtksz .. GENERATED FROM PYTHON SOURCE LINES 77-78 5.1 Plot the hydraulic head simulated in OGS with :py:mod:`ogstools.plot.contourf`. .. GENERATED FROM PYTHON SOURCE LINES 78-80 .. code-block:: Python head = ot.variables.Scalar(data_name="HEAD_OGS", data_unit="m", output_unit="m") fig = ot.plot.contourf(ogs_sim_res.slice(normal="z", origin=[50, 50, 0]), head) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_004.png :alt: plot E feflowlib H simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 81-82 6. Calculate the difference to the FEFLOW simulation and plot it. .. GENERATED FROM PYTHON SOURCE LINES 82-91 .. code-block:: Python diff = feflow_model.mesh["P_HEAD"] - ogs_sim_res["HEAD_OGS"] feflow_model.mesh["diff_HEAD"] = diff feflow_model.mesh.plot( show_edges=True, off_screen=True, scalars="diff_HEAD", cpos=[0, 1, 0.5], scalar_bar_args={"position_x": 0.1, "position_y": 0.25}, ) .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_005.png :alt: plot E feflowlib H simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_005.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /builds/ogs/tools/ogstools/docs/auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_005.vtksz .. GENERATED FROM PYTHON SOURCE LINES 92-94 6.1 Plot the differences in the hydraulic head with :py:mod:`ogstools.plot.contourf`. Slices are taken along the z-axis. .. GENERATED FROM PYTHON SOURCE LINES 94-104 .. code-block:: Python diff_head = ot.variables.Scalar( data_name="diff_HEAD", data_unit="m", output_unit="m" ) slices = np.reshape( list(feflow_model.mesh.slice_along_axis(n=4, axis="z")), (2, 2) ) fig = ot.plot.contourf(slices, diff_head) for ax, slice in zip(fig.axes, np.ravel(slices), strict=False): ax.set_title(f"z = {slice.center[2]:.1f} m", fontsize=32) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_006.png :alt: z = -99.0 m, z = -66.3 m, z = -33.7 m, z = -1.0 m :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_006.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 105-106 6.2 Slices are taken along the y-axis. .. GENERATED FROM PYTHON SOURCE LINES 106-112 .. code-block:: Python slices = np.reshape( list(feflow_model.mesh.slice_along_axis(n=4, axis="y")), (2, 2) ) fig = ot.plot.contourf(slices, diff_head) for ax, slice in zip(fig.axes, np.ravel(slices), strict=False): ax.set_title(f"y = {slice.center[1]:.1f} m", fontsize=32) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_007.png :alt: y = 1.0 m, y = 33.7 m, y = 66.3 m, y = 99.0 m :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_E_feflowlib_H_simulation_007.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 7.748 seconds) .. _sphx_glr_download_auto_examples_howto_conversions_plot_E_feflowlib_H_simulation.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_E_feflowlib_H_simulation.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_E_feflowlib_H_simulation.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_E_feflowlib_H_simulation.zip `