4. Equation and Models¶
4.1. Overview¶
DEVSIM
uses the control volume approach for assembling partialdifferential equations (PDE’s) on the simulation mesh. DEVSIM
is used to solve equations of the form:
Internally, it transforms the PDE’s into an integral form.
Equations involving the divergence operators are converted into surface integrals, while other components are integrated over the device volume.
In Fig. 4.1, 2D mesh elements are depicted. The shaded area around the center node is referred to as the node volume, and it is used for the volume integration. The lines from the center node to other nodes are referred to as edges. The flux through the edge are integrated with respect to the perpendicular bisectors (dashed lines) crossing each triangle edge.
In this form, we refer to a model integrated over the edges of triangles as edge models. Models integrated over the volume of each triangle vertex are referred to as node models. Element edge models are a special case where variables at other nodes off the edge may cause the flux to change.
There are a default set of models created in each region upon initialization of a device, and are typically based on the geometrical attributes. These are described in the following sections. Models required for describing the device behavior are created using the equation parser described in SYMDIFF. For special situations, custom matrix assembly is also available and is discussed in Custom Matrix assembly.
4.2. Bulk models¶
4.2.1. Node models¶
Node models may be specified in terms of other node models, mathematical functions, and parameters on the device. The simplest model is the node solution, and it represents the solution variables being solved for. Node models automatically created for a region are listed in Table 4.1.
In this example, we present an implementation of Shockley Read Hall recombination [MKC02].
USRH="ElectronCharge*(Electrons*Holes  n_i^2)/(taup*(Electrons + n1) \
+ taun*(Holes + p1))")
dUSRHdn="simplify(diff(%s, Electrons))" % USRH
dUSRHdp="simplify(diff(%s, Holes))" % USRH
devsim.node_model(device='MyDevice', region='MyRegion',
name="USRH", equation=USRH)
devsim.node_model(device='MyDevice', region='MyRegion',
name="USRH:Electrons", equation=dUSRHdn)
devsim.node_model(device='MyDevice', region='MyRegion',
name="USRH:Holes", equation=dUSRHdp)
The first model specified, USRH
, is the recombination model itself. The derivatives with respect to electrons and holes are USRH:Electrons
and USRH:Holes
, respectively. In this particular example Electrons
and Holes
have already been defined as solution variables. The remaining variables in the equation have already been specified as parameters.
The diff
function tells the equation parser to take the derivative of the original expression, with respect to the variable specified as the second argument. During equation assembly, these derivatives are required in order to converge upon a solution.
The simplify
function tells the expression parser to attempt to simplify the expression as much as possible.
Node Model 
Description 


Evaluates to 1 if node is a contact node, otherwise 0 

The volume of the node. Used for volume integration of node models on nodes in mesh 

The surface normal to points on the interface or contact (2D and 3D) 

The surface normal to points on the interface or contact (2D and 3D) 

The surface normal to points on the interface or contact (3D) 

The surface area of a node on interface nodes, otherwise 0 

The surface area of a node on contact nodes, otherwise 0 

Coordinate index of the node on the device 

Index of the node in the region 

x position of the node 

y position of the node 

z position of the node 
4.2.2. Edge models¶
Edge models may be specified in terms of other edge models, mathematical functions, and parameters on the device. In addition, edge models may reference node models defined on the ends of the edge. As depicted in Fig. 4.2, edge models are with respect to the two nodes on the edge, n0
and n1
.
For example, to calculate the electric field on the edges in the region, the following scheme is employed:
devsim.edge_model(device="device", region="region", name="ElectricField",
equation="(Potential@n0  Potential@n1)*EdgeInverseLength")
devsim.edge_model(device="device", region="region",
name="ElectricField:Potential@n0", equation="EdgeInverseLength")
devsim.edge_model(device="device", region="region",
name="ElectricField:Potential@n1", equation="EdgeInverseLength")
In this example, EdgeInverseLength
is a builtin model for the inverse length between nodes on an edge. Potential@n0
and Potential@n1
is the Potential
node solution on the nodes at the end of the edge. These edge quantities are created using the devsim.edge_from_node_model()
. In addition, the devsim.edge_average_model()
can be used to create edge models in terms of node model quantities.
Edge models automatically created for a region are listed in Table 4.2.
Edge Model 
Description 


The length of the perpendicular bisector of an element edge. Used to perform surface integration of edge models on edges in mesh. 

The volume for each node on an edge. Used to perform volume integration of edge models on edges in mesh. 

Inverse of the EdgeLength. 

The distance between the two nodes of an edge 

Index of the edge on the region 

x component of the unit vector along an edge 

y component of the unit vector along an edge (2D and 3D) 

