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

.. only:: html

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

        :ref:`Go to the end <sphx_glr_download_auto_examples_howto_conversions_plot_D_feflowlib_CT_simulation.py>`
        to download the full example code.

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

.. _sphx_glr_auto_examples_howto_conversions_plot_D_feflowlib_CT_simulation.py:


Workflow: Component-transport model - conversion, simulation, postprocessing
============================================================================
.. sectionauthor:: Julian Heinze (Helmholtz Centre for Environmental Research GmbH - UFZ)

In this example we show how a simple mass transport FEFLOW model can be
converted to a pyvista.UnstructuredGrid and then be simulated in OGS with the
component transport process.

.. GENERATED FROM PYTHON SOURCE LINES 12-13

0. Necessary imports

.. GENERATED FROM PYTHON SOURCE LINES 13-24

.. code-block:: Python

    import tempfile
    import xml.etree.ElementTree as ET
    from pathlib import Path

    import matplotlib.pyplot as plt

    import ogstools as ot
    from ogstools.examples import feflow_model_2D_CT_t_560

    ot.plot.setup.show_element_edges = True








.. GENERATED FROM PYTHON SOURCE LINES 25-27

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 27-44

.. code-block:: Python

    temp_dir = Path(tempfile.mkdtemp("feflow_test_simulation"))
    feflow_model = ot.FeflowModel(
        feflow_model_2D_CT_t_560, temp_dir / "2D_CT_model"
    )
    # name the feflow concentratiob result the same as in OGS for easier comparison
    feflow_model.mesh["single_species"] = feflow_model.mesh["single_species_P_CONC"]
    concentration = ot.variables.Scalar(
        data_name="single_species", output_name="concentration",
        data_unit="mg/l", output_unit="mg/l",
    )  # fmt: skip
    # The original mesh is clipped to focus on the relevant part of it, where
    # concentration is larger than 1e-9 mg/l. The rest of the mesh has concentration
    # values of 0.
    clipped_mesh = feflow_model.mesh.clip_scalar(
        scalars="single_species", invert=False, value=1.0e-9
    )
    ot.plot.contourf(clipped_mesh, concentration)



.. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_001.png
   :alt: plot D feflowlib CT simulation
   :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_001.png
   :class: sphx-glr-single-img


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

 .. code-block:: none


    <Figure size 1099.71x1080 with 2 Axes>



.. GENERATED FROM PYTHON SOURCE LINES 45-46

2. Setup a prj-file to run a OGS-simulation.

.. GENERATED FROM PYTHON SOURCE LINES 46-52

