Merge pull request festim-dev#450 from RemDelaporteMathurin/dev

Release 0.9

.github/pull_request_template.md

@@ -5,7 +5,7 @@ Describe the big picture of your changes here to communicate to the maintainers
 ## Types of changes
 
 What types of changes does your code introduce to FESTIM?
-_Put an `x` in the boxes that apply_
+<!--Put an `x` in the boxes that apply-->
 
 - [ ] Bugfix (non-breaking change which fixes an issue)
 - [ ] New feature (non-breaking change which adds functionality)
@@ -16,9 +16,9 @@ _Put an `x` in the boxes that apply_
 
 ## Checklist
 
-_Put an `x` in the boxes that apply. You can also fill these out after creating the PR. If you're unsure about any of them, don't hesitate to ask. We're here to help! This is simply a reminder of what we are going to look for before merging your code._
+<!--Put an `x` in the boxes that apply. You can also fill these out after creating the PR. If you're unsure about any of them, don't hesitate to ask. We're here to help! This is simply a reminder of what we are going to look for before merging your code.-->
 
-- [ ] Pep8 applied
+- [ ] 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)

.github/workflows/code_formatting.yml

@@ -0,0 +1,10 @@
+name: Lint
+
+on: [pull_request]
+
+jobs:
+  lint:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v2
+      - uses: psf/black@stable

CITATION.cff

@@ -0,0 +1,35 @@
+cff-version: 1.2.0
+title: FESTIM
+message: 'If you use this software, please cite it as below.'
+type: software
+authors:
+  - given-names: Rémi
+    family-names: Delaporte-Mathurin
+    email: [email protected]
+    orcid: 'https://orcid.org/0000-0003-1064-8882'
+  - given-names: James
+    family-names: Dark
+
+repository-code: 'https://github.com/RemDelaporteMathurin/FESTIM'
+version: 0.8.1
+preferred-citation:
+  authors:
+    - given-names: Rémi
+      family-names: Delaporte-Mathurin
+    - given-names: Etienne
+      family-names: Hodille
+    - given-names: Jonathan
+      family-names: Mougenot
+    - given-names: Yann
+      family-names: Charles
+    - given-names: Christian
+      family-names: Grisolia
+  title: "Finite element analysis of hydrogen retention in ITER plasma facing components using FESTIM"
+  type: article
+  year: 2019
+  doi: 10.1016/j.nme.2019.100709
+  url: "http://www.sciencedirect.com/science/article/pii/S2352179119300547"
+  issn: 2352-1791
+  volume: 21
+  pages: 100709
+  journal: "Nuclear Materials and Energy"

FESTIM/__init__.py

@@ -1,11 +1,18 @@
 import sympy as sp
-x, y, z, t = sp.symbols('x[0] x[1] x[2] t')
+
+x, y, z, t = sp.symbols("x[0] x[1] x[2] t")
 R = 8.314462618  # Gas constant J.mol-1.K-1
 k_B = 8.6173303e-5  # Boltzmann constant eV.K-1
 
-from .helpers import update_expressions, help_key, \
-    parameters_helper, read_parameters, kJmol_to_eV, \
-    extract_xdmf_labels, extract_xdmf_times
+from .helpers import (
+    update_expressions,
+    kJmol_to_eV,
+    extract_xdmf_labels,
+    extract_xdmf_times,
+    as_constant,
+    as_expression,
+    as_constant_or_expression,
+)
 
 from .meshing.mesh import Mesh
 from .meshing.mesh_1d import Mesh1D
@@ -17,10 +24,14 @@
 from .temperature.temperature_solver import HeatTransferProblem
 
 from .boundary_conditions.boundary_condition import BoundaryCondition
-from .boundary_conditions.dirichlets.dirichlet_bc import DirichletBC, \
-    BoundaryConditionTheta, BoundaryConditionExpression
+from .boundary_conditions.dirichlets.dirichlet_bc import (
+    DirichletBC,
+    BoundaryConditionTheta,
+    BoundaryConditionExpression,
+)
 from .boundary_conditions.dirichlets.dc_imp import ImplantationDirichlet
 from .boundary_conditions.dirichlets.sieverts_bc import SievertsBC
+from .boundary_conditions.dirichlets.henrys_bc import HenrysBC
 from .boundary_conditions.dirichlets.custom_dc import CustomDirichlet
 
 from .boundary_conditions.fluxes.flux_bc import FluxBC
@@ -30,7 +41,6 @@
 
 from .exports.exports import Exports
 from .exports.export import Export
-from .exports.error import Error
 from .exports.xdmf_export import XDMFExport, XDMFExports
 
 from .exports.derived_quantities.derived_quantity import DerivedQuantity
