Commit 0e34b354 authored by Lambert Theisen's avatar Lambert Theisen 🚀

Add example of bodyforce-driven Knudsen paradox in channel

- This replaces the pressure-driven example
- Matches better with Torrilhon2016 that also uses a constant dimensionless
force
- The plots, however, are equal because it's completely equivalent
parent 12f6dde1
......@@ -101,7 +101,7 @@ We provide a simple example of a flow through a finite-length channel in 2D.
.. code-block:: bash
# Move to folder:
cd examples/channel_flow
cd examples/channel_flow_force
# Create mesh:
./create_mesh.sh
# Run program with given input file:
......@@ -112,7 +112,7 @@ In the output folder the results can be post-processed to demonstrate the `Knuds
.. code-block:: bash
# Go to folder with simulation results (=casename in input.yml)
cd channel_flow
cd channel_flow_force
# Generate correlation data between Knudsen number and massflow
bash postprocessing.sh
cat table.csv
......
# General
# =======
# - output_folder: Used as output folder
output_folder: results_channel_flow_force
# Meshes
# ======
# - meshes: List of input meshes in h5 format to run simulations on
meshes:
- channel5.h5
# Numerical Parameters
# ====================
# - elements: Must contain the fields: theta, s, p, u, sigma
# - fields: List of FEM parameters (shape, degree)
# - shape: Element shape, e.g. Lagrange
# - degree: Element degree, e.g. 2
# - stabilization: Must contain cip
# - cip: Collection of Continous Interior Penalty (CIP) parameters
# - enable: Enable CIP stabilization
# - delta_1: Stabilization of grad(T)*grad(T_test) over edge
# - delta_2: Stabilization of grad(u)*grad(u_test) over edge
# - delta_3: Stabilization of grad(p)*grad(p_test) over edge
elements:
theta:
shape: Lagrange
degree: 1
s:
shape: Lagrange
degree: 1
p:
shape: Lagrange
degree: 1
u:
shape: Lagrange
degree: 1
sigma:
shape: Lagrange
degree: 1
stabilization:
cip:
enable: True
delta_1: 1.0
delta_2: 1.0
delta_3: 0.1
# Formulation Parameters
# ======================
# - nsd: Number of spatial dimensions == 2
# - mode: Formulation mode, one of heat, stress, r13
# - kn: Knudsen number
# - xi_tilde: Refaction coefficient in Maxwell accomodation model
# - heat_source: Heat source function for mode==heat||r13
# - mass_source: Mass source function for mode==stress||r13
# - body_force: Body force for mode==stress||r13
nsd: 2
mode: r13
kn: 1.0
xi_tilde: 1.0
heat_source: 0
mass_source: 0
body_force: [1,0]
# Boundary Conditions
# ===================
# - bcs: Dictionary of all boundary IDs from mesh
# - bc_id: must contain theta_w, u_t_w, u_n_w, p_w, epsilon_w
# - theta_w: Value for temperature at wall
# - u_t_w: Value for tangential velocity at wall
# - u_n_w: Value for normal velocity at wall
# - p_w: Value for pressure at wall
# - epsilon_w: Inflow-model parameter <=> Weight of pressure prescription
bcs:
3001: # bot
theta_w: 1
u_t_w: 0
u_n_w: 0
p_w: 0 # No pressure gradient, flow shall be bodyforce-induced!
epsilon_w: pow(10,-3)
3002: # right
theta_w: 1
u_t_w: 0
u_n_w: 0 # No slip
p_w: 0
epsilon_w: pow(10,+3) # No slip
3003: # top
theta_w: 1
u_t_w: 0
u_n_w: 0
p_w: 0 # No pressure gradient, flow shall be bodyforce-induced!
epsilon_w: pow(10,-3)
3004: # left
theta_w: 1
u_t_w: 0
u_n_w: 0 # No slip
p_w: 0
epsilon_w: pow(10,+3) # No slip
# Convergence Study
# =================
# - enable: Enable convergence study on given meshes
# - exact_solution: Path to exact solution in cpp-format to compare errors
# - plot: Show errors in matplotlib window. PDF output is always per default.
# - write_systemmatrix: Writes out systemmatrix (LHS) to use for analysis
# - rescale_pressure: Rescale numerical pressure result to have zero mean
# - relative_errors: Use relative errors. If exact sol. is zero, use absolute.
convergence_study:
enable: False
exact_solution: esols/01_coeffs.cpp
plot: False # to avoid error exit code due to $DISPLAY
write_systemmatrix: False
rescale_pressure: True
relative_error: True
# Postprocessing
# ==============
# - write_pdfs: Write all solution fields as PDF plot
# - massflow: List of BC IDs to compute massflow J=int_bc dot(u,n) ds
postprocessing:
write_pdfs: True
massflow: [3002]
# Parameter Study
# ==============
# - enable: Repeat simulation with different p. values (study)
# - parameter_key: Key as list, e.g. ["elemenets", "p", "degree"]
# - parameter_values: List of value for parameter, e.g. [0.01,0.1,1,10]
parameter_study:
enable: True
parameter_key: ["kn"]
parameter_values: [0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0]
\ No newline at end of file
......@@ -4,7 +4,7 @@ rm table.csv
for KNUDSEN in 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.0
do
MASSFLOW=$(cat results_channel_flow$KNUDSEN/massflow_3002)
MASSFLOW=$(cat results_channel_flow_force$KNUDSEN/massflow_3002)
echo "$KNUDSEN, $MASSFLOW" >> table.csv
done
......
table.csv
*massflow_*
\ No newline at end of file
// Command line Parameters
If(!Exists(p))
p = 5;
EndIf
// Settings
res = 100;
Mesh.CharacteristicLengthMax = 1.0 * 2^(-p);
Mesh.MshFileVersion = 2.0;
length=4;
height=1;
Point(1001) = {0, 0, 0, res};
Point(1002) = {length, 0, 0, res};
Point(1003) = {length, height, 0, res};
Point(1004) = {0, height, 0, res};
Line(2001) = {1001,1002};
Line(2002) = {1002,1003};
Line(2003) = {1003,1004};
Line(2004) = {1004,1001};
Line Loop(3001) = {2001}; Physical Curve("bot",3001) = {2001};
Line Loop(3002) = {2002}; Physical Curve("right",3002) = {2002};
Line Loop(3003) = {2003}; Physical Curve("top",3003) = {2003};
Line Loop(3004) = {2004}; Physical Curve("left",3004) = {2004};
Plane Surface(4000) = {3001:3004}; Physical Surface("mesh",4000) = {4000};
#!/bin/bash
geoToH5 channel.geo channel5.h5 "-setnumber p 5"
\ No newline at end of file
# General
# =======
# - output_folder: Used as output folder
output_folder: results_channel_flow
output_folder: results_channel_flow_pressure
# Meshes
# ======
......
import matplotlib.pyplot as plt
import numpy as np
x, y = np.loadtxt("table.csv", delimiter=",", unpack=True)
gui = False
if not gui:
plt.switch_backend("agg")
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
ax.plot(x,y, "-o", label="Channel Flow")
ax.set_xscale("log")
plt.xlabel("Knudsen number")
plt.ylabel("Dimensionless Massflow")
plt.title("Knudsen Paradox")
plt.legend()
gui = False
if gui:
plt.show()
else:
fig.savefig("fig.pdf", dpi=150)
#!/bin/bash
rm table.csv
for KNUDSEN in 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.0
do
MASSFLOW=$(cat results_channel_flow_pressure$KNUDSEN/massflow_3002)
echo "$KNUDSEN, $MASSFLOW" >> table.csv
done
cat table.csv
\ No newline at end of file
......@@ -44,11 +44,23 @@ class TestExamples(object):
self.create_meshes(working_dir)
self.run_solver("input.yml", working_dir)
def test_channel_flow(self):
@pytest.mark.skip(reason="Skip, needs too much time")
def test_channel_flow_pressure(self):
r"""
Test the channel flow case and generate table with Kn vs. massflow.
Test the pressure-induced channel flow case
and generate table with Kn vs. massflow.
"""
working_dir = "examples/channel_flow"
working_dir = "examples/channel_flow_pressure"
self.create_meshes(working_dir)
self.run_solver("input.yml", working_dir)
subprocess.check_call(["bash", "postprocessing.sh"], cwd=working_dir)
def test_channel_flow_force(self):
r"""
Test the force-induced channel flow case
and generate table with Kn vs. massflow.
"""
working_dir = "examples/channel_flow_force"
self.create_meshes(working_dir)
self.run_solver("input.yml", working_dir)
subprocess.check_call(["bash", "postprocessing.sh"], cwd=working_dir)
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment