Skip to content

Commit

Permalink
Update ver-1jb documentation and clean up
Browse files Browse the repository at this point in the history
  • Loading branch information
@simopier
simopier committed Oct 15, 2024
1 parent ecea800 commit 14a8b9b
Show file tree
Hide file tree
Showing 3 changed files with 157 additions and 26 deletions.
48 changes: 34 additions & 14 deletions doc/content/verification_and_validation/ver-1jb.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,15 +13,15 @@ This verification case is an extension of[ver-1ja](ver-1ja.md), which tests the
In ver-1jb, however, tritium decay is coupled with trapping, which was verified in several verification cases, including [ver-1d](ver-1d.md).
As [ver-1ja](ver-1ja.md), ver-1jb is based on the case published in the TMAP7 V&V suite [!citep](ambrosek2008verification).
The model assumes pre-charging of an $l=1.5$ m long slab with tritium.
As opposed to [ver-1ja](ver-1ja.md), however, traps at $C_{trap} = 0.1$% atom fraction
and $E=4.2$ eV trap energy are distributed in a normal distribution centered at the mid-plane of the slab and standard deviation of $l/4$.
As opposed to [ver-1ja](ver-1ja.md), however, traps at $C_{trap} = 0.1$% atom fraction (with a material density based on tungsten defined as 6.34 $\times 10^{28}$ atoms/m$^3$)
and $E=4.2$ eV trap energy are distributed in a normal distribution centered at the mid-plane of the slab and standard deviation of $l/4$. The traps are initially filled with trapped tritium to 50% of trap concentration.

The evolution of the mobile tritium, trapped tritium, and helium concentration, i.e.,
$C_M$, $C_T$, and $C_{He}$, respectively, is governed by

\begin{equation}
\label{eqn:diffusion_mobile}
\frac{dC_M}{dt} = - \nabla D_T \nabla C_M - \text{trap\_per\_free} \cdot \left(\frac{dC_T}{dt} - k C_T \right) - k C_M,
\frac{dC_M}{dt} = - \nabla D_T \nabla C_M - \text{trap\_per\_free} \cdot \left(\frac{dC_T}{dt} + k C_T \right) - k C_M,
\end{equation}
\begin{equation}
\label{eqn:trapped_rate}
Expand Down Expand Up @@ -63,21 +63,26 @@ $k$ is defined as $k= 0.693/t_{1/2}=1.78199 \times 10^{-9}$ 1/s instead of $1.78
Note also that TMAP7 uses a temperature of $T = 273$ K instead of $T = 300$ K.


##### TBD second set of simulations 000000 + provide model parameter values + explain dimensionless simulations
The traps are initially filled to 50% of trap concentration and the initial mobile atom concentration is set to $C_M^0 = 1$ atom/m$^3$.
We define two different initial conditions for the mobile tritium concentration $C_M^0$:

- $C_M^0 = 1$ atoms/m$^3$, which is much lower than the initial trapped tritium concentration. This case corresponds to the TMAP7 case in [!citep](ambrosek2008verification).
- $C_M^0 = 1 \times 10^{25}$ atoms/m$^3$, which makes it equivalent to the initial trapped tritium concentration. This new case demonstrate TMAP8's ability to model the decay, trapping, detrapping, and diffusion of tritium better than when the concentration of mobile tritium is negligible, as in the other case.

To limit the computational challenges related to the orders of magnitude difference in trapped and mobile tritium in the first case, we make the concentration dimensionless by dividing them by $1 \times 10^{25}$ atoms/m$^3$ and setting $\text{trap\_per\_free} = 1 \times 10^{25}$ (-). In the second case, when the concentration are equivalent, $\text{trap\_per\_free} = 1$ (-).

## Analytical Solution

The total concentration of T in atoms/m$^3$, $C_{tot} = C_M + C_T$, at any given time is given by
The total inventory of T in atoms, $I_{tot} = I_M + I_T$ where $I_M$ and $I_T$ are the mobile and trapped tritium inventories, respectively, is given at any given time by

\begin{equation}
C_{tot} = C_{tot}^0 \exp(-kt),
I_{tot} = I_{tot}^0 \exp(-kt),
\end{equation}

where $C_{tot}^0 = C_M^0 + C_T^0$ atoms/m$^3$ is the initial total concentration of tritium.
Applying a mass balance over the system, the concentration of helium in atoms/m$^3$, $C_{He}$, is given by
where $I_{tot}^0 = I_M^0 + I_T^0$ atoms/m is the initial total inventory of tritium
($I_M^0$ and $I_T^0$ are the initial mobile and trapped tritium inventories).
Applying a mass balance over the system, the inventory of helium in atoms, $I_{He}$, is given by
\begin{equation}
C_{He} = C_{tot}^0 \left[1- \exp(-kt) \right].
I_{He} = I_{tot}^0 \left[1- \exp(-kt) \right].
\end{equation}

