Fix Surface Reaction BC in Discontinuous Case

[kaelyndunnell]

Proposed changes

Fixes surface reaction BC in discontinuous cases, which was violating flux conservation in tests.

Types of changes

What types of changes does your code introduce to FESTIM?

• Bugfix (non-breaking change which fixes an issue)
• New feature (non-breaking change which adds functionality)
• Breaking change (fix or feature that would cause existing functionality to not work as expected)
• Code refactoring
• Documentation Update (if none of the other choices apply)
• New tests

Checklist

• Black formatted
• Unit tests pass locally with my changes
• I have added tests that prove my fix is effective or that my feature works
• I have added necessary documentation (if appropriate)

Further comments

Here is the mwe.py that exposed the bug (thanks @RemDelaporteMathurin):
``
model_probe = F.HydrogenTransportProblemDiscontinuous()

probe_thick = 1
breeder_thick = 1

probe = F.Material(
    D_0=1,
    E_D=0,
    K_S_0=1,
    E_K_S=0,
)
breeder = F.Material(D_0=1, E_D=0, K_S_0=2, E_K_S=0)

volume_breeder = F.VolumeSubdomain1D(
    id=1, borders=[0, breeder_thick], material=breeder
)
volume_probe = F.VolumeSubdomain1D(
    id=2,
    borders=[breeder_thick, breeder_thick + probe_thick],
    material=probe,
)

vertices_left = np.linspace(0, breeder_thick, num=1000, endpoint=False)
vertices_right = np.linspace(breeder_thick, breeder_thick + probe_thick, num=1000)
vertices = np.concatenate([vertices_left, vertices_right])
model_probe.mesh = F.Mesh1D(vertices)

inlet = F.SurfaceSubdomain1D(id=3, x=0)
outlet = F.SurfaceSubdomain1D(id=4, x=breeder_thick + probe_thick)

model_probe.subdomains = [
    volume_breeder,
    volume_probe,
    inlet,
    outlet,
]

model_probe.surface_to_volume = {
    inlet: volume_breeder,
    outlet: volume_probe,
}

model_probe.method_interface = F.InterfaceMethod.nitsche
model_probe.interfaces = [
    F.Interface(id=5, subdomains=[volume_breeder, volume_probe], penalty_term=10),
]

H = F.Species("H", subdomains=model_probe.volume_subdomains)
model_probe.species = [H]

model_probe.temperature = 600

model_probe.boundary_conditions = [
    F.FixedConcentrationBC(
        subdomain=inlet,
        value=1,
        species=H,
    ),
    F.SurfaceReactionBC(
        reactant=[H, H],
        gas_pressure=0,
        k_r0=1,
        E_kr=0,
        k_d0=0,
        E_kd=0,
        subdomain=outlet,
    ),
]
permeation_flux = F.SurfaceFlux(field=H, surface=outlet)
inlet_flux = F.SurfaceFlux(field=H, surface=inlet)
model_probe.exports = [permeation_flux, inlet_flux]

model_probe.settings = F.Settings(
    atol=1e-10,
    rtol=1e-09,
    transient=False,
)

model_probe.initialise()
model_probe.run()

assert np.isclose(permeation_flux.data[-1], -inlet_flux.data[-1], rtol=1e-5)``

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 94.90%. Comparing base (ab51de0) to head (63d0311).
⚠️ Report is 4 commits behind head on fenicsx.

Additional details and impacted files

@@             Coverage Diff             @@
##           fenicsx    #1056      +/-   ##
===========================================
+ Coverage    94.86%   94.90%   +0.03%     
===========================================
  Files           44       44              
  Lines         3099     3099              
===========================================
+ Hits          2940     2941       +1     
+ Misses         159      158       -1     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[jhdark]

Great catch!

[RemDelaporteMathurin]

Good catch @kaelyndunnell !!!