@@ -54,8 +64,8 @@
 from .sources.source import Source
 from .sources.source_implantation_flux import ImplantationFlux
 
-from .materials.materials import Materials
 from .materials.material import Material
+from .materials.materials import Materials
 
 from .concentration.concentration import Concentration
 from .initial_condition import InitialCondition
@@ -64,9 +74,10 @@
 
 from .concentration.traps.trap import Trap
 from .concentration.traps.traps import Traps
+from .concentration.traps.extrinsic_trap import ExtrinsicTrapBase
 from .concentration.traps.extrinsic_trap import ExtrinsicTrap
 from .concentration.traps.neutron_induced_trap import NeutronInducedTrap
 
 from .h_transport_problem import HTransportProblem
 
-from .generic_simulation import Simulation, run
+from .generic_simulation import Simulation

FESTIM/boundary_conditions/boundary_condition.py

@@ -4,12 +4,12 @@
 
 
 class BoundaryCondition:
-    def __init__(self, surfaces, component=0) -> None:
+    def __init__(self, surfaces, field=0) -> None:
 
         if not isinstance(surfaces, list):
             surfaces = [surfaces]
         self.surfaces = surfaces
 
-        self.component = component
+        self.field = field
         self.expression = None
         self.sub_expressions = []

FESTIM/boundary_conditions/dirichlets/custom_dc.py

@@ -14,35 +14,34 @@ def fun(T, solute, param1):
         param1=2*FESTIM.x + FESTIM.t
     )
     """
-    def __init__(self, surfaces, function, component=0, **prms) -> None:
+
+    def __init__(self, surfaces, function, field=0, **prms) -> None:
         """Inits CustomDirichlet
 
         Args:
             surfaces (list or int): the surfaces of the BC
             function (callable): the custom function
-            component (int, optional): the field the boundary condition is
+            field (int, optional): the field the boundary condition is
                 applied to. Defaults to 0.
         """
-        super().__init__(surfaces, component=component)
+        super().__init__(surfaces, field=field)
         self.function = function
         self.prms = prms
         self.convert_prms()
 
     def create_expression(self, T):
         value_BC = BoundaryConditionExpression(
-            T, self.function,
+            T,
+            self.function,
             **self.prms,
         )
         self.expression = value_BC
         self.sub_expressions = self.prms.values()
 
     def convert_prms(self):
-        """Creates Expressions or Constant for all parameters
-        """
+        """Creates Expressions or Constant for all parameters"""
         for key, value in self.prms.items():
             if isinstance(value, (int, float)):
                 self.prms[key] = f.Constant(value)
             else:
-                self.prms[key] = f.Expression(
-                    sp.printing.ccode(value),
-                    t=0, degree=1)
+                self.prms[key] = f.Expression(sp.printing.ccode(value), t=0, degree=1)

FESTIM/boundary_conditions/dirichlets/dc_imp.py

@@ -4,11 +4,11 @@
 
 
 def dc_imp(T, phi, R_p, D_0, E_D, Kr_0=None, E_Kr=None):
-    D = D_0*f.exp(-E_D/k_B/T)
-    value = phi*R_p/D
+    D = D_0 * f.exp(-E_D / k_B / T)
+    value = phi * R_p / D
     if Kr_0 is not None:
-        Kr = Kr_0*f.exp(-E_Kr/k_B/T)
-        value += (phi/Kr)**0.5
+        Kr = Kr_0 * f.exp(-E_Kr / k_B / T)
+        value += (phi / Kr) ** 0.5
 
     return value
 
@@ -22,6 +22,7 @@ class ImplantationDirichlet(DirichletBC):
     c = phi*R_p/D + (phi/Kr)**0.5
 
     """
+
     def __init__(self, surfaces, phi, R_p, D_0, E_D, Kr_0=None, E_Kr=None) -> None:
         """Inits ImplantationDirichlet
 
@@ -37,7 +38,7 @@ def __init__(self, surfaces, phi, R_p, D_0, E_D, Kr_0=None, E_Kr=None) -> None:
             E_Kr (float, optional): recombination coefficient activation
                 energy (eV). Defaults to None.
         """
-        super().__init__(surfaces, component=0)
+        super().__init__(surfaces, field=0)
         self.phi = phi
         self.R_p = R_p
         self.D_0 = D_0
@@ -46,18 +47,19 @@ def __init__(self, surfaces, phi, R_p, D_0, E_D, Kr_0=None, E_Kr=None) -> None:
         self.E_Kr = E_Kr
 
     def create_expression(self, T):
