3. Bolos API reference

This page contains the documentation of the BOLOS API.

3.1. The solver Module

This module contains the main routines to load processes, specify the physical conditions and solve the Boltzmann equation.

The data and calculations are encapsulated into the BoltzmannSolver class, which you have to instantiate with a grid.Grid instance. Use BoltzmannSolver.load_collisions() or BoltzmannSolver.add_process() to add processes with their cross-sections. Afterwards, set the density of each component with BoltzmannSolver.set_density() or BoltzmannSolver.target. The method BoltzmannSolver.maxwell() gives you a reasonable initial guess for the electron energy distribution function (EEDF) that you can then improve iteratively with BoltzmannSolver.converge(). Finally, methods such as BoltzmannSolver.rate() or BoltzmannSolver.mobility() allow you to obtain reaction rates and transport parameters for a given EEDF.

class bolos.solver.BoltzmannSolver(grid)[source]

Bases: object

Class to solve the Boltzmann equation for electrons in a gas.

This class contains the required elements to specify the conditions for the solver and obtain the equilibrium electron energy distribution function.

Parameters:

grid : grid.Grid

The grid in energies where the distribution funcition will be evaluated.

Examples

>>> import numpy as np
>>> from bolos import solver, grid
>>> grid.LinearGrid(0, 60., 400)
>>> bsolver = solver.BoltzmannSolver(grid)
>>> # Parse the cross-section file in BOSIG+ format and load it into the
>>> # solver.
>>> with open(args.input) as fp:
>>>     processes = parser.parse(fp)
>>> bsolver.load_collisions(processes)
>>> 
>>> # Set the conditions.  And initialize the solver
>>> bsolver.target['N2'].density = 0.8
>>> bsolver.target['O2'].density = 0.2
>>> bsolver.kT = 300 * co.k / co.eV
>>> bsolver.EN = 300.0 * solver.TOWNSEND
>>> bsolver.init()
>>> 
>>> # Start with Maxwell EEDF as initial guess.  Here we are starting with
>>> # with an electron temperature of 2 eV
>>> f0 = bsolver.maxwell(2.0)
>>> 
>>> # Solve the Boltzmann equation with a tolerance rtol and maxn 
>>> # iterations.
>>> f1 = bsolver.converge(f0, maxn=50, rtol=1e-5)

Attributes

benergy (array of floats) Cell boundaries of the energy grid (set automatically at initialization). Equivalent to grid.b.
benergy (array of floats) Cell lengths of the energy grid (set automatically at initialization). Equivalent to grid.d.
cenergy (array of floats) Cell centers of the energy grid (set automatically at initialization). Equivalent to grid.c.
n (int) Number of cells in the energy grid (set automatically at initialization). Equivalent to grid.n.
kT (float) Gas temperature in eV. Must be set by the user.
EN (float) Reduced electric field in Townsend (1 Td is 1e-21 V m^2). Must be set by the user.
target (dict) A dictionary with targets in the set of processes. The user needs to set the density (molar fraction) of the desired targets using this dictionary. E.g. synthetic air is represented by

Methods

add_process(**kwargs) Adds a new process to the solver.
converge(f0[, maxn, rtol, delta0, m, full]) Iterates and attempted EEDF until convergence is reached.
diffusion(F0) Calculates the diffusion coefficient from a distribution function.
init() Initializes the solver with given conditions and densities of the target species.
iter_all() Iterates over all processes.
iter_elastic() Iterates over all elastic processes.
iter_growth() Iterates over all processes that affect the growth of electron density, i.e.
iter_inelastic() Iterates over all inelastic processes.
iter_momentum()
iterate(f0[, delta]) Iterates once the EEDF.
load_collisions(dict_processes) Loads the set of collisions from the list of processes.
maxwell(kT) Calculates a Maxwell-Boltzmann distribution function.
mean_energy(F0) Calculates the mean energy from a distribution function.
mobility(F0) Calculates the reduced mobility (mobility * N) from the EEDF.
rate(F0, k[, weighted]) Calculates the rate of a process from a (usually converged) EEDF.
search(signature[, product, first]) Search for a process or a number of processes within the solver.
set_density(species, density) Sets the molar fraction of a species.
add_process(**kwargs)[source]