### Results
Expand All @@ -87,14 +92,14 @@ Applying a mass balance over the system, the concentration of helium in atoms/m$
[ver-1jb_results_comparison_analytical_time_evolution_1] shows the TMAP8 predictions and how they compare to the analytical solution
for the decay of tritium and associated growth of $^3$He in a distributed trap.
TMAP8 matches the analytical solution, with a root mean square percentage error
(RMSPE) of 3.96% and 1.02% for the $C_{tot}$ and $C_{He}$ concentration curves, respectively,
(RMSPE) of 2.85% and 0.92% for the $I_{tot}$ and $I_{He}$ concentration curves, respectively,
and can also provide the trapped and mobile tritium concentrations.

!media comparison_ver-1jb.py
image_name=ver-1jb_comparison_analytical_time_evolution.png
style=width:50%;margin-bottom:2%;margin-left:auto;margin-right:auto
id=ver-1jb_results_comparison_analytical_time_evolution_1
caption=Comparison of TMAP8 predictions against the analytical solution for the decay of tritium and associated growth of $^3$He in a distributed trap.
caption=Comparison of TMAP8 predictions against the analytical solution for the decay of tritium and associated growth of $^3$He in a distributed trap with a small concentration of mobile tritium compared to trapped tritium.

[ver-1jb_results_profile_1] shows the depth profile of the initial trapped atoms of tritium, the concentration of trapped atoms of
tritium after 45 years, and the distribution of $^3$He at the end of that time across the distributed trap as predicted by TMAP8.
Expand All @@ -107,12 +112,27 @@ including [ver-1d](ver-1d.md).
image_name=ver-1jb_profile.png
style=width:50%;margin-bottom:2%;margin-left:auto;margin-right:auto
id=ver-1jb_results_profile_1
caption=Concentration profiles of initially trapped tritium that decayed to $^3$He over 45 years.
caption=Concentration profiles of initially trapped tritium that decayed to $^3$He over 45 years with a small concentration of mobile tritium compared to trapped tritium.

#### With equivalent mobile and trapped tritium initial concentrations

[ver-1jb_results_comparison_analytical_time_evolution_2] and [ver-1jb_results_profile_2] show the results of the simulations when the initial concentrations of mobile and trapped tritium are equivalent.
[ver-1jb_results_comparison_analytical_time_evolution_2] shows TMAP8's time predictions of the inventories and how they compare to the analytical solution. The RMSPE values as low as in the previous case.

!media comparison_ver-1jb.py
image_name=ver-1jb_equivalent_concentrations_comparison_analytical_time_evolution.png
style=width:50%;margin-bottom:2%;margin-left:auto;margin-right:auto
id=ver-1jb_results_comparison_analytical_time_evolution_2
caption=Comparison of TMAP8 predictions against the analytical solution for the decay of tritium and associated growth of $^3$He in a distributed trap with equivalent initial concentrations of mobile and trapped tritium.

##### TBD 00000000000000000
[ver-1jb_results_profile_2] shows the depth profile of the initial trapped atoms of tritium, the concentration of mobile and trapped atoms of
tritium after 45 years, and the distribution of $^3$He at the end of that time across the distributed trap as predicted by TMAP8.

!media comparison_ver-1jb.py
image_name=ver-1jb_equivalent_concentrations_profile.png
style=width:50%;margin-bottom:2%;margin-left:auto;margin-right:auto
id=ver-1jb_results_profile_2
caption=Concentration profiles of initially trapped tritium that decayed to $^3$He over 45 years with equivalent initial concentrations of mobile and trapped tritium.

### Input file

Expand Down
133 changes: 122 additions & 11 deletions test/tests/ver-1jb/comparison_ver-1jb.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,10 @@
script_folder = os.path.dirname(__file__)
os.chdir(script_folder)

# ===============================================================================
# First Case - low mobile concentration
# ===============================================================================

