Coupled transient heat transfer and hydrogen transport models

[jhdark]

Proposed changes

fix for #945

With this update, users can run coupled transient heat transfer and hydrogen transport simulations

This PR introduces a new class: CoupledHeatTransferHydrogenTransport which takes arguments of a heat_problem and a hydrogen_problem.

An example of how these changes would be used:

import numpy as np

import festim as F

my_mesh = F.Mesh1D(
    vertices=np.linspace(0, 1, 1000),
)
my_mat = F.Material(
    D_0=1, E_D=0, thermal_conductivity=1, name="my_mat", density=1, heat_capacity=1
)
vol_bulk = F.VolumeSubdomain1D(id=1, borders=[0, 1], material=my_mat)
surf_left = F.SurfaceSubdomain1D(id=2, x=0)
surf_right = F.SurfaceSubdomain1D(id=3, x=1)
H = F.Species("H", mobile=True, subdomains=[vol_bulk])

concentration_func = lambda T: 10 * T
temperature_func = lambda t: t + 5

my_heat_problem = F.HeatTransferProblem(
    mesh=my_mesh,
    subdomains=[
        vol_bulk,
        surf_right,
        surf_left,
    ],
    boundary_conditions=[
        F.FixedTemperatureBC(subdomain=surf_left, value=temperature_func),
        F.FixedTemperatureBC(subdomain=surf_right, value=0),
    ],
    exports=[
        F.VTXTemperatureExport(filename="heat_problem_output_1D"),
    ],
    settings=F.Settings(
        atol=1e-10, rtol=1e-10, transient=True, final_time=10, stepsize=1.0
    ),
)

my_hydrogen_problem = F.HydrogenTransportProblem(
    mesh=my_mesh,
    subdomains=[vol_bulk, surf_right, surf_left],
    species=[H],
    boundary_conditions=[
        F.FixedConcentrationBC(
            subdomain=surf_right, value=concentration_func, species=H
        ),
        F.FixedConcentrationBC(
            subdomain=surf_left, value=concentration_func, species=H
        ),
    ],
    settings=F.Settings(atol=1e-10, rtol=1e-10, transient=True, final_time=10),
    exports=[
        F.VTXSpeciesExport(
            filename="hydrogen_problem_output_1D", subdomain=vol_bulk, field=H
        ),
    ],
)

my_hydrogen_problem.settings.stepsize = F.Stepsize(
    initial_value=1.5, growth_factor=1.5, target_nb_iterations=30
)

my_coupled_problem = F.CoupledHeatTransferHydrogenTransport(
    hydrogen_problem=my_hydrogen_problem,
    heat_problem=my_heat_problem,
)

my_coupled_problem.initialise()
my_coupled_problem.run()

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

If this is a relatively large or complex change, kick off the discussion by explaining why you chose the solution you did and what alternatives you considered, etc...

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 96.13%. Comparing base (01ab7f9) to head (a9c100f).

Additional details and impacted files

@@             Coverage Diff             @@
##           fenicsx     #948      +/-   ##
===========================================
+ Coverage    96.02%   96.13%   +0.11%     
===========================================
  Files           45       46       +1     
  Lines         2490     2564      +74     
===========================================
+ Hits          2391     2465      +74     
  Misses          99       99              

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

[RemDelaporteMathurin]

first pass on my side

[RemDelaporteMathurin]

Suggested change

	class CoupledtTransientHeatTransferHydrogenTransport:
	class CoupledTransientHeatTransferHydrogenTransport:

[RemDelaporteMathurin]

Suggested change

	exact_mobile_intial_cond = lambda x: 2 * x[0] ** 2
	exact_trapped_intial_cond = lambda x: 4 * x[0] ** 2
	exact_mobile_intial_cond = lambda x: exact_mobile_solution(x, t=0)
	exact_trapped_intial_cond = lambda x: exact_trapped_solution(x, t=0)

[RemDelaporteMathurin]

Suggested change

	exact_mobile = lambda x: 2 * x[0] ** 2 + 15 * final_time
	exact_mobile = lambda x: exact_mobile_solution(x, t=final_time)

same for the two others

[RemDelaporteMathurin]

Suggested change

	# TEST ensure L2 error below 2e07
	# TEST ensure L2 error below 2e-7

[RemDelaporteMathurin]

Could we test the whole function by computing a global error?