Documentation for validation cases 1-3 idaholab#12

doc/content/index.md

@@ -27,3 +27,9 @@ effectively changing the concentration of the trapped specie from being measured
 in #/volume to k#/volume where 'k' indicates kilo. After this transformation,
 the numerical concentrations of the trapped and mobile species are on the same
 order of magnitude.
+
+## [Validation](verification/val-list.md)
+
+Several problems originally developed for the TMAP4 code have been used for the
+validation of the TMAP8 code. These validation cases can be found
+[here](verification/val-list.md).

doc/content/verification/bibfile.bib

@@ -0,0 +1,16 @@
+@techreport{longhurst1992verification,
+  title={Verification and Validation of TMAP4},
+  author={Longhurst, GR and Harms, SL and Marwil, ES and Miller, BG},
+  year={1992},
+  institution={EG and G Idaho, Inc., Idaho Falls, ID (United States)}
+}
+@article{longhurst2005verification,
+	title={Verification and validation of the tritium transport code TMAP7},
+	author={Longhurst, Glen R and Ambrosek, James},
+	journal={Fusion science and technology},
+	volume={48},
+	number={1},
+	pages={468--471},
+	year={2005},
+	publisher={Taylor \& Francis}
+}

doc/content/verification/val-1a.md

@@ -0,0 +1,48 @@
+# val-1a
+
+# Depleting Source Problem
+
+## Test Description
+
+This validation problem is taken from [!cite](longhurst1992verification). The model consists of an enclosure containing a finite concentration of atoms which are allowed to diffuse into a SiC layer over time. No solubility or trapping effects are included. The fractional release from the outside of the shell in a depleting source model in a slab geometry is given by:
+
+\begin{equation}
+    FR = 1.0 - \sum_{n=1}^{\infty} \frac{2\ L sec \ \alpha_{n} - \exp\left(\frac{-\alpha_{n}^2 D T}{l^{2}}\right)}{L(L+1) + \alpha_n^{2}}
+\end{equation}
+
+where
+
+\begin{equation}
+    L = \frac{lA}{V \phi}
+\end{equation}
+
+\begin{equation}
+    \phi = \frac{source \ concentration}{layer \ concentration}
+\end{equation}
+
+where the layer concentration is that at the interface with the source ($\phi$ is constant in time),
+
+    $A$ = surface area
+
+    $V$ = source volume
+
+    $l$ = layer thickness
+
+and the $\alpha_n$ are the roots of
+
+\begin{equation}
+    \alpha_n = \frac{L}{tan \ \alpha_n}
+\end{equation}
+
+## Results
+
+
+[val-1a_comparison] shows the comparison of the TMAP8 calculation and the analytical solution. There is good agreement between the two plots.
+
+
+!media figures/val-1a_comparison.png
+    style=width:50%;margin-bottom:2%
+    id=val-1a_comparison
+    caption=Comparison of TMAP8 calculation with the analytical solution
+
+!bibtex bibliography

doc/content/verification/val-1b.md