# ===============================================================================
# Extract TMAP8 predictions - time evolution
if "/TMAP8/doc/" in script_folder: # if in documentation folder
Expand Down Expand Up @@ -50,19 +54,19 @@
gs = gridspec.GridSpec(1, 1)
ax = fig.add_subplot(gs[0])
tmap8_time_years = [t/conversion_years_to_s for t in tmap8_time]
ax.plot(tmap8_time_years,tmap8_tritium_trapped + tmap8_tritium_mobile,label=r"$C_{tot}$ - TMAP8",c='k', alpha=0.5)
ax.plot(tmap8_time_years,analytical_tritium,label=r"$C_{tot}$ - Analytical",c='k', linestyle='--')
ax.plot(tmap8_time_years,tmap8_tritium_trapped,label=r"$C_T$ - TMAP8",c='tab:blue', alpha=0.8, ls=':')
ax.plot(tmap8_time_years,tmap8_tritium_mobile,label=r"$C_M$ - TMAP8",c='tab:green', alpha=0.8, ls=':')
ax.plot(tmap8_time_years,tmap8_helium,label=r"$C_{He}$ - TMAP8",c='tab:red', alpha=0.5)
ax.plot(tmap8_time_years,analytical_helium,label=r"$C_{He}$ - Analytical",c='r', linestyle='--')
ax.plot(tmap8_time_years,tmap8_tritium_trapped + tmap8_tritium_mobile,label=r"$I_{tot}$ - TMAP8",c='k', alpha=0.5)
ax.plot(tmap8_time_years,analytical_tritium,label=r"$I_{tot}$ - Analytical",c='k', linestyle='--')
ax.plot(tmap8_time_years,tmap8_tritium_trapped,label=r"$I_T$ - TMAP8",c='tab:blue', alpha=0.8, ls=':')
ax.plot(tmap8_time_years,tmap8_tritium_mobile,label=r"$I_M$ - TMAP8",c='tab:green', alpha=0.8, ls=':')
ax.plot(tmap8_time_years,tmap8_helium,label=r"$I_{He}$ - TMAP8",c='tab:red', alpha=0.5)
ax.plot(tmap8_time_years,analytical_helium,label=r"$I_{He}$ - Analytical",c='r', linestyle='--')
# Root Mean Square Percentage Error calculations
RMSE_Ctot = np.linalg.norm(tmap8_tritium_trapped + tmap8_tritium_mobile-analytical_tritium)
RMSE_Ctot = np.linalg.norm(tmap8_tritium_trapped + tmap8_tritium_mobile - analytical_tritium)
RMSPE_Ctot = RMSE_Ctot*100/np.mean(analytical_tritium)
ax.text(33, 0.6e25, '$C_{tot}$ RMSPE = %.2f '%RMSPE_Ctot+'%',fontweight='bold',color = 'k')
ax.text(33, 0.6e25, '$I_{tot}$ RMSPE = %.2f '%RMSPE_Ctot+'%',fontweight='bold',color = 'k')
RMSE_CHe = np.linalg.norm(tmap8_helium-analytical_helium)
RMSPE_CHe = RMSE_CHe*100/np.mean(analytical_helium)
ax.text(52, 2.6e25, '$C_{He}$ RMSPE = %.2f '%RMSPE_CHe+'%',fontweight='bold',color = 'tab:red')
ax.text(52, 2.6e25, '$I_{He}$ RMSPE = %.2f '%RMSPE_CHe+'%',fontweight='bold',color = 'tab:red')

ax.set_xlabel(u'Time (years)')
ax.set_ylabel(r"Inventory (atoms)")
Expand All @@ -77,9 +81,9 @@
# ===============================================================================
# Extract TMAP8 predictions - concentration profile
if "/TMAP8/doc/" in script_folder: # if in documentation folder
csv_folder_profile = "../../../../test/tests/ver-1ja/gold/ver-1jb_profile_out_line_0080.csv"
csv_folder_profile = "../../../../test/tests/ver-1ja/gold/ver-1jb_profile_out_line_0068.csv"
else: # if in test folder
csv_folder_profile = "./gold/ver-1jb_profile_out_line_0080.csv"
csv_folder_profile = "./gold/ver-1jb_profile_out_line_0068.csv"
tmap8_solution = pd.read_csv(csv_folder_profile)
tmap8_position = tmap8_solution['x']
tmap8_tritium_mobile = tmap8_solution['tritium_mobile_concentration']
Expand All @@ -93,6 +97,8 @@
ax = fig.add_subplot(gs[0])
tmap8_tritium_trapped_init = trapping_sites_fraction_occupied_initial * density_material * trapping_sites_atomic_fraction_max * np.exp(-1/2*((tmap8_position-normal_center_position)/normal_standard_deviation)**2)
ax.plot(tmap8_position,tmap8_tritium_trapped_init,label=r"$C_T(t=0)$",c='tab:blue', alpha=0.5, ls=':')
tmap8_tritium_mobile_init = [tritium_mobile_concentration_initial]*len(tmap8_position)
ax.plot(tmap8_position,tmap8_tritium_mobile_init,label=r"$C_T(t=0)$",c='tab:green', alpha=0.5, ls=':')
ax.plot(tmap8_position,tmap8_tritium_trapped,label=r"$C_T(t=45 \text{years})$",c='tab:blue', alpha=1)
ax.plot(tmap8_position,tmap8_tritium_mobile,label=r"$C_M(t=45 \text{years})$",c='tab:green', alpha=1)
ax.plot(tmap8_position,tmap8_helium,label=r"$C_{He}(t=45 \text{years})$",c='tab:red', alpha=1)
Expand All @@ -106,3 +112,108 @@
ax.minorticks_on()
plt.savefig('ver-1jb_profile.png', bbox_inches='tight', dpi=300);
plt.close(fig)