.. code-block:: Python

    time_steps = list(zip([10] * 8, [8.64 * 10**i for i in range(8)], strict=False))
    feflow_model.setup_prj(end_time=int(4.8384e07), time_stepping=time_steps)
    # Save the model (mesh, subdomains and project file).
    feflow_model.save()
    # Print the prj-file as an example.
    ET.dump(ET.parse(feflow_model.mesh_path.with_suffix(".prj")))




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

 .. code-block:: none

    <OpenGeoSysProject>
        <meshes>
            <mesh>2D_CT_model.vtu</mesh>
            <mesh>single_species_P_BC_MASS.vtu</mesh>
        </meshes>
        <processes>
            <process>
                <name>CT</name>
                <type>ComponentTransport</type>
                <coupling_scheme>staggered</coupling_scheme>
                <integration_order>2</integration_order>
                <specific_body_force>0 0</specific_body_force>
                <secondary_variables>
                    <secondary_variable internal_name="darcy_velocity" output_name="v" />
                </secondary_variables>
                <process_variables>
                    <pressure>HEAD_OGS</pressure>
                    <concentration>single_species</concentration>
                </process_variables>
            </process>
        </processes>
        <media>
            <medium id="0">
                <phases>
                    <phase>
                        <type>AqueousLiquid</type>
                        <properties>
                            <property>
                                <name>viscosity</name>
                                <type>Constant</type>
                                <value>1</value>
                            </property>
                            <property>
                                <name>density</name>
                                <type>Constant</type>
                                <value>1</value>
                            </property>
                        </properties>
                        <components>
                            <component>
                                <name>single_species</name>
                                <properties>
                                    <property>
                                        <name>decay_rate</name>
                                        <type>Constant</type>
                                        <value>0.0</value>
                                    </property>
                                    <property>
                                        <name>pore_diffusion</name>
                                        <type>Constant</type>
                                        <value>3.5999998241701783e-10</value>
                                    </property>
                                    <property>
                                        <name>retardation_factor</name>
                                        <type>Constant</type>
                                        <value>16441.72737282367</value>
                                    </property>
                                </properties>
                            </component>
                        </components>
                    </phase>
                </phases>
                <properties>
                    <property>
                        <name>porosity</name>
                        <type>Constant</type>
                        <value>0.10999999940395355</value>
                    </property>
                    <property>
                        <name>longitudinal_dispersivity</name>
                        <type>Constant</type>
                        <value>0.0</value>
                    </property>
                    <property>
                        <name>transversal_dispersivity</name>
                        <type>Constant</type>
                        <value>0.0</value>
                    </property>
                    <property>
                        <name>permeability</name>
                        <type>Constant</type>
                        <value>1.1574074074074073e-05</value>
                    </property>
                </properties>
            </medium>
        </media>
        <time_loop>
            <processes>
                <process ref="CT">
                    <nonlinear_solver>basic_picard</nonlinear_solver>
                    <convergence_criterion>
                        <type>DeltaX</type>
                        <norm_type>NORM2</norm_type>
                        <reltol>1e-6</reltol>
                    </convergence_criterion>
                    <time_discretization>
                        <type>BackwardEuler</type>
                    </time_discretization>
                    <time_stepping>
                        <type>FixedTimeStepping</type>
                        <t_initial>0</t_initial>
                        <t_end>48384000</t_end>
                        <timesteps>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>8.64</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>86.4</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>864.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>8640.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>86400.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>864000.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>8640000.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>86400000.0</delta_t>
                            </pair>
                        </timesteps>
                    </time_stepping>
                </process>
                <process ref="CT">
                    <nonlinear_solver>basic_picard</nonlinear_solver>
                    <convergence_criterion>
                        <type>DeltaX</type>
                        <norm_type>NORM2</norm_type>
                        <reltol>1e-6</reltol>
                    </convergence_criterion>
                    <time_discretization>
                        <type>BackwardEuler</type>
                    </time_discretization>
                    <time_stepping>
                        <type>FixedTimeStepping</type>
                        <t_initial>0</t_initial>
                        <t_end>48384000</t_end>
                        <timesteps>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>8.64</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>86.4</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>864.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>8640.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>86400.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>864000.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>8640000.0</delta_t>
                            </pair>
                            <pair>
                                <repeat>10</repeat>
                                <delta_t>86400000.0</delta_t>
                            </pair>
                        </timesteps>
                    </time_stepping>
                </process>
            </processes>
            <output>
                <type>VTK</type>
                <prefix>/tmp/tmpeu4i8cwufeflow_test_simulation/2D_CT_model</prefix>
                <timesteps>
                    <pair>
                        <repeat>1</repeat>
                        <each_steps>1</each_steps>
                    </pair>
                </timesteps>
                <variables>
                    <variable>single_species</variable>
                    <variable>HEAD_OGS</variable>
                </variables>
                <fixed_output_times>48384000</fixed_output_times>
            </output>
            <global_process_coupling>
                <max_iter>1</max_iter>
                <convergence_criteria>
                    <convergence_criterion>
                        <type>DeltaX</type>
                        <norm_type>NORM2</norm_type>
                        <reltol>1e-10</reltol>
                    </convergence_criterion>
                    <convergence_criterion>
                        <type>DeltaX</type>
                        <norm_type>NORM2</norm_type>
                        <reltol>1e-10</reltol>
                    </convergence_criterion>
                </convergence_criteria>
            </global_process_coupling>
        </time_loop>
        <parameters>
            <parameter>
                <name>C0</name>
                <type>Constant</type>
                <value>0</value>
            </parameter>
            <parameter>
                <name>p0</name>
                <type>Constant</type>
                <value>0</value>
            </parameter>
            <parameter>
                <name>single_species_P_BC_MASS</name>
                <type>MeshNode</type>
                <mesh>single_species_P_BC_MASS</mesh>
                <field_name>single_species_P_BC_MASS</field_name>
            </parameter>
        </parameters>
        <process_variables>
            <process_variable>
                <name>single_species</name>
                <components>1</components>
                <order>1</order>
                <initial_condition>C0</initial_condition>
                <boundary_conditions>
                    <boundary_condition>
                        <type>Dirichlet</type>
                        <mesh>single_species_P_BC_MASS</mesh>
                        <parameter>single_species_P_BC_MASS</parameter>
                    </boundary_condition>
                </boundary_conditions>
            </process_variable>
            <process_variable>
                <name>HEAD_OGS</name>
                <components>1</components>
                <order>1</order>
                <initial_condition>p0</initial_condition>
            </process_variable>
        </process_variables>
        <nonlinear_solvers>
            <nonlinear_solver>
                <name>basic_picard</name>
                <type>Picard</type>
                <max_iter>100</max_iter>
                <linear_solver>general_linear_solver</linear_solver>
            </nonlinear_solver>
        </nonlinear_solvers>
        <linear_solvers>
            <linear_solver>
                <name>general_linear_solver</name>
                <lis>-i cg -p jacobi -tol 1e-10 -maxiter 100000</lis>
                <eigen>
                    <solver_type>CG</solver_type>
                    <precon_type>DIAGONAL</precon_type>
                    <max_iteration_step>100000</max_iteration_step>
                    <error_tolerance>1e-10</error_tolerance>
                </eigen>
            </linear_solver>
        </linear_solvers>
    </OpenGeoSysProject>