@@ -0,0 +1,72 @@
+# val-1b
+
+# Diffusion Problem with Constant Source Boundary Condition
+
+This validation problem is taken from [!cite](longhurst1992verification). Diffusion of tritium through a semi-infinite SiC layer is modeled with a constant
+source located on one boundary. No solubility or traping is included. The
+concentration as a function of time and position is given by:
+
+\begin{equation}
+C = C_o \; erfc \left(\frac{x}{2\sqrt{Dt}}\right)
+\end{equation}
+
+Comparison of the TMAP8 results and the analytical solution is shown in
+[val-1b_comparison_time] as a function of time at
+x = 0.2 mm. For simplicity, both the diffusion coefficient and the initial
+concentration were set to unity. The TMAP8 code predictions match very well with
+the analytical solution.
+
+!media figures/val-1b_comparison_time.png
+    style=width:50%;margin-bottom:2%
+    id=val-1b_comparison_time
+    caption=Comparison of concentration as function of time at x\=0.2m calculated
+     through TMAP8 and analytically
+
+As a second check, the concentration as a function of position at a given time
+t = 25s, from TMAP8 was compared with the analytical solution as shown in
+[val-1b_comparison_dist]. The predicted concentration profile from TMAP8 is in
+good agreement with the analytical solution.
+
+!media figures/val-1b_comparison_dist.png
+    style=width:50%;margin-bottom:2%
+    id=val-1b_comparison_dist
+    caption=Comparison of concentration as function of distance from the source
+    at t\=25sec calculated through TMAP8 and analytically
+
+Finally, the diffusive flux ($J$) was compared with the analytic solution where the
+flux is proportional to the derivative of the concentration with respect to x and
+is given by:
+
+\begin{equation}
+J = C_o \; \sqrt{\frac{D}{t\pi}} \; exp \left(\frac{x}{2\sqrt{Dt}}\right)
+\end{equation}
+
+The flux as given by Equation (?) is compared with values calculated by TMAP8 in
+Table ?. The diffusivity, D, and the initial concentration, C$_o$, were both
+taken as unity, and the distance, x, was taken as 0.5 in this comparison.
+TMAP8 initially under predicts but the results match well subsequently. Comparison
+results are shown in []
+
+!media figures/val-1b_comparison_flux.png
+    style=width:50%;margin-bottom:2%
+    id=val-1b_comparison_flux
+    caption=Comparison of flux as function of time at x\=0.5m calculated through
+    TMAP8 and analytically
+
+### Notes
+
+The trapping test features some oscillations in the solution for whatever
+reason. In order for the oscillations to not take over the simulation, it seems
+that the ratio of the **inverse of the Fourier number** must be kept
+sufficiently high, e.g. `h^2 / (D * dt)`. Included in this directory are three
+`png` files that show the permeation for different `h` and `dt` values. They are
+summarized below:
+
+- `nx-80.png`: `nx = 80` and `dt = .0625`
+- `nx-40.png`: `nx = 40` and `dt = .25`
+- `nx-20.png`: `nx = 20` and `dt = 1`
+
+The oscillations in the permeation graph go away with increasing fineness in the
+mesh and in `dt`.
+
+!bibtex bibliography

doc/content/verification/val-1c.md

@@ -0,0 +1,39 @@
+# val-1c
+
+# Diffusion Problem with Partially Preloaded Slab
+
+This validation problem is taken from [!cite](longhurst1992verification). Diffusion of tritium through a semi-infinite SiC layer is modeled with an initial
+loading of 1 atom/m{^3} in the first 10 m of a 2275-m slab. Solubility is unity
+and no traping is included. The analytical solution is given by:
+
+\begin{equation}
+C = \frac{C_o}{2} \left( erf \left( \frac{h-x}{2} \sqrt{Dt} \right) + erf \left( \frac{h+x}{2\sqrt{Dt}}) \right) \right)
+\end{equation}
+
+
+where h is the thickness of the pre-loaded portion of the layer.
+
+At the surface (x = 0) the concentration is given by:
+
+\begin{equation}
+C = \frac{C_o}{2} \; erf \left( \frac{h}{2\sqrt{Dt}} \right)
+\end{equation}
+
+while at x = h its value is described by
+
+\begin{equation}
+C = \frac{C_o}{2} \; erf \left( \frac{h}{\sqrt{Dt}} \right)
+\end{equation}
+
+A comparison of the mobile species concentration values at x = 0 m, 10 m and
+12.5 m calculated through TMAP8 and analytically is shown in
+[val-1c_comparison_time]. The TMAP8 calculations are found to be in good agreement
+with the analytical solution.
+
+!media figures/val-1c_comparison_time.png
+    style=width:50%;margin-bottom:2%
+    id=val-1c_comparison_time
+    caption=Comparison of concentration as function of time at x\=0 m, 10 m and 12 m
+    calculated through TMAP8 and analytically
+
+!bibtex bibliography