# ===============================================================================
# Second Case - equivalent mobile and tritium concentrations
# ===============================================================================

# ===============================================================================
# Extract TMAP8 predictions - time evolution
if "/TMAP8/doc/" in script_folder: # if in documentation folder
csv_folder = "../../../../test/tests/ver-1ja/gold/ver-1jb_equivalent_concentrations_time_dependent_out.csv"
else: # if in test folder
csv_folder = "./gold/ver-1jb_equivalent_concentrations_time_dependent_out.csv"
tmap8_solution = pd.read_csv(csv_folder)
tmap8_time = tmap8_solution['time']
tmap8_tritium_mobile = tmap8_solution['tritium_mobile_inventory']
tmap8_tritium_trapped = tmap8_solution['tritium_trapped_inventory']
tmap8_helium = tmap8_solution['helium_inventory']

# ===============================================================================
# Case set up
slab_length = 1.5 # m
slab_height = 1 # m - assumed for integrations
slab_width = 1 # m - assumed for integrations
trapping_sites_atomic_fraction_max = 0.001 # at.frac.
trapping_sites_fraction_occupied_initial = 0.5 # (-)
density_material = 6.34e28 # atoms/m^3 # for tungsten
normal_center_position = slab_length/2 # m
normal_standard_deviation = slab_length/4 # m

# ===============================================================================
# Calculate the analytical solution
tritium_mobile_concentration_initial = 1e25 # atoms/m^3
tritium_trapped_inventory_initial = 2.8438315780556e25 * slab_height * slab_width # atoms -- integral of the normal distribution over the slab length * area
tritium_mobile_inventory_initial = tritium_mobile_concentration_initial * slab_length * slab_height * slab_width # atoms
half_life = 12.3232 # years
conversion_years_to_s = 365.25*24*60*60
half_life_s = half_life*conversion_years_to_s # s
decay_rate_constant = 0.693/half_life_s # 1/s

analytical_tritium = (tritium_mobile_inventory_initial + tritium_trapped_inventory_initial) * np.exp(- decay_rate_constant * tmap8_time)
analytical_helium = (tritium_mobile_inventory_initial + tritium_trapped_inventory_initial) * ( 1. - np.exp(- decay_rate_constant * tmap8_time))

# ===============================================================================
# Plot figure for verification of tritium decay - time evolution
fig = plt.figure(figsize=[6.5, 5.5])
gs = gridspec.GridSpec(1, 1)
ax = fig.add_subplot(gs[0])
tmap8_time_years = [t/conversion_years_to_s for t in tmap8_time]
ax.plot(tmap8_time_years,tmap8_tritium_trapped + tmap8_tritium_mobile,label=r"$I_{tot}$ - TMAP8",c='k', alpha=0.5)
ax.plot(tmap8_time_years,analytical_tritium,label=r"$I_{tot}$ - Analytical",c='k', linestyle='--')
ax.plot(tmap8_time_years,tmap8_tritium_trapped,label=r"$I_T$ - TMAP8",c='tab:blue', alpha=0.8, ls=':')
ax.plot(tmap8_time_years,tmap8_tritium_mobile,label=r"$I_M$ - TMAP8",c='tab:green', alpha=0.8, ls=':')
ax.plot(tmap8_time_years,tmap8_helium,label=r"$I_{He}$ - TMAP8",c='tab:red', alpha=0.5)
ax.plot(tmap8_time_years,analytical_helium,label=r"$I_{He}$ - Analytical",c='r', linestyle='--')
# Root Mean Square Percentage Error calculations
RMSE_Ctot = np.linalg.norm(tmap8_tritium_trapped + tmap8_tritium_mobile - analytical_tritium)
RMSPE_Ctot = RMSE_Ctot*100/np.mean(analytical_tritium)
ax.text(33, 0.75e25, '$I_{tot}$ RMSPE = %.2f '%RMSPE_Ctot+'%',fontweight='bold',color = 'k')
RMSE_CHe = np.linalg.norm(tmap8_helium-analytical_helium)
RMSPE_CHe = RMSE_CHe*100/np.mean(analytical_helium)
ax.text(52, 3.8e25, '$I_{He}$ RMSPE = %.2f '%RMSPE_CHe+'%',fontweight='bold',color = 'tab:red')

