.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/howto_conversions/plot_F_feflowlib_HT_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_F_feflowlib_HT_simulation.py: Workflow: Hydro-thermal model - conversion, simulation and post-processing ========================================================================== .. sectionauthor:: Julian Heinze (Helmholtz Centre for Environmental Research GmbH - UFZ) In this example we show how a simple hydro thermal FEFLOW model can be converted to a pyvista.UnstructuredGrid and then be simulated in ot. .. GENERATED FROM PYTHON SOURCE LINES 11-12 0. Necessary imports .. GENERATED FROM PYTHON SOURCE LINES 12-20 .. code-block:: Python import tempfile from pathlib import Path import pyvista as pv import ogstools as ot from ogstools.examples import feflow_model_2D_HT .. GENERATED FROM PYTHON SOURCE LINES 21-23 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 23-25 .. code-block:: Python temp_dir = Path(tempfile.mkdtemp("feflow_test_simulation")) feflow_model = ot.FeflowModel(feflow_model_2D_HT, temp_dir / "2D_HT_model.vtu") .. GENERATED FROM PYTHON SOURCE LINES 26-27 2. Plot the temperature simulated in FEFLOW on the mesh, and print information about the mesh. .. GENERATED FROM PYTHON SOURCE LINES 27-30 .. code-block:: Python feflow_temperature = ot.variables.temperature.replace(data_name="P_TEMP") ot.plot.contourf(feflow_model.mesh, feflow_temperature) print(feflow_model.mesh) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_001.png :alt: plot F feflowlib HT simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none UnstructuredGrid (0x780b32531e40) N Cells: 6260 N Points: 3228 X Bounds: -2.500e+02, 9.000e+02 Y Bounds: 3.500e+02, 6.500e+02 Z Bounds: 0.000e+00, 0.000e+00 N Arrays: 34 .. GENERATED FROM PYTHON SOURCE LINES 31-32 3. Extract the subdomains. .. GENERATED FROM PYTHON SOURCE LINES 32-42 .. 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.view_xy() plotter.show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_002.png :alt: plot F feflowlib HT simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_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_F_feflowlib_HT_simulation_002.vtksz .. GENERATED FROM PYTHON SOURCE LINES 43-45 4. Setup a prj-file to run a OGS-simulation. Get the ogs6py model to create a prj-file and run the simulation. .. GENERATED FROM PYTHON SOURCE LINES 45-46 .. code-block:: Python feflow_model.setup_prj(end_time=1e11, time_stepping=[(1, 1e10)]) .. GENERATED FROM PYTHON SOURCE LINES 47-48 5. Run the model. .. GENERATED FROM PYTHON SOURCE LINES 48-49 .. code-block:: Python feflow_model.run() .. rst-class:: sphx-glr-script-out .. code-block:: none OGS finished with project file /tmp/tmpda5d0hb6feflow_test_simulation/2D_HT_model.prj. Execution took 2.0455374717712402 s Project file written to output. .. GENERATED FROM PYTHON SOURCE LINES 50-51 6. Read the results and plot them. .. GENERATED FROM PYTHON SOURCE LINES 51-65 .. code-block:: Python ms = ot.MeshSeries(temp_dir / "2D_HT_model.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 / "sim_2D_HT_model_ts_10_t_100000000000.000000.vtu" ) """ # Plot the hydraulic head/height, which was simulated in OGS. hydraulic_head = ot.variables.Scalar( data_name="HEAD_OGS", data_unit="m", output_unit="m" ) ot.plot.contourf(ogs_sim_res, hydraulic_head) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_003.png :alt: plot F feflowlib HT simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none
.. GENERATED FROM PYTHON SOURCE LINES 66-67 Plot the temperature, which was simulated in OGS. .. GENERATED FROM PYTHON SOURCE LINES 67-68 .. code-block:: Python ot.plot.contourf(ogs_sim_res, ot.variables.temperature) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_004.png :alt: plot F feflowlib HT simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none
.. GENERATED FROM PYTHON SOURCE LINES 69-70 7. Plot the difference between the FEFLOW and OGS simulation. .. GENERATED FROM PYTHON SOURCE LINES 70-78 .. code-block:: Python feflow_model.mesh["HEAD"] = feflow_model.mesh["P_HEAD"] ogs_sim_res["HEAD"] = ogs_sim_res["HEAD_OGS"] # Plot differences in hydraulic head/height. diff_mesh = ot.meshlib.difference(feflow_model.mesh, ogs_sim_res, "HEAD") hydraulic_head_diff = ot.variables.Scalar( data_name="HEAD_difference", data_unit="m", output_unit="m" ) ot.plot.contourf(diff_mesh, hydraulic_head_diff, vmin=-1.5e-9, vmax=1.5e-9) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_005.png :alt: plot F feflowlib HT simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_005.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none
.. GENERATED FROM PYTHON SOURCE LINES 79-80 Plot differences in temperature. .. GENERATED FROM PYTHON SOURCE LINES 80-88 .. code-block:: Python feflow_model.mesh["temperature"] = feflow_model.mesh["P_TEMP"] # Plot differences in temperature. diff_mesh = ot.meshlib.difference( feflow_model.mesh, ogs_sim_res, ot.variables.temperature ) ot.plot.contourf( diff_mesh, ot.variables.temperature.difference, vmin=-8.7e-9, vmax=8.7e-9 ) .. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_006.png :alt: plot F feflowlib HT simulation :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_F_feflowlib_HT_simulation_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none
.. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 4.392 seconds) .. _sphx_glr_download_auto_examples_howto_conversions_plot_F_feflowlib_HT_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_F_feflowlib_HT_simulation.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_F_feflowlib_HT_simulation.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_F_feflowlib_HT_simulation.zip `