Merge pull request #1 from LIBRA-project/parametrisation

Parametrisation

Author: RemDelaporteMathurin

README.md

@@ -1 +1,21 @@
-# libra_sparging
+# libra_sparging
+
+To run the model:
+
+1. Create the right conda environment
+
+> [!NOTE]
+> Requires and conda to be installed
+
+> [!NOTE]
+> This uses `dolfinx` which doesn't run on Linux. For windows users, consider using Windows Subsystem for Linux (WSL)
+
+```
+conda env create -f environment.yml
+conda activate libra_sparging
+```
+
+
+```
+python model.py
+```

model.py

@@ -11,7 +11,36 @@
 
 log.set_log_level(log.LogLevel.INFO)
 
-mesh = dolfinx.mesh.create_interval(MPI.COMM_WORLD, 1000, points=[0, 1])
+avogadro_number = 6.022e23  # 1/mol
+
+tank_height = 1  # m
+tank_diameter = 0.5  # m
+K_S = 1e20 / avogadro_number  # solubility of breeder mol/m3/Pa
+diffusivity = 1e-9  # breeder diffusivity  m2/s
+
+
+source_term = 0.001  # mol/m3/s generation term
+
+u_b = 0.3  # m/s bubble velocity
+d_b = 0.002  # m bubble diameter
+
+P = 151988  # total gas pressure  # TODO should be 7 PSIG - differential at the top
+h_l = (
+    (diffusivity * u_b) / (ufl.pi * d_b)
+) ** 0.5  # mass transport coefficient Higbie penetration model
+
+R = 8.314  # J/mol/K
+T = 900  # K temperature
+epsilon_g = 0.03  # gas void fraction  # TODO from correlations
+epsilon_l = 1 - epsilon_g  # liquid void fraction
+a = 6 * epsilon_g / d_b  # specific interfacial area
+
+# FIXME is this homogeneous?
+E_g = 0.2 * tank_diameter**2 * u_b  # gas phase diffusivity (dispersion coefficient)
+E_l = diffusivity  # liquid phase diffusivity  # FIXME
+
+# MESH AND FUNCTION SPACES
+mesh = dolfinx.mesh.create_interval(MPI.COMM_WORLD, 10000, points=[0, tank_height])
 fdim = mesh.topology.dim - 1
 cg_el = basix.ufl.element("Lagrange", mesh.basix_cell(), degree=1, shape=(2,))
 
@@ -27,47 +56,49 @@
 
 dt = 0.2
 
-vel_x = 1.0
+vel_x = u_b
 vel = dolfinx.fem.Constant(mesh, PETSc.ScalarType([vel_x]))
 
-K_S = 1  # solubility of breeder
-D = 0.001
-P = 1  # total gas pressure
-h_l = 0.1  # mass transport coefficient
+
 EPS = 1e-16
 
-gen = dolfinx.fem.Constant(mesh, PETSc.ScalarType(1.0))  # generation term (neutrons)
+gen = dolfinx.fem.Constant(
+    mesh, PETSc.ScalarType(source_term)
+)  # generation term (neutrons)
+
+
+# VARIATIONAL FORMULATION
 
 # mass transfer rate
-J = h_l * (c_T - K_S * (P * y_T2 + EPS) ** 0.5)
+J = a * h_l * (c_T - K_S * (P * y_T2 + EPS))
 
 F = 0  # variational formulation
 
 # transient terms
-F += ((c_T - c_T_n) / dt) * v_c * ufl.dx + ((y_T2 - y_T2_n) / dt) * v_y * ufl.dx
-
+F += epsilon_l * ((c_T - c_T_n) / dt) * v_c * ufl.dx
+F += epsilon_g * 1 / (R * T) * (P * (y_T2 - y_T2_n) / dt) * v_y * ufl.dx
 
 # diffusion/dispersion terms
-F += D * ufl.dot(ufl.grad(c_T), ufl.grad(v_c)) * ufl.dx
-F += ufl.dot(ufl.grad(y_T2), ufl.grad(v_y)) * ufl.dx
+F += epsilon_l * E_l * ufl.dot(ufl.grad(c_T), ufl.grad(v_c)) * ufl.dx
+F += epsilon_g * E_g * ufl.dot(ufl.grad(P * y_T2), ufl.grad(v_y)) * ufl.dx
 
 
 # mass exchange (coupling term)
-F += J * v_c * ufl.dx - 0.5 * J * v_y * ufl.dx
+F += J * v_c * ufl.dx - J * v_y * ufl.dx
 
 # Generation term in the breeder
 F += -gen * v_c * ufl.dx
 
 # advection of gas
-F += ufl.inner(ufl.dot(ufl.grad(y_T2), vel), v_y) * ufl.dx
+F += 1 / (R * T) * ufl.inner(ufl.dot(ufl.grad(P * y_T2), vel), v_y) * ufl.dx
 
 
 # BOUNDARY CONDITIONS
 gas_inlet_facets = dolfinx.mesh.locate_entities_boundary(
     mesh, fdim, lambda x: np.isclose(x[0], 0.0)
 )
 gas_outlet_facets = dolfinx.mesh.locate_entities_boundary(
-    mesh, fdim, lambda x: np.isclose(x[0], 1.0)
+    mesh, fdim, lambda x: np.isclose(x[0], tank_height)
 )
 bc1 = dolfinx.fem.dirichletbc(
     dolfinx.fem.Constant(mesh, 0.0),