Adding val-1e case

Ref: idaholab#12

doc/content/verification/ver-1e.md

@@ -0,0 +1,61 @@
+# val-1e
+
+# Diffusion in a composite layer
+
+# Depleting Source Problem
+
+## Test Description
+
+This verification problem is taken from [!cite](ambrosek2008verification). In this problem a composite structure of PyC and SiC is modeled with a constant concentration boundary condition of the free surface of PyC and zero concentration boundary condition on the free surface of the SiC. The steady state solution for the PyC is given as:
+
+\begin{equation}
+\label{eqn:steady_state_pyc}
+    C = C_o \left[1 + \frac{x}{l}  \left(\frac{a D_{PyC}}{a D_{PyC} + l D_{SiC}} - 1 \right) \right]
+\end{equation}
+
+while the concentration profile for the SiC layer is given as:
+
+\begin{equation}
+\label{eqn:steady_state_sic}
+    C = C_o \left(\frac{a+l-x}{l} \right) \left(\frac{a D_{PyC}}{a D_{PyC} + l D_{SiC}} \right)
+\end{equation}
+
+where
+
+    $x$ = distance from free surface of PyC
+
+    $a$ = thickness of the PyC layer (33 $\mu m$)
+
+    $l$ = thickness of the SiC layer (66 $\mu m$)
+
+    $C_o$ = concentration at the PyC free surface (3.0537 x 10$^{25}$ atoms/m$^3$)
+
+    $D_{PyC}$ = diffusivity in PyC (1.274 x 10$^{-7}$ m$^2$/sec)
+
+    $D_{SiC}$ = diffusivity in SiC (2.622 x 10$^{-11}$ m$^2$/sec)
+
+The analytical transient solution for the concentration in the SiC side of the composite slab is given as:
+
+\begin{equation}
+\label{eqn:transient}
+    C = C_o \Bigg\{ \frac{D_{PyC}(l-x)}{l D_{PyC} + a D_{SiC}} - 2 \sum_{n=1}^{\infty} \frac{\sin(a \lambda_n) \sin(k l \lambda_n) \sin \left[k (l-x) \lambda_n \right]}{\lambda_n \left[ a \sin^2(k l \lambda_n) + l \sin^2 (a \lambda_n) \right]} \exp(-D_{PyC} \lambda_n^2 t) \Bigg\}
+\end{equation}
+
+
+## Results
+
+[ver-1e_comparison_dist] shows the comparison of the TMAP8 calculation and the analytical solution for concentration after steady-state is reached.
+
+!media figures/ver-1e_comparison_dist.png
+    style=width:50%;margin-bottom:2%
+    id=ver-1e_comparison_dist
+    caption=Comparison of TMAP8 calculation with the analytical solution
+
+For transient solution comparison, the concentation at a point, which is 15.75 $\mu m$ away from the PyC-SiC interface into the SiC layer, is obtained using the TMAP code as well as analytically. [ver-1e_comparison_time] shows comparison of the TMAP calculation with the analytical solution for this transient case. There is good agreement between TMAP and the analytical solution for both steady state as well as transient cases.
+
+!media figures/ver-1e_comparison_time.png
+    style=width:50%;margin-bottom:2%
+    id=ver-1e_comparison_time
+    caption=Comparison of TMAP8 calculation with the analytical solution
+
+!bibtex bibliography

doc/content/verification/ver-list.md

@@ -6,5 +6,6 @@
 | ver-1b  | [Diffusion Problem with Constant Source Boundary Condition](ver-1b.md) |
 | ver-1c  | [Diffusion Problem with Partially Preloaded Slab](ver-1c.md)           |
 | ver-1d  | [Permeation Problem with Trapping](ver-1d.md)                          |
+| ver-1e  | [Diffusion in Composite Material Layers](ver-1e.md)                    |
 | ver-1fa | [Heat Conduction with Heat Generation](ver-1fa.md)                     |
 | ver-1fb | [Thermal Transient](ver-1fb.md)                                        |

test/tests/ver-1e/comparison_ver-1e.py

