fedm package#
Submodules#
fedm.file_io module#
- class fedm.file_io.Files#
Bases:
object
The Files class provides access to the files and directories used by FEDM. It should normally be accessed via the global instance fedm.files.
It has the following attributes:
- – file_input: Directory containing input files, default ‘./file_input’. Can be
assigned to a new path. Throws an error if the new path doesn’t exist.
- – output_folder_path: Directory which will contain output files, default ‘./output’.
Can be assigned to a new path. This will create a new directory if it doesn’t already exist.
- – error_file: File storing relative errors. Cannot be assigned, and will always
return ‘output_folder_path/relative error.log’. The first time the user accesses it per run, the file will be created/truncated.
- – model_log: File storing logs describing the run. Cannot be assigned, and will
always return ‘output_folder_path/model.log’. The first time the user accesses it per run, the file will be created/truncated.
Examples:
– Assign a new value to file_input (folder must exist) files.file_input = “my_file_input” – Get the current value of output_folder_path output = files.output_folder_path – Get error file (first access will truncate the file) error_file = files.error_file
- property error_file: Path#
- property file_input: Path#
- property model_log: Path#
- property output_folder_path: Path#
- fedm.file_io.decomment(lines: List[str]) str #
Removes comment at end of each line, denoted with ‘#’. If the line is blank, or the line starts with a ‘#’, returns an empty string. Works as a generator function. Code snippet addapted from: https://stackoverflow.com/questions/14158868/python-skip-comment-lines-marked-with-in-csv-dictreader
- fedm.file_io.file_output(t, t_old, t_out, step, t_out_list, step_list, file_type, output_file_list, particle_name, u_old, u_old1, unit='s')#
Writing value of desired variable to file. Value is calculated by linear interpolation. Input arguments are current time step, previous time step, current time step length, previous time step length, list of output times, list of steps (which needs have same size), list of output files, list of variable values in current and previous time steps, and (optional) units of time.
- fedm.file_io.flatten(input_list: List[List[Any]]) List[Any] #
Reduces 2D list to 1D list.
- fedm.file_io.flatten_float(input_list: List[List[Any]]) List[float] #
Reduces 2D list to 1D list and converts elements to float.
- fedm.file_io.log(log_type, log_file_name, *args)#
The function is used to log model data and its progress. The input arguments vary depending on the type of logging. For type = properties, input are gas, model, particle species names, particle mass, and charge. For type = conditions, input arguments are time step length, working voltage, pressure, gap length, gas number density, and gas temperature. For type = matrices, input arguments are gain, loss and power matrices. For type = initial_time, input argument is time. For type = time, input argument is time. For type = mesh, input argument is mesh.
- fedm.file_io.mesh_statistics(mesh: Mesh) None #
Returns mesh size and, maximum and minimum element size. Input is mesh.
- fedm.file_io.no_convert(x: Any) Any #
This utility function is used in functions that may optionally try to convert inputs to a given type. ‘no_convert’ may be used in place of ‘float’, ‘str’, etc to prevent the function from converting types.
- fedm.file_io.numpy_2d_array_to_str(x: ndarray) str #
- fedm.file_io.output_files(file_type: str, type_of_output: str, output_file_names: List[str]) List[Any] #
Creates desired number of output files.
- Parameters:
file_type (str) – ‘pvd’ or ‘xdmf’
type_of_output (str) – Name of folder where files should be saved
output_file_names – List of files to create
- fedm.file_io.print_time(t)#
Prints time.
- fedm.file_io.print_time_step(dt)#
Prints time step.
- fedm.file_io.rate_coefficient_file_names(path)#
Reads names of reaction coefficient files from “reacscheme.cfg” file. Input parameter is the path to the folder.
- fedm.file_io.reaction_matrices(path: str, species: List[str])#
Reads the reaction scheme from “reacscheme.cfg” file and creates power, loss and gain matrices.
- fedm.file_io.read_and_decomment(file_name: str) List[str] #
Reads file, returns list of strings. Comment lines and blank lines are removed.
- fedm.file_io.read_dependence(file_name: str)#
Reads dependence of rate coefficients from the corresponding file.
- fedm.file_io.read_dependences(file_names: List[str], zero_if_file_missing=False)#
Reads dependence of rate coefficients from a list of corresponding files.
- fedm.file_io.read_energy_loss(path)#
Reads energy loss values from “reacscheme.cfg” file. Input argument is the path to the folder.
- fedm.file_io.read_particle_properties(file_names: List[str], model: str)#
Reads particle properties from input file. Input are file names and model.
- fedm.file_io.read_rate_coefficients(rc_file_names: List[str], k_dependences: List[str])#
Reading rate coefficients from files. Input are file names and dependences.
- fedm.file_io.read_single_float(file_name, convert=<function no_convert>)#
Reads file containing only constant. Input parameter is file name.
- fedm.file_io.read_single_string(file_name)#
Reads file containing one column. Input parameter is file name.
- fedm.file_io.read_single_value(file_name, convert=<function no_convert>)#
Reads file containing only constant. Input parameter is file name.
- fedm.file_io.read_speclist(file_path)#
Function reads list of species from ‘speclist.cfg’ file and extracts species names and species properties filename.
- fedm.file_io.read_transport_coefficients(particle_names: List[str], transport_type: str, model: str)#
Reading transport coefricients from files. Input are particle names, type of transport coefficient (diffusion or mobility) and model.
- fedm.file_io.read_two_columns(file_name)#
Reads two column files. Input parameter is file name.
- fedm.file_io.truncate_file(path: Path) None #
Deletes contents of file, leaving it as an empty file. If the file does not exist, creates the empty file and any directories up to it.
fedm.functions module#
- fedm.functions.BoundaryGradient(var, zeroDomain, source_term, ds_extract, epsilon=8.854187817e-12)#
The function is an adaptation of the code snippet by D. Kamensky from https://fenicsproject.discourse.group/t/compute-gradient-of-scalar-field-on-boundarymesh/1172/2. It is used for the accurate calculation of the flux (in this case the electric field) across the specific boundary. Input parameters are variables whose gradient needs to be determined, marker of the whole domain except the boundary on which flux is calculated, the source term of the equation, list of ds of the specific boundaries (irrelevant boundaries should be marked as zero).
- fedm.functions.Boundary_flux(bc_type: str, equation_type: str, particle_type: str, sign: float, mu: Any, E: Any, normal: Any, u: Any, gamma: Any, v: Any, ds_temp: Measure, r: float = 0.15915494309189535, vth: float = 0.0, ref: float = 1.0, Ion_flux: float = 0.0)#
Function defines boundary conditions for different equations.
- Parameters:
flux' (Type of boundary condition. Options are 'zero) –
'Neumann' ('flux source' or) –
- equation_typestr
Choices are ‘reaction’, ‘diffusion-reaction’ or ‘drift-diffusion-reaction’
- particle_typestr
Choices are ‘Heavy’ or ‘electrons’.
- sign: float
Particle charge sign
- mu: df.Function
Mobility
- E
Electric field, Dolfin ComponentTensor
- normal
Dolfin FacetNormal
- u
Trial function
- gamma
Secondary electron emission coefficient
- v
Test function
- ds_tempdf.Measure
ds, surface element used to build integrals (??)
- r: float, default 0.5/pi
r coordinate
- vth: float, default 0.0
Thermal velocity
- ref: float, default 1.0
Reflection coefficient for specified particles species and boundary.
- Ion_flux: float, default 0.0
Flux of ions
- Returns:
df.Form – If bc_type is ‘flux source’ and there is a diffusion term in equation_type. Also if the user combines Neumann boundaries with a drift-diffusion-reaction equation.
float – Otherwise.
- Raises:
ValueError – If bc_type is not recognised. If bc_type is not ‘zero flux’, also raised if equation_type is not recognised. Furthermore, if equation_type is ‘drift-diffusion-reaction’, raises if particle_type is not recognised.
- class fedm.functions.CircleSubDomain(center_z: float, center_r: float, radius: float, gap_length: float, submesh: bool = False, tol: float = 1e-08)#
Bases:
SubDomain
- inside(self: dolfin.cpp.mesh.SubDomain, arg0: numpy.ndarray[numpy.float64[m, 1]], arg1: bool) bool #
- fedm.functions.Energy_Source_term(coupling: str, p_matrix: ndarray, l_matrix: ndarray, g_matrix: ndarray, k_coeffs: List[Function], u_loss: List[float], mean_energy: Any, N0: float, n: Any, Ei=0) List[Any] #
Defines energy source term for LMEA approximation. Function arguments are power, loss and gain matrices, rate coeffiients, energy losses for specific process, mean electron energy, gas number density and particle number density variable.
- exception fedm.functions.ErrorGreaterThanTTOL#
Bases:
Exception
- fedm.functions.Flux(sign, u, D, mu, E, grad_diffusion=True, logarithm_representation=True)#
Defines particle flux using drift-diffusion approximation. Input arguments are particle charge, number density, diffusion coefficient, mobility and electric field. The additional arguments are used to specify if the gradient of the diffusion should be considered and whether logarithmic representation is used.
- fedm.functions.Function_definition(function_space: FunctionSpace, function_type: str, eq_number: int = 1) List[Any] #
Defines list of desired function type (TrialFunction, TestFunction or ordinary Function). Input arguments are function space, type of desired function and number of equations, where the default value is one.
- fedm.functions.Function_space_list(number_of_equations: int, function_space: FunctionSpace) List[FunctionSpace] #
Defines list of function spaces. Input arguments are number of equations and function space.
- class fedm.functions.LineSubDomain(r_range: Tuple[float, float], z_range: Tuple[float, float])#
Bases:
SubDomain
- inside(self: dolfin.cpp.mesh.SubDomain, arg0: numpy.ndarray[numpy.float64[m, 1]], arg1: bool) bool #
- fedm.functions.Marking_boundaries(mesh: Mesh, boundaries: List[List[Any]], submesh: bool = False, gap_length: float = 0.01) MeshFunction #
Marking boundaries of a provided mesh. Currently, straight-line and circular boundaries are supported. First argument is the mesh, the second argument is a list of boundary properties (boundary type and coordinates).
- fedm.functions.Max(a, b)#
Returns maximum value of a and b.
- fedm.functions.Min(a, b)#
Returns minimum value of a and b.
- fedm.functions.Mixed_element_list(number_of_equations: int, element: FiniteElement) List[FiniteElement] #
Defines list of mixed elements. Input arguments are number of equations and element type.
- fedm.functions.Normal_vector(mesh: Mesh)#
- fedm.functions.Poisson_solver(A, L, b, bcs, u, solver_type='mumps', preconditioner='hypre_amg')#
- class fedm.functions.Problem(J, F, bcs)#
Bases:
NonlinearProblem
Nonlinear problem definition. Input parameters are F - weak form of the equation J - Jacobian and bcs -Dirichlet boundary conditions. The part of code is provided by user Nate Sime from FEniCS forum: https://fenicsproject.discourse.group/t/set-krylov-linear-solver-paramters-in-newton-solver/1070/3
- F(b, x)#
Linear form assembly
- J(A, x)#
Bilinear form assembly
- fedm.functions.Rate_coefficient_interpolation(status: str, dependences: List[str], k_coeffs: List[Function], kxs: List[Any], kys: List[Any], energy: Function, redfield: Function, Te: float = 300.0, Tgas: float = 300.0) None #
Function for linear interpolation of rate coefficients.
WARNING: When using the dependences ‘fun:Te,Tgas’ and ‘fun:Tgas’, this function will evaluate the corresponding contents of ‘ky’ using ‘eval’. It is the user’s responsibility to ensure no malicious code is injected here.
- Parameters:
status (str) – Possible values are ‘initial’ or ‘update’. The former should be used when determining the initial definition.
dependences (List[str]) – Possible values are ‘const’, ‘Umean’, ‘E/N’, ‘ESR’, or ‘Tgas’. Can also be set to zero to skip this k_coeff.
k_coeff (List[df.Function]) – List of coefficient variables
kx (List[Any]) – Look-up table data
ky (List[Any]) – Look-up table data
energy (df.Function) – Energy function
redfield (df.Function) – Reduced electric field.
Tgas (float, default 300.0) – Gas temperature. Not used directly, but may be used by user-supplied function.
Te (float, default 0.0) – Electron temperature. Not used directly, but may be used by user-supplied function.
- Return type:
None
- Raises:
ValueError – If given incorrect status or dependences, or if the provided lists are not all of the same length.
Exception – The user may supply their own functions using the ‘fun:Te,Tgas’ or ‘fun:Tgas’ dependences. These are evaluated using ‘eval’, meaning anything could happen.
- fedm.functions.Source_term(coupling: str, approx: str, p_matrix: ndarray, l_matrix: ndarray, g_matrix: ndarray, k_coeffs: List[Function], N0: float, u: Any) List[Any] #
Defines source term for coupled or uncoupled approach, with LFA (counting particles from 0) or LMEA approximation (counting particle from 1 and the zeroth equation is energy).
- Parameters:
coupling (str) – Either ‘coupled’ or ‘uncoupled’
approx (str) – Either ‘LFA’ or ‘LMEA’
p_matrix (np.ndarray) – Power matrix
l_matrix (np.ndarray) – Loss matrix
g_matrix (np.ndarray) – Gain matrix
k_coeffs (List[df.Function]) – Rate coefficient, energy loss and mean energy
N0 (float) – Gas number density
u – Trial function
- Returns:
Source terms
- Return type:
List[Any]
- Raises:
ValueError – If coupling or approx are invalid. Could also be raised if the inputs have the wrong shapes.
- fedm.functions.Transport_coefficient_interpolation(status: str, dependences: List[str], N0: float, Tgas: float, k_coeffs: List[Function], kxs: List[Any], kys: List[Any], energy: Function, redfield: Function, mus: List[Function] | None = None) None #
Function for linear interpolation of transport coefficients. Modifies k_coeffs.
- Parameters:
status (str) – Possible values are ‘initial’ or ‘update’. The former should be used when determining the initial definition.
dependences (List[str]) – Possible values are ‘const’, ‘Umean’, ‘E/N’, ‘ESR’, or ‘Tgas’. Can also be set to 0 to skip this k_coeff.
N0 (float) – Gas number density
Tgas (float) – Gas temperature
k_coeff (List[df.Function]) – List of coefficient variables
kx (List[Any]) – Look-up table data
ky (List[Any]) – Look-up table data
energy (df.Function) – Energy function
redfield (df.Function) – Reduced electric field.
mus (List[df.Function], default None) – Mobilities
- Return type:
None
- Raises:
ValueError – If given incorrect status or dependences, or if using ESR dependence without providing mus, or if the provided lists are not all of the same length.
- fedm.functions.adaptive_solver(nonlinear_solver: PETScSNESSolver | NonlinearVariationalSolver, problem: Problem, t: float, dt: Expression, dt_old: Expression, u_new: Function, u_old: Function, var_list_new: List[Any], var_list_old: List[Any], assigner: FunctionAssigner, error: List[float], error_file: Path, max_error: List[float], ttol: float, dt_min: float, time_dependent_arguments: List[Any] | None = None, approximation: str = 'LMEA') float #
This function is used for solving the problem when adaptive time stepping is used.
- Parameters:
nonlinear_solver (Union[df.PETScSNESSolver, df.NonlinearVariationalSolver]) – The Dolfin solver to use.
problem (Problem) – Nonlinear problem definition
t (float) – Time
dt (df.Expression) – Time step size
dt_old (df.Expression) – Previous time step size
u_new (df.Function) – New variables defined on mixed function space
u_old (df.Function) – Old variables defined on mixed function space
var_list_new (List[Any]) – List of new variables for postprocessing
var_list_old (List[Any]) – List of old variables for postprocessing
assigner (df.FunctionAssigner) – Assigns values between variables
error (List[float]) – Record of errors
error_file (Path) – Error file name
max_error (List[float]) – Maximum errors
ttol (float) – Timestepping tolerance
dt_min (float) – Minimum timestep allowed
time_dependent_arguments (List[Any], default None) – List of functions that need to be updated with time
approximation (str, default 'LMEA') – Type of approximation to use, options are ‘LMEA’ or ‘LFA’.
- Returns:
The new time ‘t’
- Return type:
float
- Raises:
SystemExit – If the time step size is reduced below dt_min
- fedm.functions.adaptive_timestep(dt, error, tol=0.0001, dt_min=1e-13, dt_max=1e-09)#
Function calculates new time step based on a PID controller M. Moeller, Time stepping methods, ATHENS course: Introductioninto Finite Elements, Delft Institute of Applied Mathematics, TU Delft(2015). Input arguments are time step, error, tolerance for time step controller, minimal and maximal time step.
- fedm.functions.adaptive_timestep_H211b(dt, dt_old, error, tol=0.0001, dt_min=1e-13, dt_max=1e-09)#
Function calculates new time step size using H211b controller (G. Soederlind, Acm. T. Math. Software 29: 1-26, 2003). Input arguments are time step size, previous time step size, error, tolerance for time step controller, minimal and maximal time step.
- fedm.functions.adaptive_timestep_PI34(dt, error, tol=0.0001, dt_min=1e-13, dt_max=1e-09)#
Function calculates new time step size using PI.3.4 controller (G. Soederlind Numerical Algorithms 31: 281-310, 2002). Input arguments are time-step size, error, tolerance for time step controller, minimal and maximal time-step size.
- fedm.functions.modify_approximation_vars(approximation_type: str, number_of_species: int, particle_species: List[str], masses: List[float], charges: List[float]) Tuple[int, int, List[str], List[float], List[float]] #
Depending on approximation used, the number of equations, charge and mass variables are modified. Returns number species, number of equations, particle species, masses, and charges. particle_species, masses, and charges may be modified.
- fedm.functions.semi_implicit_coefficients(dependences: List[str], mean_energy_new: Any, mean_energy_old: Function, coefficients: List[Function], coefficient_diffs: List[Function]) List[Function] #
- fedm.functions.weak_form_Poisson_equation(dx, u, v, f, r=0.15915494309189535) Form #
Returns a weak form of Poisson equation. Input arguments are dV, trial function, test function, source term and r coordinate.
- Parameters:
dx (df.Measure) – dV, used to build integrals. Recommended to use dolfin.dx.
u – Trial function
v – Test function
f – Source term
r – r coordinate
- Return type:
df.Form
- fedm.functions.weak_form_balance_equation(equation_type: str, dt: Expression, dt_old: Expression, dx: Measure, u: Any, u_old: Any, u_old1: Any, v: Any, f: Function, Gamma: Function, r: float = 0.15915494309189535, D: Function | None = None, log_representation: bool = False) Form #
Returns the weak form of the particle balance equations.
If log_representation is True, solves:
\[2 \pi \int_{\Omega}{} e^{u_\mathrm{k+1}} \frac{1+2 \omega_k}{1 + \omega_k} \Bigl( u_{k+1} - u_{k} \frac{1 + \omega_k^2}{1 + 2 \omega_k} + u_{k-2} \frac{\omega_k^2}{1 + 2 \omega_k}) \Bigr) \frac{v}{\Delta t_k} r \mathrm{d}x\]otherwise, solves:
\[2 \pi \int_{\Omega}{} \frac{1+2 \omega_k}{1 + \omega_k} \Bigl( u_{k+1} - u_{k} \frac{1 + \omega_k^2}{1 + 2 \omega_k} + u_{k-2} \frac{\omega_k^2}{1 + 2 \omega_k}) \Bigr) \frac{v}{\Delta t_k} r \mathrm{d}x\]where \(\omega_\mathrm{k} = \frac{\Delta t_\mathrm{k}}{\Delta t_\mathrm{k-1}}\).
If solving diffusion-reaction equation, also includes the term:
for log representation
\[-2 \pi \int_{\Omega}{} (- \nabla(D_k e^u_{k+1}) \cdot \nabla v) r \mathrm{d}x\]otherwise
\[-2 \pi \int_{\Omega}{} (- \nabla(D_k u_{k+1}) \cdot \nabla v) r \mathrm{d}x\]If solving drift-diffusion-reaction, instead includes the term:
\[-2 \pi \int_{\Omega}{} \Gamma_{k+1} \cdot \nabla v r \mathrm{d}x\]In all cases, also includes the source term:
\[- 2 \pi \int_{\Omega}{} f v r \mathrm{d}x\]- Parameters:
equation_type (str) – Type of equation to solve. Options are ‘reaction’, ‘diffusion-reaction’, or ‘drift-diffusion-reaction’.
dt (df.Expression) – Current time-step size
dt_old (df.Expression) – Previous time-step size
dx (df.Measure) – dV, used to build integrals. Recommended to use dolfin.dx.
u – Trial function
u_old – Value of variable in current time step
u_old1 – Value of variable in previous time step
v – Test function
f (df.Function) – Source term
Gamma (df.Function) – particle flux
r (float, default 0.5 / pi) – r coordinate
D (df.Function, default None) – diffusion coefficient, only required for the diffusion equation.
log_representation (bool, default False) – Use logarithmic representation.
- Return type:
df.Form
- Raises:
ValueError – If equation_type is not recognised, or if D is not supplied when solving the diffusion-reaction equation.
- fedm.functions.weak_form_balance_equation_log_representation(*args, **kwargs) Form #
Convenience function, calls weak_form_balance_equation with log_representation set to True.
fedm.physical_constants module#
This file collects physical constants used by other FEDM modules.
fedm.utils module#
- fedm.utils.comma_separated(strings: List[str]) str #
Utility function, takes a list of strings and returns a single string containing each of those strings in single quotes and separated by ‘, ‘
- fedm.utils.mesh_info(mesh: Mesh) str #
- fedm.utils.print_rank_0(*args, **kwargs) None #
Utility function, print to terminal if MPI rank 0, otherwise do nothing.