doc/content/verification/val-1d.md

@@ -1,13 +1,94 @@
 # val-1d
 
-This reproduces the val-1d verification test from TMAP4. We examine two
-different regimes, one where diffusion is the rate-limiting step, and one where
-trapping is the rate-limiting step.
+# Permeation Problem with Trapping
+
+## Test Description
+
+This validation problem is taken from [!cite](longhurst1992verification). It models permeation through a membrane with a constant source in which traps are operative. The breakthrough time may have one of two limiting values depending on whether the trapping is in the effective diffusivity or strong-trapping regime. A trapping parameter is defined by:
+
+\begin{equation}
+  \label{eqn:zeta}
+    \zeta = \frac{\lambda^2 \nu}{\rho D_o} exp \left( \frac{E_d - \epsilon}{kT} \right) + \frac{c}{\rho}
+\end{equation}
+
+where
+
+$\lambda$ = lattice parameter
+
+$\nu$ = Debye frequency ($\approx$ $10^{13} \; s^{-1}$)
+
+$\rho$ = trapping site fraction
+
+$D_o$ = diffusivity pre-exponential
+
+$E_d$ = diffusion activation energy
+
+$\epsilon$ = trap energy
+
+$k$ = Boltzmann's constant
+
+$T$ = temperature
+
+$c$ = dissolved gas atom fraction
+
+The discriminant for which regime is dominant is the ratio of $\zeta$ to c/$\rho$. If $\zeta$ > c/$\rho$ then the effective diffusivity regime applies, and the permeation transient is identical to the standard diffusion transient but with the diffusivity replaced by an effective diffusivity.
+
+\begin{equation}
+\label{eqn:Deff}
+    D_eff = \frac{D}{1 + \frac{1}{\zeta}}
+\end{equation}
+
+In this limit, the breakthrough time, defined as the intersection of the steepest tangent to the diffusion transient with the time axis, will be
+
+\begin{equation}
+\label{eqn:tau_be}
+    \tau_{b_e} = \frac{l^2}{2 \; \pi^2 \; D_eff}
+\end{equation}
+
+where $l$ is the thickness of the slab and D is the diffusivity of the gas through the material. The permeation transient is then given by
+
+
+\begin{equation}
+\label{eqn:Jp}
+    J_p = \frac{c_o D}{l} \Bigg\{ 1 + 2 \sum_{m=1}^{\infty} \left[ (-1)^m \exp \left( -m^2 \frac{t}{2 \; \tau_{b_e}} \right) \right] \Bigg\}
+\end{equation}
+
+
+[!cite](longhurst2005verification) where $\tau_{b_e}$ is defined in [eqn:tau_be]
+
+In the deep-trapping limit, $\zeta$ < c/$\rho$, and no permeation occurs until essentially all the traps have been filled. Then permeation rapidly turns on to its state value. The breakthrough time is given by
+
+\begin{equation}
+\label{eqn:tau_bd}
+    \tau_{b_d} = \frac{l^2 \rho}{2 \; c_o \; D}
+\end{equation}
+
+where $c_o$ is the steady dissolved gas concentration at the upstream (x = 0) side.
+
+Using TMAP8 we examine these two different regimes, one where diffusion is the rate-limiting step, and one where trapping is the rate-limiting step. The upstream-side starting concentration of 0.0001 atom fraction, a diffusivity of 1 $m^2$/s, a trapping site fraction of 0.1, $\lambda^2 = 10^{-15} \; m^2$, and a temperature of 1000 K is considered.
+
 
 ## Diffusion-limited
 
