Organic Electronic Device Simulator (oedes) documentation¶
Introduction¶
oedes
(Open Electronic DEvice Simulator) is a Python software package for modeling of electronic devices. It is primarily focused on emerging electronic devices. It is applicable to organic electronic, electrochemical, bioelectronic, and Peroskvite based devices.
oedes was written to take into account both the special aspects of non-conventional electronic devices, and modern trends in software development. The result a small but powerful package, which is tightly integrated with standard scientific software stack. This simplifies both running and sharing the simulations.
oedes mission is to:
- enable open science research in emerging electronic devices, by enabling full disclosure of numerical simulations;
- allow reuse of established device models, which were usually published without a runnanble implementation;
- accelerate research by providing a platform for parameter extraction from measurements and for device optimization driven by simulation;
- simplify development of new simulation models and tools, by providing a library of robust components;
- provide a standard way of running and sharing device simulations.
oedes is released as open-source under GNU Affero General Public License v3. It is free to use and distributed with complete source code.
Conventions¶
Code which is intended to execute in system command-line starts with $
, for example:
$ python
All code which is indendent to be run in Python starts with >>>
or ...
, for example:
>>> print('hello world')
Sequences of characters which can be copied into code in current context are emphasized as print
or 1+2+3
. oedes
formatted in that way refers to name of Python software package, or to Python symbol in current example.
General concepts are emphasized as params. oedes written in that way refers to the software project.
Features¶
General
- Pure Python implementation
- Reliable fully implicit solver
- Flexible precision
Electric
- Cell centered finite-volume with arbitrary mesh spacing
- Scharfetter-Gummel discretization
- Conservation of charge guaranteed in a converged solution
- Arbitrary number of layers
- Ramo-Shockley calculation of terminal current
- Arbitrary temperature
- Support for multidimensional discretizations
- Arbitrary number of transported species
Optical
- Transfer matrix approach
- Arbitrary illumination
Analysis
- DC
- Transient
- AC small signal
- Sensitivity analysis
Models
- Gaussian DOS
- Extended Gaussian Disorder Model
- Generalized Einstein’s relation
- Transient traps
- Shockley-Read-Hall recombination
- Langevin recombination
- Image-force barrier lowering
- Doping (also position dependent)
- Generation (also position dependent)
- Onsager-Braun model of exciton dissociation
- User supplied mobility model
- User supplied DOS
Outputs
- Current-voltage
- Current-voltage-light
- SCL
- Transient SCL
- Capacitance-voltage
- Impedance spectroscopy
- Efficiency (photovoltaic)
- Fill factor (photovoltaic)
- Dark current (photodetector)
- Responsivity (photodetector)
- Light output (light emitting diodes)
Installation¶
Requirements¶
- Python 2.7 or Python 3.4+
- sparsegrad
- Python scientific stack
Installation from PyPI¶
It is recommended to use Python Package Index (PyPI) to install oedes
package. This is done using command-line by pip
program, which is normally installed together with Python. Using this method, all dependencies are resolved automatically.
Two variants of the installation are possible:
- system wide installation:
$ pip install oedes
- local installation not requiring administrator’s rights:
$ pip install oedes --user
In the case of local installation, oedes is installed inside user’s home directory. In Linux, this defaults to $HOME/.local
.
Verifying the installation¶
After installing, it is advised to run the test suite to ensure that oedes works correctly on your system:
>>> import oedes
>>> oedes.test()
Running unit tests for oedes...
OK
<nose.result.TextTestResult run=15 errors=0 failures=0>
If any errors are found, oedes is not compatible with your system. Either your Python scientific stack is too old, or there is a bug.
oedes is evolving, and backward compatibility is not yet offered. It is recommended to check which version is in use by running:
>>> import oedes
>>> oedes.version
'0.0.18'
Upgrading/downgrading¶
By default, pip
does not upgrade packages unless required. To change this behavior, option --upgrade
should be used. For example, oedes is updated to the most recent version by running
$ pip install oedes --upgrade
It is also possible to upgrade/downgrade oedes by specifying an exact version to be installed. For example, to ensure that installed version is 0.0.18
, run
$ pip install oedes==0.0.18
Again, --user
option can be given to restrict changes to users’ home directory.
Development installation (advanced)¶
Current development version of sparsegrad can be installed from the development repository by running
$ git clone https://github.com/mzszym/oedes.git
$ cd oedes
$ pip install -e .
The option -e
tells that oedes
code should be loaded from git
controlled directory, instead of being copied to the Python libraries directory. As with the regular installation, --user
option should be appended for local installation.
Tutorial: PN junction¶
In this tutorial a simple model of PN junction will be constructed. The same model is used in example included in oedes distribution in file examples/interactive/pn.ipynb
.
oedes
contains models of typical devices. However, in this tutorial, the model will be constructed equation-by-equation to demonstrate basic functionality.
The drift-diffusion system¶
For simulating transport, a model of temperature distribution must be assumed. In the simplest case of isothermal simulation, the same effective temperature is assumed equal everywhere inside the device. The model is created as follows:
>>> temperature = oedes.models.ConstTemperature()
The specification of model is separated from the specification of parameter values. The actual value of temperature, in Kelvins, will be given later.
The next step is to create Poisson’s equation of electrostatics
>>> poisson = oedes.models.PoissonEquation()
Constructed objects poisson
and temperature
must are passed as arguments when creating transport equations:
>>> electron = oedes.models.equations.BandTransport(poisson=poisson, name='electron', z=-1, thermal=temperature)
>>> hole = oedes.models.equations.BandTransport(poisson=poisson, name='hole', z=1, thermal=temperature)
Above creates two conventional drift-diffusion equations for electrons and holes respectively. By default, the mobility is assumed constant and the DOS is modeled by using the Boltzmann approximation. name
arguments are names which are used to identify parameters and outputs. z
are charges of species expressed in units of elementary charge. oedes allows arbitrary number of species, and arbitrary values of name
and z
. This allows to construct complicated models with, for example, mixed ionic-electronic transport.
The doping profile¶
To model the PN junction, a doping profile must be defined. In the example, left half of the device is doped with ionized donor concentration given by parameter Nd
, and right half of device is doped with doped with ionized acceptor concentration given by parameter Na
.
The pn_doping
function is called during model evaluation and is given as parameters the mesh, the evaluation context object ctx
, and the discretized equation object eq
. In the example, it uses ctx
object to access parameters values of the dopant concentations ('Na'
, ''Nd'
).
>>> def pn_doping(mesh, ctx, eq):
... return oedes.ad.where(mesh.x<mesh.length*0.5,ctx.param(eq, 'Nd'),- ctx.param(eq, 'Na'))
>>> doping = oedes.models.FixedCharge(poisson, density=pn_doping)
The dopants are assumed to be fully ionized and therefore it is modeled as fixed charge. The FixedCharge
adds calculated doping profile to the previously created Poisson’s equation poisson
.
Above code uses a specialized version of function where
is used instead of version from numpy
. This is required for support of sensitivity analysis with respect to parameter values.
The Ohmic contacts¶
To keep the example simple, Ohmic contacts are assumed on both sides of the device. They are created as follows:
>>> semiconductor = oedes.models.Electroneutrality([electron, hole, doping],name='semiconductor')
>>> anode = oedes.models.OhmicContact(poisson, semiconductor, 'electrode0')
>>> cathode = oedes.models.OhmicContact(poisson, semiconductor, 'electrode1')
Ohmic contacts require knowledge of equilibrium charge carrier concentrations in semiconductor. This is calculated by Electroneutrality
. Note that since concentrations in doped semiconductor are of interest, all charged species are passed to Electroneutrality. 'electrode0'
and 'electrode1'
refers to names of boundaries in the mesh.
Putting all together¶
To avoid divergence of the simulation due to infinitely large lifetime of electrons and holes, recombination should be added. Duirecrecombination model is created by
>>> recombination = oedes.models.DirectRecombination(semiconductor)
The calculation of terminal current is a non-trivial post-processing step. It is recommended to use Ramo-Shockley current calculation in most cases, which is created by
>>> current = oedes.models.RamoShockleyCurrentCalculation([poisson])
The discrete model is constructed and initialized by calling oedes.fvm.discretize
. It takes two arguments: the system of equations and terms to solve, and the specification of domain. Below oedes.fvm.mesh1d
creates a 1-D domain with length specified as argument.
>>> all_equations_and_terms = [ poisson, temperature, electron, hole, doping, current, semiconductor, anode, cathode, recombination ]
>>> domain = oedes.fvm.mesh1d(100e-9)
>>> model = oedes.fvm.discretize(all_equations_and_terms, domain)
Parameters¶
The physical parameters are provided as dict
.
>>> params={
... 'T':300,
... 'epsilon_r':12,
... 'Na':1e24,
... 'Nd':1e24,
... 'hole.mu':1,
... 'electron.mu':1,
... 'hole.energy':-1.1,
... 'electron.energy':0,
... 'electrode0.voltage':0,
... 'electrode1.voltage':0,
... 'hole.N0':1e27,
... 'electron.N0':1e27,
... 'beta':1e-9
... }
Above, 'T'
key is used to specify temperature in Kelvins. It is used by ConstTemperature
object. 'epsilon_r'
specifies the relative dielectric permittivity. It is used by discretized PoissonEquation
object. 'Na'
and 'Nd'
are parameters accessed by pn_doping
function, the concentrations of dopants. 'beta'
is used by DirectRecombination
.
Other parameters are in form name.parameter
. name
is passed to the equation, and they can be nested. For example, if a transport equation were created as
something = oedes.models.BandTransport(name='zzz',...)
then the corresponding mobility parameter would be identified by key 'zzz.mu'
.
The mobilities electron.mu
and hole.mu
are given in , therefore are equal to 1000 each. In the example above, instead of specifying electron affinity and band-gap, the energies of both bands are specified directly by energy parameters, in eV. The voltages are applied to Ohmic contacts are specified by 'electrode0.voltage'
and 'electrode1.voltage'
, in Volts. N0
denotes the total density of states, in .
params¶
By convention, values of physical parameters are specified in dict
object named params
, with string keys, and float values. All values of given in SI base units, except for small energies which are specified in eV.
oedes currently does not assume default values of parameters. If any necessary parameter is not specified in params
, exception KeyError
is raised.
Solving¶
oedes.context
objects binds models with their parameters and solutions. It also provides convenience functions for solving, post-processing and plotting the data.
The following calculates soltuion for parameters specified in dict params
>>> c = oedes.context(model)
>>> c.solve(params)
Examining output¶
The solution can be investigated by calling output function, which returns a dict of available outputs:
>>> out=c.output()
>>> print(sorted(out.keys()))
['.meta', 'D', 'Dt', 'E', 'Et', 'J', 'R', 'c', 'charge', 'electrode0.J', 'electrode1.J', 'electron.Eband', 'electron.Ef', 'electron.J', 'electron.Jdiff', 'electron.Jdrift', 'electron.c', 'electron.charge', 'electron.ct', 'electron.j', 'electron.jdiff', 'electron.jdrift', 'electron.phi_band', 'electron.phi_f', 'hole.Eband', 'hole.Ef', 'hole.J', 'hole.Jdiff', 'hole.Jdrift', 'hole.c', 'hole.charge', 'hole.ct', 'hole.j', 'hole.jdiff', 'hole.jdrift', 'hole.phi_band', 'hole.phi_f', 'potential', 'semiconductor.Ef', 'semiconductor.electron.c', 'semiconductor.hole.c', 'semiconductor.phi', 'total_charge_density']
The outputs are numpy arrays. For example, the electrostatic potential is
>>> print(out['potential'])
[-0.17857923 -0.17858006 -0.17858115 -0.17858258 -0.17858445 -0.17858693
-0.17859023 -0.17859468 -0.1785994 -0.17860442 -0.17860982 -0.17861567
...
-0.92141696 -0.92141772]
To access additional information about output (such as its mesh), use .meta
subdictionary.
>>> out['.meta']['potential']
OutputMeta(mesh=<oedes.fvm.mesh.mesh1d object at ...>, face=False, unit='V')
outputs¶
Most useful outputs are given below. Just as for params, all values are in SI base units, except for small energies in eV. * denotes prefix identifying the equation, such as electron
or hole
.
*.c
: concentration of particles, in*.j
: flux of particles, in*.Ef
: quasi Fermi level, in eV*.Eband
: band energy, in eVR
: recombination density, inJ
: total electric current density, inE
: electric field, in V/mpotential
: electrostatic potential, V
Plotting¶
oedes.context
object simplifies plotting results using matplotlib
. For example, bands and quasi Fermi levels are plotted as
>>> import matplotlib.pylab as plt
>>> fig,ax = plt.subplots()
>>> p=c.mpl(fig, ax)
>>> p.plot(['electron.Eband'],label='$E_c$')
>>> p.plot(['hole.Eband'],label='$E_v$')
>>> p.plot(['electron.Ef'],linestyle='--',label='$E_{Fn}$')
>>> p.plot(['hole.Ef'],linestyle='-.',label='$E_{Fp}$')
>>> p.apply_settings({'xunit':'n','xlabel':'nm'})
Physical models¶
Concentrations of charges in thermal equilibrium¶
Probability that an electronic state with energy is occupied is given by the Fermi-Dirac distribution:
with denoting the Fermi energy, denoting the Boltzmann constant and denoting the temperature.
The states in the conduction band are distributed in energy according to the density of states function . The total concentration of electrons is given by integral
In the valence band, almost all states are occupied by the electrons. It is therefore useful to track unoccupied states, holes, instead of the occupied states. The concentration of holes is given by:
Noting that
and changing the integration variable, both concentrations can be written in a common form as:
(1)¶
Practical note: The SI base unit of energy is Joule (J), however the energies such as are very small and should be expressed in electronvolts (eV). The SI basic unit of concentration is , although is often encountered. The value of at the room temperature is approximately 26 meV. The SI unit of temperature is Kelvin, .
Band energies¶
In an idealized case, the energies and of the conduction and valence bands are
(2)¶
where is the electron affinity energy, and is the bandgap energy. is the electrostatic potential. q is the elementary charge.
Electrostatic potential¶
The electric field is related to the elestrostatic potential as
(3)¶
In linear, isotropic, homogeneous medium the electric displacement field is
with permittivity
where is the vacuum permittivity, and is the relative permittivity of the material.
The electric displacement field satisfies the electric the Gauss’s equation
(4)¶
where is the density of free charge
with denoting the elementary charge. Above, denotes other charges, such as ionized dopants.
Combining the above equations gives the usual Poisson’s equation for electrostatics:
(5)¶
The SI unit of electrostatic potential is Volt, and the unit of electric field is Volt/meter. The unit of permittivity is Farad/meter. The unit of charge density is .
Approximation for low concentrations¶
If the concentration of charge carriers is low enough, only states on the edge of band gap are important. In such case, the density of states can be assumed as a sharp energetic level,
in case of electrons and
in case of holes. Substiting into (1) gives
At low charge carrier concentrations, Fermi-Dirac distribution is simplified as
The approximation is considered valid when .
Approximate charge carrier concentrations are
(6)¶
Gaussian density of states¶
In the case of Gaussian DOS, the density of states shape function is the Gaussian distribution function scaled by total density of states :
Concentrations of species are given by integral
Conservation equation¶
The conservation equation is:
where denotes the concentration, is time, and is the flux density. denotes source term, which is positive for generating particles, and negative for sinking particles of type i. The SI unit of source term is .
The conservation equation must be satisfied for each species separately. In the case of transport of electrons and holes, this gives
(7)¶
where the source S term contains for example generation and recombination terms
The conservation of electric charge must be satisfied everywhere. Therefore, the source terms acting at given point must not create a net electric charge. In the case of system of electron and holes, this requires
Current density¶
Current density is related to the density flux by the charge of single particle . Obviously, for electrons and for holes , therefore
(8)¶
Note that a convention is adopted to denote the electric current with uppercase letter , and the flux density with lowercase letter . The SI unit of density flux is , while the unit of electric current density is .
Equilibrium conditions¶
In the equilibrium conditions, Fermi level energy has the same value everywhere. The electrostatic potential can vary, and the density of free charge does not need to be zero. Equations (1), (2), (4) are satisfied simultaneously. The current flux, the source terms, and the time dependence are all zeros, so conservation (7) is trivially satisfied.
Nonequilibrium conditions¶
In the non-equilibrium conditions, the transport is introduced as a perturbation from equilibrium. The Fermi energy level is replaced with quasi Fermi level, which is different for each species. In (1), the equilibrium Fermi level for electrons is replaced with a quasi Fermi level . Similarly,, the equilibrium Fermi level for holes is replaced wuth quasi Fermi level for holes , giving
(9)¶
Quasi Fermi levels have associated quasi Fermi potential according to the formula for energy of an electron in electrostatic field :
(10)¶
The transport is modeled by approximating electric current density as
(11)¶
where denotes the respective mobilities. The SI unit of mobility is , although is often used.
Equations (2), (5), (9), (11), (7) are simultaneously satisfied in non-equilibrium conditions.
Drift-diffusion system¶
Standard form of density fluxes in the drift-diffusion system is
(12)¶
or more generally, allowing arbitrary charge per particle
(13)¶
is the diffusion coefficient, with SI unit .
Drift-diffusion system: low concentration limit¶
To obtain the conventional drift-diffusion formulation (12), the the low concentration approximation (6) should be used. After introducing quasi Fermi levels, as it is done in (9), one obtains
From that, the quasi Fermi energies are calculated as
Using (2), and assuming constant ionization potential , bandgap , constant total densities of states , and constant temperature , substituting into (11), and using (3)
In terms of density flux (8), this reads
(14)¶
where thermal voltage
Einstein’s relation¶
Equation (14) is written in the standard drift-diffusion form (12) when the diffusion coefficient satisfies
(15)¶
This is called Einstein’s relation.
Drift-diffusion system: general case¶
Using functions defined in (1), bands (2) and approximation (9)
current densities under assumptions are
(16)¶
Generalized Einstein’s relation¶
In equation (16), assuming
In order to express equation (16) in the standard drift-diffusion form (12), the diffusion coefficient must satisfy
This is so called generalized Einstein’s relation .
Intrinsic concentrations¶
Intrinsic concentrations , , and intrinsic Fermi level satisfy electric neutrality conditions
Unidimensional form¶
By substituting and , the equations (5), (7), (12) of the basic drift-diffusion device model are
Total electric current density¶
Total electric current is a sum of currents due to transport of each species and the displacement current
Total electric current satisfies the conservation law
This can be verified by taking time derivative (4), using (7) and considering that the sum of all charge created by the source terms must be zero.
Metal¶
In metal, the relation between the electrostatic potential , the workfunction energy and the Fermi level is
On the other hand, the Fermi potential corresponds to the applied voltage
This leads to electrostatic potential at metal surface
Ohmic contact¶
Ohmic contact is an idealization assuming that there is no charge accumulation at the contact, and the applied voltage is equal to quasi Fermi potentials (10) of charged species
Above three conditions uniquely determine the charge concentrations , , and the electrostatic potential at the contact.
Electrochemical transport¶
Electrochemical potential for ionic species is
It should be noticed that so defined “potential” has the unit of energy, unlike the electrostatic potential and quasi Fermi potentials. Above denote corrections, for example due to steric interactions. Electrochemical potential should not be confused with mobility .
Density flux is approximated as
yielding the standard form (13) using Einstein’s relation (15).
Electrochemical species should be included in Poisson’s equation, by including proper source terms of form . A variant of Poisson’s equation (5) where are free charges are ions can be written as
Steric corrections¶
To account for finite size of ions, the electrochemical potential in the form introduced in [LE13] is useful
where denotes volume of particle of type . is the unoccupied fraction of space
where summing is taken over all species occupying space, including solvent.
[LE13] | Jinn-Liang Liu and Bob Eisenberg. Correlated Ions in a Calcium Channel Model: A Poisson–Fermi Theory. The Journal of Physical Chemistry B, 117(40):12051–12058, 2013. PMID: 24024558. URL: https://doi.org/10.1021/jp408330f, arXiv:https://doi.org/10.1021/jp408330f, doi:10.1021/jp408330f. |
Optical models¶
Coherent transfer matrix method¶
Transfer matrix method a convenient way of modeling thin film stacks. It is assumed that layers are stacked along axis, with being interface between layer and layer . Optical properties of each layer are specified by wavelength dependent complex refraction coefficient .
Optical field inside layer at given point along axis is specified by column vector , with being complex amplitude of forward traveling wave, and being complex amplitude of backward traveling wave.
Snell law is determines angles of propagation in each layer
where index 0 refers to medium before first layer. is angle of illuminating wave. All angles can be complex numbers. Since is multivalued function, angle of forward traveling wave is found from conditions that forward wave has forward pointing Poynting vector, or alternatively, that the amplitude of forward wave decays in absorbing medium.
In this convention, interface between layers is described by matrix as
with entries of matrix specified as
where transmission coefficient and reflection coefficient are given by Fresnel equations for complex amplitudes of light passing from layer i to layer i+1. Coefficients for backward propagating wave and are eliminated using Stokes relations.
For s-polarized wave:
For p-polarized wave:
Propagation inside layer is described by matrix as
with phase-shift
Light entering layer , on side of layer has vector of complex amplitudes
with vector denoting light leaving the device on the side opposite to illumination, with being complex amplitude of transmitted wave.
Applying above to whole device gives
with amplitude of illuminating wave set arbitrarily to and being complex amplitude of reflected wave.
When analyzing stack, firstly, solution , is found. Then intensity of light anywhere inside the device is calculated using found vectors and propagation matrices . Total intensity is found by applying Poynting formula. Absorbed energy is found by differentiating with respect to .
Incoherent light¶
Incoherent light is described by spectrum . Absorption of incoherent light is calculated as
where is calculated using coherent transfer matrix method.
Examples¶
oedes comes with examples demonstrating typical use. Each example shows usage of oedes and is used as a test suite to check for correctness. Many examples additionally partially or completely reproduce published papers. They therefore allow to experiment with published models.
The examples provided as files with extension .ipynb
run in Jupyter Notebook. Before running the examples, oedes should be installed.
The examples are included in oedes source distribution, and can also be downloaded from oedes repository on github. Also, they can be browsed on oedes website.