Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Ver 1ka #183

Merged
merged 24 commits into from
Oct 14, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion doc/content/verification_and_validation/index.md
simopier marked this conversation as resolved.
Show resolved Hide resolved
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,7 @@ TMAP8 also contains [example cases](examples/tmap_index.md), which showcase how
| ver-1ha | [Convective Gas Outflow Problem](ver-1ha.md) |
| ver-1ja | [Radioactive Decay of Mobile Tritium in a Slab](ver-1ja.md) |
| ver-1jb | [Radioactive Decay of Mobile Tritium in a Slab with a Distributed Trap Concentration](ver-1jb.md) |

| ver-1ka | [Simple Volumetric Source](ver-1ka.md) |

# List of benchmarking cases

Expand Down
45 changes: 45 additions & 0 deletions doc/content/verification_and_validation/ver-1ka.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,45 @@
# ver-1ka

# Simple Volumetric Source

## General Case Description

This problem involves two enclosures connected by a diffusive membrane that follows Sieverts law for diffusion. Both enclosures contain hydrogen (H$_2$), tritium (T$_2$), and tritium gas (HT). In the first enclosure, there is an initial inventory of only hydrogen (H$_2$) along with a constant volumetric source rate of tritium (T$_2$). The second enclosure starts out empty.

## Case Set Up

This verification problem is taken from [!cite](ambrosek2008verification).
The rise in pressure of T$_2$ molecules in the first enclosure can be monitored by not enabling T$_2$ molecules to traverse the membrane between the two enclosures (no tritium flux). Consequently, the rate of pressure increase in the first enclosure can be expressed as:

\begin{equation} \label{eq:source_term}
\frac{dP_{T_2}}{dt} = \frac{S}{V} kT,
\end{equation}

where $S$ represents the volumetric T$_2$ source rate, $V$ is the volume of the enclosure, $k$ is the Boltzmann constant, and $T$ is the temperature of the enclosure.
In this case, $S$ is set to 10$^{20}$ molecules/m$^{-3}$/s, $V = 1$ m$^3$, and the temperature of the enclosure is constant at $T = 500$ K.

## Analytical Solution

The analytical solution for [eq:source_term] is simply

\begin{equation}
P_{T_2}(t) = \frac{S}{V} kT t.
\end{equation}

## Results
Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved

Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved
Comparison of the TMAP8 results and the analytical solution is shown in
Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved
[ver-1ka_comparison_time] as a function of time. The TMAP8 code predictions match very well with the analytical solution with a root mean squared percentage error of RMSPE $= 0.00$ %.

!media comparison_ver-1ka.py
Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved
image_name=ver-1ka_comparison_time.png
style=width:50%;margin-bottom:2%;margin-left:auto;margin-right:auto
id=ver-1ka_comparison_time
caption=Comparison of T$_2$ partial pressure in an enclosure with no loss pathways as function of time calculated through TMAP8 and analytically. TMAP8 matches the analytical solution.

## Input files

!style halign=left
The input file for this case can be found at [/ver-1ka.i], which is also used as tests in TMAP8 at [/ver-1ka/tests].

!bibtex bibliography
64 changes: 64 additions & 0 deletions test/tests/ver-1ka/comparison_ver-1ka.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
import numpy as np
Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved
import pandas as pd
import os
from matplotlib import gridspec
import matplotlib.pyplot as plt

# Changes working directory to script directory (for consistent MooseDocs usage)
script_folder = os.path.dirname(__file__)
os.chdir(script_folder)

# Constants
S = 1e20 # Source term (m^-3 * s^-1)
V = 1 # Volume (m^3)
kb = 1.380649e-23 # Boltzmann constant (J/K)
T = 500 # Temperature (K)

# Extract time and pressure data
if "/TMAP8/doc/" in script_folder: # if in documentation folder
csv_folder = "../../../../test/tests/ver-1ka/gold/ver-1ka_out.csv"
else: # if in test folder
csv_folder = "./gold/ver-1ka_out.csv"
expt_data = pd.read_csv(csv_folder)
TMAP8_time = expt_data['time']
TMAP8_pressure = expt_data['pressure']

# Calculate the theoretical expression for pressure
analytical_pressure = (S / V) * kb * T * TMAP8_time

fig = plt.figure(figsize=[6.5, 5.5])
gs = gridspec.GridSpec(1, 1)
ax = fig.add_subplot(gs[0])