@@ -0,0 +1,61 @@
+import csv
+import matplotlib.pyplot as plt
+import numpy as np
+from matplotlib import gridspec
+import pandas as pd
+from scipy import special
+
+# ========= Comparison of concentration as a function of time ===================
+
+fig = plt.figure(figsize=[6.5, 5.5])
+gs = gridspec.GridSpec(1, 1)
+ax = fig.add_subplot(gs[0])
+
+tmap_sol = pd.read_csv("./ver-1e_csv.csv")
+tmap_time = tmap_sol['time']
+tmap_conc = tmap_sol['conc_point1']*1e25
+
+analytical_sol = pd.read_csv("./analytical_time.csv")
+analytical_time = analytical_sol['t']
+analytical_conc = analytical_sol['u']
+ax.plot(tmap_time, tmap_conc, label=r"TMAP8", c='tab:gray')
+ax.plot(analytical_time, analytical_conc,
+        label=r"Analytical", c='k', linestyle='--')
+
+
+ax.set_xlabel(u'Time (s)')
+ax.set_ylabel(r"Concentration (atoms/m$^3$)")
+ax.legend(loc="best")
+plt.grid(b=True, which='major', color='0.65', linestyle='--', alpha=0.3)
+
+ax.minorticks_on()
+plt.savefig('ver-1e_comparison_time.png', bbox_inches='tight')
+plt.close(fig)
+
+
+# ============ Comparison of concentration as a function of distance ============
+
+fig = plt.figure(figsize=[6.5, 5.5])
+gs = gridspec.GridSpec(1, 1)
+ax = fig.add_subplot(gs[0])
+
+tmap_sol = pd.read_csv("./ver-1e_u_vs_x_steadyState.csv")
+tmap_distance = tmap_sol['x']*1e6
+tmap_conc = tmap_sol['u']*1e25
+
+analytical_sol = pd.read_csv("./analytical_u_vs_x_steadyState.csv")
+analytical_distance = analytical_sol['x']*1e6
+analytical_conc = analytical_sol['u']*1e25
+
+ax.plot(tmap_distance, tmap_conc, label=r"TMAP8", c='tab:gray')
+ax.plot(analytical_distance, analytical_conc,
+        label=r"Analytical", c='k', linestyle='--')
+
+ax.set_xlabel(u'Distance ($\mu$m)')
+ax.set_ylabel(r"Concentration (atoms/m$^3$)")
+ax.legend(loc="best")
+plt.grid(b=True, which='major', color='0.65', linestyle='--', alpha=0.3)
+
+ax.minorticks_on()
+plt.savefig('ver-1e_comparison_dist.png', bbox_inches='tight')
+plt.close(fig)

test/tests/ver-1e/gold/ver-1e_csv.csv

@@ -0,0 +1,28 @@
+time,conc_point1
+0,0
+1,0.039724446516061
+3,0.44866390368448
+5,0.96674184367197
+7,1.2185269732843
+9,1.3956518479133
+11,1.5482494533886
+13,1.6469639683069
+15,1.7436396808399
+17,1.8125093853855
+19,1.8777838413292
+21,1.9300801469721
+23,1.9767579931054
+25,2.0175755192723
+27,2.0523991113006
+29,2.0844649390419
+31,2.1111764970046
+33,2.1363570592144
+35,2.1571964316116
+37,2.1769290931904
+39,2.1933473835442
+41,2.2087806453039
+43,2.2217850617483
+45,2.233839989968
+47,2.2441672957771
+49,2.2535778858794
+50,2.2577929195466

test/tests/ver-1e/tests

@@ -0,0 +1,10 @@
+[Tests]
+  design = 'Diffusion.md TimeDerivative.md DirichletBC.md'
+  issues = '#12'
+  [csv]
+    type = CSVDiff
+    input = ver-1e.i
+    csvdiff = ver-1e_csv.csv
+    requirement = 'The system shall be able to model transient diffusion through a composite slab with a constant concentration boundary conditions.'
+  []
+[]

test/tests/ver-1e/ver-1e.i

@@ -0,0 +1,93 @@
+[Mesh]
+  type = GeneratedMesh
+  dim = 1
+  nx = 1000
+  xmax = 99e-6
+[]
+
+[Variables]
+  [u]
+  []
+[]
+
+[Functions]
+  [diffusivity_value]
+    type = ParsedFunction
+    value = 'if(x<33e-6, 1.274e-7, 2.622e-11)'
+  []
+[]
+
+[Kernels]
+  [diff]
+    type = FunctionDiffusion
+    variable = u
+    function = diffusivity_value
+  []
+  [time]
+    type = TimeDerivative
+    variable = u
+  []
+[]
+
+[BCs]
+  [left]
+    type = DirichletBC
+    variable = u
+    boundary = left
+    value = 3.0537
+  []
+  [right]
+    type = DirichletBC
+    variable = u
+    boundary = right
+    value = 0
+  []
+[]
+
+# Used while obtaining steady-state solution
+#
+# [VectorPostprocessors]
+#   [line]
+#     type = LineValueSampler
+#     start_point = '0 0 0'
+#     end_point = '99e-6 0 0'
+#     num_points = 199
+#     sort_by = 'x'
+#     variable = u
+#   []
+# []
+
+[Postprocessors]
+  [conc_point1]
+    type = PointValue
+    variable = u
+    point = '48.75e-6 0 0'
+  []
+[]
+
+[Executioner]
+  type = Transient
+  # end_time = 5000 # for obtaining steady-state solution
+  # dtmax = 200.0 # for obtaining steady-state solution
+  end_time = 50
+  dtmax = 2.0
+  solve_type = NEWTON
+  petsc_options_iname = '-pc_type -pc_hypre_type'
+  petsc_options_value = 'hypre boomeramg'
+  scheme = 'crank-nicolson'
+
+  [TimeStepper]
+    type = IterationAdaptiveDT
+    dt = 1.0
+    optimal_iterations = 4
+  []
+[]
+
+[Outputs]
+  exodus = true
+  [csv]
+    # interval = 10 # used while obtaining steady-state solution
+    type = CSV
+  []
+  perf_graph = true
+[]