Adds a new process to the solver.

Adds a new process to the solver. The process data is passed with keyword arguments.

Parameters:

type : string

one of “EFFECTIVE”, “MOMENTUM”, “EXCITATION”, “IONIZATION” or “ATTACHMENT”.

target : string

the target species of the process (e.g. “O”, “O2”...).

ratio : float

the ratio of the electron mass to the mass of the target (for elastic/momentum reactions only).

threshold : float

the energy threshold of the process in eV (only for inelastic reactions).

data : array or array-like

cross-section of the process array with two columns: column 0 must contain energies in eV, column 1 contains the cross-section in square meters for each of these energies.

Returns:

process : process.Process

The process that has been added.

See also

load_collisions
Add a set of collisions.

Examples

>>> import numpy as np
>>> from bolos import solver, grid
>>> grid.LinearGrid(0, 60., 400)
>>> solver = BoltzmannSolver(grid)
>>> # This is an example cross-section that decays exponentially
>>> energy = np.linspace(0, 10)
>>> cross_section = 1e-20 * np.exp(-energy)
>>> solver.add_process(type="EXCITATION", target="Kriptonite", 
>>>                    ratio=1e-5, threshold=10, 
>>>                    data=np.c_[energy, cross_section])
converge(f0, maxn=100, rtol=1e-05, delta0=100000000000000.0, m=4.0, full=False, **kwargs)[source]

Iterates and attempted EEDF until convergence is reached.

Parameters:

f0 : array of floats

Initial EEDF.

maxn : int

Maximum number of iteration until the convergence is declared as failed (default: 100).

rtol : float

Target tolerance for the convergence. The iteration is stopped when the difference between EEDFs is smaller than rtol in L1 norm (default: 1e-5).

delta0 : float

Initial value of the iteration parameter. This parameter is adapted in succesive iterations to improve convergence. (default: 1e14)

m : float

Attempted reduction in the error for each iteration. The Richardson extrapolation attempts to reduce the error by a factor m in each iteration. Larger m means faster convergence but also possible instabilities and non-decreasing errors. (default: 4)

full : boolean

If true returns convergence information besides the EEDF.

Returns:

f1 : array of floats

Final EEDF

iters : int (returned only if full is True)

Number of iterations required to reach convergence.

err : float (returned only if full is True)

Final error estimation of the EEDF (must me smaller than rtol).

Notes

If convergence is not achieved after maxn iterations, an exception of type ConvergenceError is raised.

diffusion(F0)[source]

Calculates the diffusion coefficient from a distribution function.

Parameters:

F0 : array of floats

The EEDF used to compute the diffusion coefficient.

Returns:

diffn : float

The reduced diffusion coefficient of electrons in SI units..

See also

mobility
Find the reduced mobility from the EEDF.
grid
init()[source]

Initializes the solver with given conditions and densities of the target species.

This method does all the work previous to the actual iterations. It has to be called whenever the densities, the gas temperature or the electric field are changed.

Notes

The most expensive calculations in this method are cached so they are not repeated in each call. Therefore the execution time may vary wildly in different calls. It takes very long whenever you change the solver’s grid; therefore is is strongly recommended not to change the grid if is not strictly neccesary.

iter_all()[source]

Iterates over all processes.

Returns:An iterator over (target, process) tuples.
iter_elastic()[source]

Iterates over all elastic processes.

Returns:An iterator over (target, process) tuples.
iter_growth()[source]

Iterates over all processes that affect the growth of electron density, i.e. ionization and attachment.