z component of the unit vector along an edge (3D only) 
4.2.3. Element edge models¶
Element edge models are used when the edge quantitites cannot be specified entirely in terms of the quantities on both nodes of the edge, such as when the carrier mobility is dependent on the normal electric field. In 2D, element edge models are evaluated on each triangle edge. As depicted in Fig. 4.3, edge models are with respect to the three nodes on each triangle edge and are denoted as en0
, en1
, and en2
. Derivatives are with respect to each node on the triangle.
In 3D, element edge models are evaluated on each tetrahedron edge. Derivatives are with respect to the nodes on both triangles on the tetrahedron edge. Element edge models automatically created for a region are listed in Table 4.3.
As an alternative to treating integrating the element edge model with respect to ElementEdgeCouple
, the integration may be performed with respect to ElementNodeVolume
. See devsim.equation()
for more information.
Element Edge Model 
Description 


The length of the perpendicular bisector of an edge. Used to perform surface integration of element edge model on element edge in the mesh. 

The node volume at either end of each element edge. 
4.2.4. Model derivatives¶
To converge upon the solution, derivatives are required with respect to each of the solution variables in the system. DEVSIM
will look for the required derivatives. For a model model
, the derivatives with respect to solution variable variable
are presented in Table 4.4.
Model Type 
Derivatives Required 

Node Model 

Edge Model 

Element Edge Model 

4.2.5. Conversions between model types¶
The devsim.edge_from_node_model()
is used to create edge models referring to the nodes connecting the edge. For example, the edge models Potential@n0
and Potential@n1
refer to the Potential
node model on each end of the edge.
The devsim.edge_average_model()
creates an edge model which is either the arithmetic mean, geometric mean, gradient, or negative of the gradient of the node model on each edge.
When an edge model is referred to in an element edge model expression, the edge values are implicity converted into element edge values during expression evaluation. In addition, derivatives of the edge model with respect to the nodes of an element edge are required, they are converted as well. For example, edgemodel:variable@n0
and edgemodel:variable@n1
are implicitly converted to edgemodel:variable@en0
and edgemodel:variable@en1
, respectively.
The devsim.element_from_edge_model()
is used to create directional components of an edge model over an entire element. The derivative
option is used with this command to create the derivatives with respect to a specific node model. The devsim.element_from_node_model()
is used to create element edge models referring to each node on the element of the element edge.
4.2.6. Equation assembly¶
Bulk equations are specified in terms of the node, edge, and element edge models using the devsim.equation()
. Node models are integrated with respect to the node volume. Edge models are integrated with the perpendicular bisectors along the edge onto the nodes on either end.
Element edge models are treated as flux terms and are integrated with respect to ElementEdgeCouple
using the element_model
option. Alternatively, they may be treated as source terms and are integrated with respect to ElementNodeVolume
using the volume_node0_model
and volume_node1_model
option.
In this example, we are specifying the Potential Equation in the region to consist of a flux term named PotentialEdgeFlux
and to not have any node volume terms.
devsim.equation(device="device", region="region", name="PotentialEquation",
variable_name="Potential", edge_model="PotentialEdgeFlux",
variable_update="log_damp" )
In addition, the solution variable coupled with this equation is Potential
and it will be updated using logarithmic damping.
Model Type 
Model Name 
Derivatives Required 

Node Model (region 0) 


Node Model (region 1) 


Interface Node Model 


