From d37c008e0dd8df73b9a341cbaf53b50d709e6756 Mon Sep 17 00:00:00 2001 From: mivaia <59896845+mivaia@users.noreply.github.com> Date: Thu, 21 Nov 2024 09:43:15 -0500 Subject: [PATCH] Review warming-up-a-viscous-fluid example. Post-processing python file added. Modified the doc to take the new boundaries default parameters into account (#1376) Results are identical. Added a post-processing file for the example. Modified the boundary condition section in the documentation. ALL the boundaries are now described. Beta is also set to 0 because flow is entering in that kind of problem. --- .../warming-up-a-viscous-fluid.rst | 80 ++++++++++------ ...warming-up-viscous-fluid-postprocessing.py | 93 +++++++++++++++++++ .../warming-up-viscous-fluid.prm | 11 +-- 3 files changed, 147 insertions(+), 37 deletions(-) create mode 100644 examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid-postprocessing.py diff --git a/doc/source/examples/multiphysics/warming-up-a-viscous-fluid/warming-up-a-viscous-fluid.rst b/doc/source/examples/multiphysics/warming-up-a-viscous-fluid/warming-up-a-viscous-fluid.rst index 83115674d2..bcc4040e43 100644 --- a/doc/source/examples/multiphysics/warming-up-a-viscous-fluid/warming-up-a-viscous-fluid.rst +++ b/doc/source/examples/multiphysics/warming-up-a-viscous-fluid/warming-up-a-viscous-fluid.rst @@ -118,7 +118,7 @@ The ``mesh`` considered is a very basic rectangle, using the ``dealii`` grid typ Multiphysics ~~~~~~~~~~~~~~ -The ``multiphysics`` subsection enable to turn on (``true``) and off (``false``) the physics of interest. Here ``heat transfer`` and ``viscous dissipation`` must be set (see Bonuses for results without viscous dissipation). +The ``multiphysics`` subsection enables to turn on (``true``) and off (``false``) the physics of interest. Here ``heat transfer`` and ``viscous dissipation`` must be set (see Bonuses for results without viscous dissipation). .. code-block:: text @@ -156,13 +156,13 @@ Boundary Conditions The ``boundary conditions`` are set for: -* the fluid dynamic in ``subsection boundary conditions``, with ``noslip`` at the left wall (``bc 0``) and a velocity of ``2`` in the y-direction at the right wall (``bc 1``), -* the heat transfer in ``subsection boundary conditions heat transfer``, with a ``convection`` imposed at the left wall (``bc 0``) with a heat transfer coefficient ``h = 0`` to represent an insulation condition, and an imposed ``temperature`` of ``80`` at the right wall. +* the fluid dynamic in ``subsection boundary conditions``, with ``noslip`` at the left wall (``bc 0``) and a velocity of ``2`` in the y-direction at the right wall (``bc 1``). The other walls (``bc 2`` and ``bc 3``) are set as ``outlet`` with a ``beta = 0`` to represent an open boundary condition. +* the heat transfer in ``subsection boundary conditions heat transfer``, with an imposed ``temperature`` of ``80`` at the right wall. All the other walls are set as ``noflux``. .. code-block:: text subsection boundary conditions - set number = 2 + set number = 4 subsection bc 0 set id = 0 set type = noslip @@ -177,19 +177,23 @@ The ``boundary conditions`` are set for: set Function expression = 2 end end + subsection bc 2 + set id = 2 + set type = outlet + set beta = 0 + end + subsection bc 3 + set id = 3 + set type = outlet + set beta = 0 + end end subsection boundary conditions heat transfer - set number = 2 + set number = 4 subsection bc 0 set id = 0 - set type = convection-radiation-flux - subsection h - set Function expression = 0 - end - subsection Tinf - set Function expression = 0 - end + set type = noflux end subsection bc 1 set id = 1 @@ -198,8 +202,18 @@ The ``boundary conditions`` are set for: set Function expression = 80 end end + subsection bc 2 + set id = 2 + set type = noflux + end + subsection bc 3 + set id = 3 + set type = noflux + end end +.. note:: + In lethe, beta = 1 is the default value for the ``outlet`` boundary condition. Beta acts as a penalizing factor for fluid entering the domain. Because the fluid will be forced to enter through one of the outlets, beta is set to 0. ----------------------- Running the Simulation @@ -336,7 +350,7 @@ Several adjustments have to be made in the `.prm` to turn the domain clockwise, .. code-block:: text subsection boundary conditions - set number = 2 + set number = 4 subsection bc 0 set id = 2 set type = noslip @@ -345,41 +359,47 @@ Several adjustments have to be made in the `.prm` to turn the domain clockwise, set id = 3 set type = function subsection u - set Function expression = 2 + set Function expression = 0 end subsection v - set Function expression = 0 + set Function expression = 2 end end + subsection bc 2 + set id = 0 + set type = outlet + set beta = 0 + end + subsection bc 3 + set id = 1 + set type = outlet + set beta = 0 + end end subsection boundary conditions heat transfer set number = 4 - subsection bc 2 + subsection bc 0 set id = 2 - set type = convection-radiation-flux - subsection h - set Function expression = 0 - end - subsection Tinf - set Function expression = 0 - end + set type = noflux end - subsection bc 3 + subsection bc 1 set id = 3 set type = temperature subsection value set Function expression = 80 end end + subsection bc 2 + set id = 0 + set type = noflux + end + subsection bc 3 + set id = 1 + set type = noflux + end end -.. important:: - For the fluid ``boundary conditions``, we use ``set number = 2``, whereas for ``boundary conditions heat transfer`` we use ``set number = 4``. These two notations are perfectly equivalent, as the boundary conditions are ``none`` by default (or ``noflux`` in the case of heat transfer, see :doc:`../../../parameters/cfd/boundary_conditions_multiphysics`). However, it is important to make sure that: - - * the index in ``subsection bc ..`` is coherent with the ``number`` set (if ``number = 2``, ``bc 0`` and ``bc 1`` are created but ``bc 2`` does not exist), - * the index in ``set id = ..`` is coherent with the ``id`` of the boundary in the mesh (here, the deal.II generated mesh). - ---------------------------- Possibilities for Extension diff --git a/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid-postprocessing.py b/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid-postprocessing.py new file mode 100644 index 0000000000..c58dafb2ea --- /dev/null +++ b/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid-postprocessing.py @@ -0,0 +1,93 @@ +# SPDX-FileCopyrightText: Copyright (c) 2023 The Lethe Authors +# SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception OR LGPL-2.1-or-later + +""" +IMPORTS +""" +from lethe_pyvista_tools import * +import matplotlib.pyplot as plt +import numpy as np +import pyvista as pv + +""" +FUNCTION +""" +def process_temperature(params, prm_folder, prm_file_name, output_folder, pvd_file_name): + + temperature = lethe_pyvista_tools(prm_folder, prm_file_name, pvd_file_name) + line = pv.Line(np.array([0, 0.5, 0]), np.array([0.5, 0.5, 0]), resolution = 100) + + time_list = temperature.time_list + + max_temp_list = np.zeros(len(temperature.list_vtu)) + min_temp_list = np.zeros(len(temperature.list_vtu)) + mean_temp_list = np.zeros(len(temperature.list_vtu)) + + sim_error = np.zeros(len(temperature.list_vtu)) + + for i in range(len(temperature.list_vtu)): + sim = pv.read(output_folder + "/" + temperature.list_vtu[i]) + + max_temp_list[i] = max(sim["temperature"]) + min_temp_list[i] = min(sim["temperature"]) + mean_temp_list[i] = np.mean(sim["temperature"]) + + position_on_sampled_line = line.sample(sim)["Distance"] + temperature_on_sampled_line = line.sample(sim)["temperature"] + + # Calculate the analytical solution + analytical_temperature = params.Tw+(((params.rho*params.nu)*params.v*params.v)/(2*params.K))*(1-(position_on_sampled_line/params.B)*(position_on_sampled_line/params.B)) + sim_error[i] = np.sqrt(np.sum((analytical_temperature - temperature_on_sampled_line)**2)) + + return time_list, max_temp_list, min_temp_list, mean_temp_list, sim_error + +""" +PARAMETERS +""" +prm_folder = "." +prm_file_name = "warming-up-viscous-fluid" + +output_folder = "output" +pvd_file_name = "warming-up.pvd" + +class params: + rho = 0.9 + nu = 0.5 + K = 0.12 + Tw = 80 + v = 2 + B = 0.5 + +""" +MAIN +""" +time_list, max_temp_list, min_temp_list, mean_temp_list, sim_error = process_temperature(params, prm_folder, prm_file_name, output_folder, pvd_file_name) + +""" +PLOTS FOR TEMPERATURE +""" +plt.figure(figsize=(8, 6)) +plt.plot(time_list, mean_temp_list, color='green', label='Mean Temperature') +plt.fill_between(time_list, min_temp_list, max_temp_list, color='green', alpha=0.3, label='Temperature envelope') + +# Labeling +plt.title('Temperature envelope and mean value', fontsize=14) +plt.xlabel('Time (s)', fontsize=12) +plt.ylabel('Temperature (C)', fontsize=12) +plt.ylim(0, 90) +plt.xlim(0, 7) +plt.grid(True, linestyle='--', alpha=0.5) +plt.legend() + +""" +PLOTS FOR ERROR +""" +plt.figure(figsize=(8, 6)) +plt.title('Error on the temperature according to time', fontsize=14) +plt.plot(time_list, sim_error, color='red', label='L2 Error') +plt.xlabel('Time (s)', fontsize=12) +plt.ylabel('L2 norm of the error on temperature (-)', fontsize=12) +plt.xlim(0, 7) +plt.grid(True, linestyle='--', alpha=0.5) + +plt.show() \ No newline at end of file diff --git a/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid.prm b/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid.prm index e5424ca48b..f7dee4132c 100644 --- a/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid.prm +++ b/examples/multiphysics/warming-up-viscous-fluid/warming-up-viscous-fluid.prm @@ -15,6 +15,7 @@ subsection simulation control set time step = 0.05 set time end = 7.0 set output name = warming-up + set output path = ./output/ set output frequency = 1 end @@ -128,10 +129,12 @@ subsection boundary conditions subsection bc 2 set id = 2 set type = outlet + set beta = 0 end subsection bc 3 set id = 3 set type = outlet + set beta = 0 end end @@ -139,13 +142,7 @@ subsection boundary conditions heat transfer set number = 4 subsection bc 0 set id = 0 - set type = convection-radiation-flux - subsection h - set Function expression = 0 - end - subsection Tinf - set Function expression = 0 - end + set type = noflux end subsection bc 1 set id = 1