Returns:An iterator over (target, process) tuples.
iter_inelastic()[source]

Iterates over all inelastic processes.

Returns:An iterator over (target, process) tuples.
iter_momentum()[source]
iterate(f0, delta=100000000000000.0)[source]

Iterates once the EEDF.

Parameters:

f0 : array of floats

The previous EEDF

delta : float

The convergence parameter. Generally a larger delta leads to faster convergence but a too large value may lead to instabilities or slower convergence.

Returns:

f1 : array of floats

A new value of the distribution function.

Notes

This is a low-level routine not intended for normal uses. The standard entry point for the iterative solution of the EEDF is the BoltzmannSolver.converge() method.

load_collisions(dict_processes)[source]

Loads the set of collisions from the list of processes.

Loads a list of dictionaries containing processes.

Parameters:

dict_processes : List of dictionary or dictionary-like elements.

The processes to add to this solver class. See :method:`solver.add_process` for the required fields of each of the dictionaries.

Returns:

processes : list

A list of all added processes, as process.Process instances.

See also

add_process
Add a single process, with its cross-sections, to this solver.
maxwell(kT)[source]

Calculates a Maxwell-Boltzmann distribution function.

Parameters:

kT : float

The electron temperature in eV.

Returns:

f : array of floats

A normalized Boltzmann-Maxwell EEDF with the given temperature.

Notes

This is often useful to give a starting value for the EEDF.

mean_energy(F0)[source]

Calculates the mean energy from a distribution function.

Parameters:

F0 : array of floats

The EEDF used to compute the diffusion coefficient.

Returns:

energy : float

The mean energy of electrons in the EEDF.

mobility(F0)[source]

Calculates the reduced mobility (mobility * N) from the EEDF.

Parameters:

F0 : array of floats

The EEDF used to compute the mobility.

Returns:

mun : float

The reduced mobility (mu * n) of the electrons in SI units (V / m / s).

See also

diffusion
Find the reduced diffusion rate from the EEDF.

Examples

>>> mun = bsolver.mobility(F0)
rate(F0, k, weighted=False)[source]

Calculates the rate of a process from a (usually converged) EEDF.

Parameters:

F0 : array of floats

Distribution function.

k : process.Process or string

The process whose rate we want to calculate. If k is a string, it is passed to search() to obtain a process instance.

weighted : boolean, optional

If true, the rate is multiplied by the density of the target.

Returns:

rate : float

The rate of the given process according to F0.

See also

search
Find a process that matches a given signature.

Examples

>>> k_ionization = bsolver.rate(F0, "N2 -> N2^+")
search(signature, product=None, first=True)[source]

Search for a process or a number of processes within the solver.

Parameters:

signature : string

Signature of the process to search for. It must be in the form “TARGET -> RESULT [+ RESULT2]...”.

product : string

If present, the first parameter is interpreted as TARGET and the second parameter is the PRODUCT.

first : boolean

If true returns only the first process matching the search; if false returns a list of them, even if there is only one result.

Returns:

processes : list or process.Process instance.

If first was true, returns the first process matching the search. Otherwise returns a (possibly empty) list of matches.

Examples

>>> ionization = solver.search("N2 -> N2^+")[0]
>>> ionization = solver.search("N2", "N2^+", first=True)
set_density(species, density)[source]

Sets the molar fraction of a species.

Parameters:

species : str

The species whose density you want to set.

density : float

New value of the density.

Examples

These are two equivalent ways to set densities for synthetic air:

Using set_density():

bsolver.set_density('N2', 0.8)
bsolver.set_density('O2', 0.2)

Using bsolver.target:

bsolver.target['N2'].density = 0.8
bsolver.target['O2'].density = 0.2
exception bolos.solver.ConvergenceError[source]

Bases: exceptions.Exception

3.2. The grid Module

Routines to handle different kinds of grids (linear, quadratic, logarithmic)