4.3. Interface¶
4.3.1. Interface models¶
Fig. 4.4 depicts an interface in DEVSIM
. It is a collection of overlapping nodes existing in two regions, r0
and r1
.
Interface models are node models specific to the interface being considered. They are unique from bulk node models, in the sense that they may refer to node models on both sides of the interface. They are specified using the devsim.interface_model()
. Interface models may refer to node models or parameters on either side of the interface using the syntax nodemodel@r0
and nodemodel@r1
to refer to the node model in the first and second regions of the interface. The naming convention for node models, interface node models, and their derivatives are shown in Table 4.5.
devsim.interface_model(device="device", interface="interface",
name="continuousPotential", equation="Potential@r0Potential@r1")
4.3.2. Interface model derivatives¶
For a given interface model, model
, the derivatives with respect to the variable variable
in the regions are
model:variable@r0
model:variable@r1
devsim.interface_model(device="device", interface="interface",
name="continuousPotential:Potential@r0", equation="1")
devsim.interface_model(device="device", interface="interface",
name="continuousPotential:Potential@r1", equation="1")
4.3.3. Interface equation assembly¶
There are three types of interface equations considered in DEVSIM
. They are both activated using the devsim.interface_equation()
.
In the first form, continuous
, the equations for the nodes on both sides of the interface are integrated with respect to their volumes and added into the same equation. An additional equation is then specified to relate the variables on both sides. In this example, continuity in the potential solution across the interface is enforced, using the continuousPotential
model defined in the previous section.
devsim.interface_equation(device="device", interface="interface", name="PotentialEquation",
variable_name="Potential", interface_model="continuousPotential",
type="continuous")
In the second form, fluxterm
, a flux term is integrated over the surface area of the interface and added to the first region, and subtracted from the second.
In the third form, hybrid
, equations for nodes on both sides of the interface are added into the equation for the node in the first region. The equation for the node on the second interface is integrated in the second region, and the fluxterm is subracted in the second region.
4.4. Contact¶
4.4.1. Contact models¶
Fig. 4.5 depicts how a contact is treated in a simulation. It is a collection of nodes on a region. During assembly, the specified models form an equation, which replaces the equation applied to these nodes for a bulk node.
Contact models are equivalent to node and edge models, and are specified using the devsim.contact_node_model()
and the devsim.contact_edge_model()
, respectively. The key difference is that the models are only evaluated on the contact nodes for the contact specified.
4.4.2. Contact model derivatives¶
The derivatives are equivalent to the discussion in Model derivatives. If external circuit boundary conditions are being used, the model model
derivative with respect to the circuit node node
name should be specified as model:node
.
4.4.3. Contact equation assembly¶
The devsim.contact_equation()
is used to specify the boundary conditions on the contact nodes. The models specified replace the models specified for bulk equations of the same name. For example, the node model specified for the contact equation is assembled on the contact nodes, instead of the node model specified for the bulk equation. Contact equation models not specified are not assembled, even if the model exists on the bulk equation for the region attached to the contact.
As an example
devsim.contact_equation(device="device", contact="contact", name="PotentialEquation",
variable_name="Potential", node_model="contact_bc",
edge_charge_model="DField")
Current models refer to the instantaneous current flowing into the device. Charge models refer to the instantaneous charge at the contact.
During a transient, smallsignal or ac simulation, the time derivative is taken so that the net current into a circuit node is
where \(i\) is the integrated current and \(q\) is the integrated charge.
4.5. Custom Matrix assembly¶
The devsim.custom_equation()
command is used to register callbacks to be called during matrix and right hand side assembly. The Python
procedure should expect to receive two arguments and return two lists and a boolean value. For example a procedure named myassemble
registered with
devsim.custom_equation(name="test1", procedure="myassemble")
expects two arguments
def myassemble(what, timemode):
.
.
.
return rcv, rv, True
where what
may be passed as one of
MATRIXONLYRHSMATRIXANDRHS
and timemode
may be passed as one of
DCTIME
When timemode
is DC
, the timeindependent part of the equation is returned.
When timemode
is TIME
, the timederivative part of the equation is returned. The simulator will scale the timederivative terms with the proper frequency or time scale.
The return value from the procedure must return two lists and a boolean value of the form
[1 1 1.0 2 2 1.0 1 2 1.0 2 1 1.0 2 2 1.0], [1 1.0 2 1.0 2 1.0], True
where the length of the first list is divisible by 3 and contains the row, column, and value to be assembled into the matrix. The second list is divisible by 2 and contains the right hand side entries. Either list may be empty.
The boolean value denotes whether the matrix and right hand side entries should be row permutated. A value of True
should be used for assembling bulk equations, and a value of False
should be used for aseembling contact and interface boundary conditions.
The devsim.get_circuit_equation_number()
may be used to get the equation numbers corresponding to circuit node names. The devsim.get_equation_numbers()
may be used to find the equation number corresponding to each node index in a region.
The matrix and right hand side entries should be scaled by the NodeVolume
if they are assembled into locations in a device region as volume integration.
4.6. Cylindrical Coordinate Systems¶
In 2D, models representing the edge couples, surface areas and node volumes may be generated using the following commands:
In order to change the integration from the default models to cylindrical models, the following parameters may be set
set_parameter(name="node_volume_model",
value="CylindricalNodeVolume")
set_parameter(name="edge_couple_model",
value="CylindricalEdgeCouple")
set_parameter(name="edge_node0_volume_model",
value="CylindricalEdgeNodeVolume@n0")
set_parameter(name="edge_node1_volume_model",
value="CylindricalEdgeNodeVolume@n1")
set_parameter(name="element_edge_couple_model",
value="ElementCylindricalEdgeCouple")
set_parameter(name="element_node0_volume_model",
value="ElementCylindricalNodeVolume@en0")
set_parameter(name="element_node1_volume_model",
value="ElementCylindricalNodeVolume@en1")