# Plot the experimental data
ax.plot(TMAP8_time/3600, TMAP8_pressure, label=r"TMAP8", c='tab:gray',linestyle='-')

# Plot the selected theoretical data
ax.plot(TMAP8_time/3600, analytical_pressure, label=r"Analytical",c='k', linestyle='--')

RMSE_pressure = np.linalg.norm(TMAP8_pressure-analytical_pressure)
err_percent_pressure = RMSE_pressure*100/np.mean(analytical_pressure)

# Add text annotation for RMSPE on the plot
ax.text(0.05, 0.95, '$P_{T_2}$ RMSPE = %.2f ' % err_percent_pressure + '%', fontweight='bold')

# Format the y-axis to use scientific notation
plt.gca().yaxis.set_major_formatter(plt.FuncFormatter(lambda val, pos: '{:.1e}'.format(val)))

# Label the axes
ax.set_xlabel('Time (hr)')
ax.set_ylabel('Pressure (Pa)')
Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved

# define axis range
ax.set_xlim(left=0)
ax.set_xlim(right=3)
ax.set_ylim(bottom=0)

# Add a legend
ax.legend(loc="best")

# Add a grid
plt.grid(which='major', color='0.65', linestyle='--', alpha=0.3)

# Save the plot as a PNG file
plt.savefig('ver-1ka_comparison_time.png', bbox_inches='tight')
22 changes: 22 additions & 0 deletions test/tests/ver-1ka/gold/ver-1ka_out.csv
Original file line number Diff line number Diff line change
@@ -0,0 +1,22 @@
time,pressure
0,0
500,345.16225000001
1000,690.32450270146
1500,1035.4867536019
2000,1380.6490039021
2500,1725.8112540022
3000,2070.9735040355
3500,2416.1357540466
4000,2761.2980040503
4500,3106.4602540515
5000,3451.6225040519
5500,3796.784754052
6000,4141.9470040521
6500,4487.1092540521
7000,4832.2715040521
7500,5177.4337540521
8000,5522.5960040519
8500,5867.7582540518
9000,6212.9205040517
9500,6558.0827540518
10000,6903.2450040518
17 changes: 17 additions & 0 deletions test/tests/ver-1ka/tests
Original file line number Diff line number Diff line change
@@ -0,0 +1,17 @@
[Tests]
design = 'ODETimeDerivative.md ParsedODEKernel.md'
issues = '#12'
verification = 'ver-1ka.md'
[ver-1ka_csv]
type = CSVDiff
input = ver-1ka.i
csvdiff = ver-1ka_out.csv
requirement = 'The system shall be able to model a tritium volumetric source in one enclosure'
[]
[ver-1ka_comparison]
type = RunCommand
command = 'python3 comparison_ver-1ka.py'
requirement = 'The system shall be able to generate comparison plots between the analytical solution and simulated solution of verification case 1ka, modeling a tritium volumetric source in one enclosure.'
required_python_packages = 'matplotlib numpy pandas os git'
[]
[]
Lee01Atom marked this conversation as resolved.
Show resolved Hide resolved
48 changes: 48 additions & 0 deletions test/tests/ver-1ka/ver-1ka.i
Original file line number Diff line number Diff line change
@@ -0,0 +1,48 @@
initial_pressure = ${units 0 Pa} # initial internal pressure
kb = ${units 1.380649e-23 J/K} # Boltzmann constant J/K - from PhysicalConstants.h
T = ${units 500 K} # Temperature
S = ${units 1e20 1/m^3/s} # Source term
V = ${units 1 m^3} # Volume
end_time = ${units 10000 s}
time_step = ${units 500 s}

[Mesh]
type = GeneratedMesh
dim = 2
nx = 10
ny = 10
xmax = 1
ymax = 1
[]

[Variables]
[pressure]
family = SCALAR
order = FIRST
initial_condition = '${fparse initial_pressure}'
[]
[]

[ScalarKernels]
[time]
type = ODETimeDerivative
variable = pressure
[]
[source]
type = ParsedODEKernel
variable = pressure
expression = '${fparse - S/V * kb * T}'
[]
[]

[Executioner]
type = Transient
dt = ${time_step}
end_time = ${end_time}
scheme = 'bdf2'
[]

[Outputs]
file_base = 'ver-1ka_out'
csv = true
[]