.. GENERATED FROM PYTHON SOURCE LINES 53-54

3. Run the model.

.. GENERATED FROM PYTHON SOURCE LINES 54-55

.. code-block:: Python

    feflow_model.run()




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

 .. code-block:: none

    OGS finished with project file /tmp/tmpeu4i8cwufeflow_test_simulation/2D_CT_model.prj.
    Execution took 3.749140739440918 s
    Project file written to output.




.. GENERATED FROM PYTHON SOURCE LINES 56-58

4. Read the last timestep and plot the results along a line on the upper edge
of the mesh parallel to the x-axis.

.. GENERATED FROM PYTHON SOURCE LINES 58-68

.. code-block:: Python

    ogs_sim_res = ot.MeshSeries(temp_dir / "2D_CT_model.pvd")[-1]
    fig, ax = plt.subplots(1, 1, figsize=(16, 10))
    pts = [[0.038 + 1.0e-8, 0.005, 0], [0.045, 0.005, 0]]
    for i, mesh in enumerate([ogs_sim_res, feflow_model.mesh]):
        sample = mesh.sample_over_line(*pts)
        label = ["OGS", "FEFLOW"][i]
        ot.plot.line(
            sample, concentration, ax=ax, color="kr"[i], label=label, ls="-:"[i]
        )
    fig.tight_layout()



.. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_002.png
   :alt: plot D feflowlib CT simulation
   :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_002.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 69-70

5. Concentration difference plotted on the mesh.

.. GENERATED FROM PYTHON SOURCE LINES 70-74

.. code-block:: Python


    diff = ot.meshlib.difference(feflow_model.mesh, ogs_sim_res, concentration)
    diff_clipped = diff.clip_box([0.038, 0.045, 0, 0.01, 0, 0], invert=False)
    fig = ot.plot.contourf(diff_clipped, concentration.difference, fontsize=20)



.. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_003.png
   :alt: plot D feflowlib CT simulation
   :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_003.png
   :class: sphx-glr-single-img





.. GENERATED FROM PYTHON SOURCE LINES 75-76

5.1 Concentration difference plotted along a line.

.. GENERATED FROM PYTHON SOURCE LINES 76-81

.. code-block:: Python

    diff_sample = diff.sample_over_line(*pts)
    fig = ot.plot.line(
        diff_sample, concentration.difference, label="Difference FEFLOW-OGS"
    )
    fig.tight_layout()



.. image-sg:: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_004.png
   :alt: plot D feflowlib CT simulation
   :srcset: /auto_examples/howto_conversions/images/sphx_glr_plot_D_feflowlib_CT_simulation_004.png
   :class: sphx-glr-single-img






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

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


.. _sphx_glr_download_auto_examples_howto_conversions_plot_D_feflowlib_CT_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_D_feflowlib_CT_simulation.ipynb <plot_D_feflowlib_CT_simulation.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: plot_D_feflowlib_CT_simulation.py <plot_D_feflowlib_CT_simulation.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: plot_D_feflowlib_CT_simulation.zip <plot_D_feflowlib_CT_simulation.zip>`