class bolos.grid.AutomaticGrid(grid, f0, delta=0.0001)[source]

Bases: bolos.grid.Grid

A grid set automatically using a previous estimation of the EEDF to fix a peak energy.

Methods

cell(x) Returns the cell index containing the value x.
interpolate(f, other) Interpolates into this grid an eedf defined in another grid.
class bolos.grid.GeometricGrid(x0, x1, n, r=1.1)[source]

Bases: bolos.grid.Grid

A grid with geometrically progressing spacing. To be more precise, here the length of cell i+1 is r times the length of cell i.

Methods

cell(x) Returns the cell index containing the value x.
f(x)
finv(w)
interpolate(f, other) Interpolates into this grid an eedf defined in another grid.
f(x)[source]
finv(w)[source]
class bolos.grid.Grid(x0, x1, n)[source]

Bases: object

Class to define energy grids.

This class encapsulates the information about an energy grid.

Parameters:

x0 : float

Lowest boundary energy.

x1 : float

Highest energy boundary.

n : float

Number of cells

See also

LinearGrid
A grid with linear spacings (constant cell length).
QuadraticGrid
A grid with quadratic spacings (linearly increasing cell length).
GeometricGrid
A grid with geometrically increasing cell lengths.
LogGrid
A logarithmic grid.

Notes

This is a base class and you usually do not want to instantiate it directly. You can define new grid classes by subclassing this class and then defining an f method that maps energy to a new variable y that is divided uniformly.

Methods

cell(x) Returns the cell index containing the value x.
interpolate(f, other) Interpolates into this grid an eedf defined in another grid.
cell(x)[source]

Returns the cell index containing the value x.

Parameters:

x : float

The value x which you want to localize.

Returns:

index : int

The index to the cell containing x

interpolate(f, other)[source]

Interpolates into this grid an eedf defined in another grid.

Parameters:

f : array or array-like

The original EEDF

other : Grid

The old grid, where f is defined.

Returns:

fnew : array or array-like

An EEDF defined in our grid.

class bolos.grid.LinearGrid(x0, x1, n)[source]

Bases: bolos.grid.Grid

A grid with linear spacing.

Methods

cell(x) Returns the cell index containing the value x.
f(x)
finv(w)
interpolate(f, other) Interpolates into this grid an eedf defined in another grid.
f(x)[source]
finv(w)[source]
class bolos.grid.LogGrid(x0, x1, n, s=10.0)[source]

Bases: bolos.grid.Grid

A pseudo-logarithmic grid. We add a certain s to the variable to avoid log(0) = -inf. The grid is actually logarithmic only for x >> s.

Methods

cell(x) Returns the cell index containing the value x.
f(x)
finv(w)
interpolate(f, other) Interpolates into this grid an eedf defined in another grid.
f(x)[source]
finv(w)[source]
class bolos.grid.QuadraticGrid(x0, x1, n)[source]

Bases: bolos.grid.Grid

A grid with quadratic spacing.

Methods

cell(x) Returns the cell index containing the value x.
f(x)
finv(w)
interpolate(f, other) Interpolates into this grid an eedf defined in another grid.
f(x)[source]
finv(w)[source]
bolos.grid.mkgrid(kind, *args, **kwargs)[source]

Builds and returns a grid of class kind. Possible values are ‘linear’, ‘lin’, ‘quadratic’, ‘quad’, ‘logarithmic’, ‘log’.

3.3. The parser Module

This module contains the code required to parse BOLSIG+-compatible files. To make the code re-usabe in other projects it is independent from the rest of the BOLOS code.

Most user would only use the method parse() in this module, which is documented below.

bolos.parser.parse(fp)[source]

Parses a BOLSIG+ cross-sections file.

Parameters:

fp : file-like

A file object pointing to a Bolsig+-compatible cross-sections file.

Returns:

processes : list of dictionaries

A list with all processes, in dictionary form, included in the file.