Issues with multi-material Henry's law

User Support
1

Having some issues running a multi-material simulation with a salt-metal interface.

I have a model here with half salt, half metal, with a time-dependent fixed concentration at the salt boundary and insulated on the end of the metal.

It seems the hydrogen enters the metal, but then cant leave despite there being no traps.

Here is a gif to illustrate:

[animate output image]

here is the script to reproduce (runs on latest version of fenicsx branch commit - # d308842)

import festim as F
from dolfinx.log import set_log_level, LogLevel
import numpy as np

my_model = F.HydrogenTransportProblemDiscontinuous()

vertices = np.concatenate(
    [
        np.linspace(0, 0.03, num=300),
        np.linspace(0.03, 0.06, num=300),
    ]
)

my_mesh = F.Mesh1D(vertices=vertices)
my_model.mesh = my_mesh

mat_flibe = F.Material(
    D_0=9.3e-7,
    E_D=0.435,
    K_S_0=4.76e22,
    E_K_S=0.363,
    solubility_law="henry",
)
mat_inconel = F.Material(
    D_0=1.75e-6,
    E_D=0.545,
    K_S_0=1.26e23,
    E_K_S=0.109,
    solubility_law="sievert",
)

vol_flibe = F.VolumeSubdomain1D(id=1, material=mat_flibe, borders=[0, 0.03])
vol_inconel = F.VolumeSubdomain1D(id=2, material=mat_inconel, borders=[0.03, 0.06])

inner_flibe_surface = F.SurfaceSubdomain1D(id=3, x=0)
outer_radius_surface = F.SurfaceSubdomain1D(id=4, x=0.06)

my_model.subdomains = [
    vol_inconel,
    vol_flibe,
    inner_flibe_surface,
    outer_radius_surface,
]

salt_metal_interface = F.Interface(
    id=5, subdomains=(vol_flibe, vol_inconel), penalty_term=1e25
)

my_model.interfaces = [salt_metal_interface]
my_model.surface_to_volume = {
    inner_flibe_surface: vol_flibe,
    outer_radius_surface: vol_inconel,
}

T = F.Species("T", mobile=True, subdomains=[vol_flibe, vol_inconel])
my_model.species = [T]

irradiation_time = 12 * 3600


def my_value(t):
    if t < irradiation_time:
        return 5e4
    else:
        return 0


my_model.boundary_conditions = [
    F.FixedConcentrationBC(
        subdomain=inner_flibe_surface, species=T, value=lambda t: my_value(t)
    ),
]

my_model.temperature = 923  # 650 degC

dt = F.Stepsize(
    100,
    growth_factor=1.1,
    cutback_factor=0.9,
    target_nb_iterations=4,
    milestones=[irradiation_time],
)

my_model.settings = F.Settings(
    transient=True,
    atol=1e-10,
    rtol=1e-12,
    final_time=30 * 24 * 3600,
    stepsize=dt,
)

my_model.exports = [
    F.VTXSpeciesExport(
        "results/tritium_concentration_flibe_1D.bp", field=T, subdomain=vol_flibe
    ),
    F.VTXSpeciesExport(
        "results/tritium_concentration_inconel_1D.bp", field=T, subdomain=vol_inconel
    ),
]

set_log_level(LogLevel.INFO)

my_model.initialise()

my_model.run()
2

Here it seems that it’s not going down because of solving in zero iterations