Ability to export trap densities to xdmf and derived quantites

[jhdark]

I think it could be fine if we have the ability to export the trap densities as both .xdmf and in derived quantities for analysing density variation with a trap creation model

[RemDelaporteMathurin]

@jhdark This can already be done. See:

import FESTIM as F
import fenics as f

class CustomSim(F.Simulation):
    def iterate(self):  # overide iterate
        super().iterate()

        # export trap density
        trap = self.traps[0]  # change accordingly
        density = trap.density

        # if density is a fenics.Expression, project it first
        density_as_function = f.project(density, self.V_DG1)

        # write function to XDMF
        density_file.write_checkpoint(
            density_as_function, 'density', self.t,
            f.XDMFFile.Encoding.HDF5,
            append=append)

        append = True

append = False  # initialise at False
density_file = f.XDMFile("trap_density.xdmf")

# run simulation
my_sim = CustomSim(...)

my_sim.initialise()
my_sim.run()

[RemDelaporteMathurin]

@jhdark this is much better tbf

import FESTIM as F
import fenics as f

class DensityXDMF(F.XDMFExport):
    def __init__(self, trap, **kwargs) -> None:
        """Inits DensityXDMF

        Args:
            trap (FESTIM.Trap): the trap we want to export (density)
        """
        super().__init__(
            field="1", **kwargs
        )  # field is "1" just to make the code not crash

        self.trap = trap

    def write(self, t):
        functionspace = self.function.function_space().collapse()
        density_as_function = f.project(self.trap.density[0], functionspace)
        self.function = density_as_function
        super().write(t)

# Create a Simulation
my_sim = F.Simulation()

mesh = F.MeshFromVertices([0, 1, 2, 3, 4])
my_sim.mesh = mesh

my_sim.materials = F.Materials([F.Material(1, 1, 0)])

trap_1 = F.Trap(0, 0, 0, 0, [1], 1 + F.x + F.t)
my_sim.traps = F.Traps([trap_1])

# use our custom XDMFExport class
my_sim.exports = F.Exports(
    [DensityXDMF(trap=trap_1, label="density1", folder="results", checkpoint=False)]
)

my_sim.boundary_conditions = []
my_sim.T = F.Temperature(100)

my_sim.settings = F.Settings(1e-10, 1e-10, final_time=10)
my_sim.dt = F.Stepsize(1)

# run simulation
my_sim.initialise()
my_sim.run()

[RemDelaporteMathurin]

@jhdark this now works on dev with NeutronInducedTrap

import FESTIM as F
import fenics as f
import numpy as np

class DensityXDMF(F.XDMFExport):
    def __init__(self, trap, **kwargs) -> None:
        """Inits DensityXDMF

        Args:
            trap (FESTIM.Trap): the trap we want to export (density)
        """
        super().__init__(
            field="1", **kwargs
        )  # field is "1" just to make the code not crash

        self.trap = trap

    def write(self, t):
        functionspace = self.function.function_space().collapse()
        density_as_function = f.project(self.trap.density[0], functionspace)
        self.function = density_as_function
        super().write(t)

# Create a Simulation
my_sim = F.Simulation()

mesh = F.MeshFromVertices(np.linspace(0, 10, num=100))
my_sim.mesh = mesh

my_sim.materials = F.Materials([F.Material(1, 1, 0)])

trap_1 = F.NeutronInducedTrap(
    0, 0, 0, 0, materials=[1], phi=F.x**2, K=1, n_max=5, A_0=0, E_A=0
)
my_sim.traps = F.Traps([trap_1])

# use our custom XDMFExport class
my_sim.exports = F.Exports(
    [DensityXDMF(trap=trap_1, label="density1", folder="results", checkpoint=False)]
)

my_sim.boundary_conditions = []
my_sim.T = F.Temperature(100)

my_sim.settings = F.Settings(1e-10, 1e-10, final_time=100)
my_sim.dt = F.Stepsize(1)

# run simulation
my_sim.initialise()
my_sim.run()

[jhdark]

Very nice

[RemDelaporteMathurin]

@jhdark do you think this should go in the next release v0.10 ?

[jhdark]

As in make it an actual class within the exports?

[RemDelaporteMathurin]

Yeah as in implement this class in FESTIM