+For the effective diffusivity limit, we selected $\epsilon/k = 100 K$ to give $\zeta = 90.48 c/\rho$. The comparison results are presented in [val-1d_comparison_diffusion].
+
+!media figures/val-1d_comparison_diffusion.png
+    style=width:50%;margin-bottom:2%
+    id=val-1d_comparison_diffusion
+    caption=Permeation history of a slab subject to effective-diffusivity limit trapping.
+
 ## Trapping-limited
 
+For the deep trapping limit we took $\epsilon/k = 10000 K$ to give $\zeta = 0.04533 c/\rho$.  The comparison results are presented in [val-1d_comparison_trapping].
+
+!media figures/val-1d_comparison_trapping.png
+    style=width:50%;margin-bottom:2%
+    id=val-1d_comparison_trapping
+    caption=Permeation transient in a slab subject to strong trapping.
+
+
+
+
 ### Notes
 
 The trapping test features some oscillations in the solution for whatever
@@ -23,3 +104,5 @@ summarized below:
 
 The oscillations in the permeation graph go away with increasing fineness in the
 mesh and in `dt`.
+
+!bibtex bibliography

doc/content/verification/val-list.md

@@ -0,0 +1,8 @@
+# List of validation cases
+
+| Case | Title |
+| - | - |
+| val-1a | [Depleting Source Problem](val-1a.md) |
+| val-1b | [Diffusion Problem with Constant Source Boundary Condition](val-1b.md) |
+| val-1c | [Diffusion Problem with Partially Preloaded Slab](val-1c.md) |
+| val-1d | [Permeation Problem with Trapping](val-1d.md) |

test/tests/val-1a/analytical.csv

@@ -0,0 +1,47 @@
+time(s),frac_rel
+0,0
+1,6.00E-07
+2,0.000262
+3,0.002399
+4,0.007873
+5,0.016781
+6,0.028554
+7,0.042499
+8,0.06218
+9,0.074568
+10,0.091823
+11,0.109482
+12,0.127337
+13,0.145233
+14,0.163059
+15,0.180733
+16,0.198197
+17,0.21541
+18,0.232343
+19,0.248977
+20,0.265301
+21,0.281307
+22,0.296991
+23,0.312353
+24,0.327395
+25,0.342119
+26,0.356529
+27,0.37063
+28,0.384426
+29,0.397923
+30,0.411126
+31,0.424042
+32,0.436677
+33,0.449035
+34,0.461122
+35,0.472945
+36,0.484509
+37,0.49582
+38,0.506883
+39,0.517703
+40,0.528286
+41,0.538637
+42,0.54876
+43,0.558662
+44,0.568346
+45,0.577818

test/tests/val-1a/comparison_val-1a.py

@@ -0,0 +1,34 @@
+import csv
+import matplotlib.pyplot as plt
+import numpy as np
+from matplotlib import gridspec
+import pandas as pd
+
+#===============================================================================
+
+fig = plt.figure(figsize=[6.5,5.5])
+gs = gridspec.GridSpec(1,1)
+ax = fig.add_subplot(gs[0])
+
+tmap_sol = pd.read_csv("./val-1a_csv.csv")
+analytical_sol = pd.read_csv("./analytical.csv")
+tmap_time = tmap_sol['time']
+tmap_fr = tmap_sol['rhs_release']
+analytical_time = analytical_sol['time(s)']
+analytical_fr = analytical_sol['frac_rel']
+
+ax.plot(tmap_time,tmap_fr,label=r"TMAP8",c='tab:gray')
+ax.plot(analytical_time,analytical_fr,label=r"Analytical",c='k', linestyle='--')
+
+
+ax.set_xlabel(u'Time(s)')
+ax.set_ylabel(r"Fractional release")
+ax.legend(loc="best")
+ax.set_xlim(left=0)
+ax.set_xlim(right=45)
+ax.set_ylim(bottom=0)
+plt.grid(b=True, which='major', color='0.65', linestyle='--', alpha=0.3)
+
+ax.minorticks_on()
+plt.savefig('comparison.png', bbox_inches='tight');
+plt.close(fig)