-        phi = f.Expression(
-            sp.printing.ccode(self.phi),
-            t=0, degree=1)
-        R_p = f.Expression(
-            sp.printing.ccode(self.R_p),
-            t=0, degree=1)
+        phi = f.Expression(sp.printing.ccode(self.phi), t=0, degree=1)
+        R_p = f.Expression(sp.printing.ccode(self.R_p), t=0, degree=1)
         sub_expressions = [phi, R_p]
 
         value_BC = BoundaryConditionExpression(
-            T, dc_imp,
-            phi=phi, R_p=R_p, D_0=self.D_0, E_D=self.E_D,
-            Kr_0=self.Kr_0, E_Kr=self.E_Kr,
-            )
+            T,
+            dc_imp,
+            phi=phi,
+            R_p=R_p,
+            D_0=self.D_0,
+            E_D=self.E_D,
+            Kr_0=self.Kr_0,
+            E_Kr=self.E_Kr,
+        )
         self.expression = value_BC
         self.sub_expressions = sub_expressions

FESTIM/boundary_conditions/dirichlets/dirichlet_bc.py

@@ -4,19 +4,19 @@
 
 
 class DirichletBC(BoundaryCondition):
-    """Class to enforce the solution on boundaries.
-    """
-    def __init__(self, surfaces, value=None, component=0) -> None:
+    """Class to enforce the solution on boundaries."""
+
+    def __init__(self, surfaces, value=None, field=0) -> None:
         """Inits DirichletBC
 
         Args:
             surfaces (list or int): the surfaces of the BC
             value (float or sp.Expr, optional): the value of the boundary
                 condition. Defaults to None.
-            component (int, optional): the field the boundary condition is
+            field (int, optional): the field the boundary condition is
                 applied to. Defaults to 0.
         """
-        super().__init__(surfaces, component=component)
+        super().__init__(surfaces, field=field)
         self.value = value
         self.dirichlet_bc = []
 
@@ -47,14 +47,19 @@ def normalise_by_solubility(self, materials, volume_markers, T):
         self.sub_expressions.append(self.expression)
         # create modified BC based on solubility
         expression_BC = BoundaryConditionTheta(
-                            self.expression,
-                            materials,
-                            volume_markers, T)
+            self.expression, materials, volume_markers, T
+        )
         self.expression = expression_BC
 
     def create_dirichletbc(
-            self, V, T, surface_markers, chemical_pot=False, materials=None,
-            volume_markers=None):
+        self,
+        V,
+        T,
+        surface_markers,
+        chemical_pot=False,
+        materials=None,
+        volume_markers=None,
+    ):
         """creates a list of fenics.DirichletBC and stores it in
         self.dirichlet_bc
 
@@ -73,18 +78,17 @@ def create_dirichletbc(
         self.dirichlet_bc = []
         self.create_expression(T)
         # TODO: this should be more generic
-        mobile_components = [0, "0", "solute"]
-        if self.component in mobile_components and chemical_pot:
+        mobile_fields = [0, "0", "solute"]
+        if self.field in mobile_fields and chemical_pot:
             self.normalise_by_solubility(materials, volume_markers, T)
 
         # create a DirichletBC and add it to bcs
         if V.num_sub_spaces() == 0:
             funspace = V
-        else:  # if only one component, use subspace
-            funspace = V.sub(self.component)
+        else:  # if only one field, use subspace
+            funspace = V.sub(self.field)
         for surface in self.surfaces:
-            bci = f.DirichletBC(funspace, self.expression,
-                                surface_markers, surface)
+            bci = f.DirichletBC(funspace, self.expression, surface_markers, surface)
             self.dirichlet_bc.append(bci)
 
 
@@ -95,6 +99,7 @@ class BoundaryConditionTheta(f.UserExpression):
     Args:
         UserExpression (fenics.UserExpression):
     """
+
     def __init__(self, bci, materials, vm, T, **kwargs):
         """initialisation
 
@@ -116,18 +121,22 @@ def eval_cell(self, value, x, ufc_cell):
         cell = f.Cell(self._mesh, ufc_cell.index)
         subdomain_id = self._vm[cell]
         material = self._materials.find_material_from_id(subdomain_id)
-        # TODO this requires changes for Henry's law
         S_0 = material.S_0
         E_S = material.E_S
-        value[0] = self._bci(x)/(S_0*f.exp(-E_S/k_B/self._T(x)))
+        c = self._bci(x)
+        S = S_0 * f.exp(-E_S / k_B / self._T(x))
+        if material.solubility_law == "sievert":
+            value[0] = c / S
+        elif material.solubility_law == "henry":
+            value[0] = (c / S + f.DOLFIN_EPS) ** 0.5
 
     def value_shape(self):
         return ()
 
 
 class BoundaryConditionExpression(f.UserExpression):
     def __init__(self, T, eval_function, **kwargs):
-        """"[summary]"
+        """ "[summary]"
 
         Args:
             T (fenics.Function): the temperature