ax.set_xlabel(u'Time (years)')
ax.set_ylabel(r"Inventory (atoms)")
ax.legend(loc="best")
ax.set_xlim(left=0)
ax.set_ylim(bottom=0)
plt.grid(which='major', color='0.65', linestyle='--', alpha=0.3)
ax.minorticks_on()
plt.savefig('ver-1jb_equivalent_concentrations_comparison_analytical_time_evolution.png', bbox_inches='tight', dpi=300);
plt.close(fig)

# ===============================================================================
# Extract TMAP8 predictions - concentration profile
if "/TMAP8/doc/" in script_folder: # if in documentation folder
csv_folder_profile = "../../../../test/tests/ver-1ja/gold/ver-1jb_equivalent_concentrations_profile_out_line_0068.csv"
else: # if in test folder
csv_folder_profile = "./gold/ver-1jb_equivalent_concentrations_profile_out_line_0068.csv"
tmap8_solution = pd.read_csv(csv_folder_profile)
tmap8_position = tmap8_solution['x']
tmap8_tritium_mobile = tmap8_solution['tritium_mobile_concentration']
tmap8_tritium_trapped = tmap8_solution['tritium_trapped_concentration']
tmap8_helium = tmap8_solution['helium_concentration']

# ===============================================================================
# Plot figure for verification of tritium decay - concentration profile
fig = plt.figure(figsize=[6.5, 5.5])
gs = gridspec.GridSpec(1, 1)
ax = fig.add_subplot(gs[0])
tmap8_tritium_trapped_init = trapping_sites_fraction_occupied_initial * density_material * trapping_sites_atomic_fraction_max * np.exp(-1/2*((tmap8_position-normal_center_position)/normal_standard_deviation)**2)
ax.plot(tmap8_position,tmap8_tritium_trapped_init,label=r"$C_T(t=0)$",c='tab:blue', alpha=0.5, ls=':')
tmap8_tritium_mobile_init = [tritium_mobile_concentration_initial]*len(tmap8_position)
ax.plot(tmap8_position,tmap8_tritium_mobile_init,label=r"$C_T(t=0)$",c='tab:green', alpha=0.5, ls=':')
ax.plot(tmap8_position,tmap8_tritium_trapped,label=r"$C_T(t=45 \text{years})$",c='tab:blue', alpha=1)
ax.plot(tmap8_position,tmap8_tritium_mobile,label=r"$C_M(t=45 \text{years})$",c='tab:green', alpha=1)
ax.plot(tmap8_position,tmap8_helium,label=r"$C_{He}(t=45 \text{years})$",c='tab:red', alpha=1)
ax.set_xlabel(u'Position (m)')
ax.set_ylabel(r"Concentration (atoms/m$^3$)")
ax.legend(loc="best")
ax.set_xlim(left=0)
ax.set_xlim(right=slab_length)
ax.set_ylim(bottom=0)
plt.grid(which='major', color='0.65', linestyle='--', alpha=0.3)
ax.minorticks_on()
plt.savefig('ver-1jb_equivalent_concentrations_profile.png', bbox_inches='tight', dpi=300);
plt.close(fig)
2 changes: 1 addition & 1 deletion test/tests/ver-1jb/ver-1jb.i
View file
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@ trapping_sites_fraction_occupied_initial = 0.5 # (-)
normal_center_position = ${fparse slab_length/2} # m
normal_standard_deviation = ${fparse slab_length/4} # m
density_material = ${units 6.34e28 atoms/m^3} # for tungsten
density_scalar = ${units 1e25 atoms/m^3} # used to scale variables and use a dimentionless system
density_scalar = ${units 1e25 atoms/m^3} # used to scale variables and use a dimensionless system
temperature = ${units 300 K} # assumed (TMAP7's input file lists 273 K)
tritium_diffusivity_prefactor = ${units 1.58e-4 m^2/s} # from TMAP7 V&V input file
tritium_diffusivity_energy = ${units 308000.0 J/K} # from TMAP7 V&V input file
Expand Down

0 comments on commit 14a8b9b

Please sign in to comment.