# Equation and Models¶

## Overview¶

`DEVSIM`

uses the control volume approach for assembling partial-differential 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. 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.

## Bulk models¶

### 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 2.

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 |
---|---|

`AtContactNode` |
Evaluates to 1 if node is a contact node, otherwise 0 |

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

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

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

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

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

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

`coordinate_index` |
Coordinate index of the node on the device |

`node_index` |
Index of the node in the region |

`x` |
x position of the node |

`y` |
y position of the node |

`z` |
z position of the node |

### 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. 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 built-in 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 3.

Edge Model | Description |
---|---|

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

`EdgeInverseLength` |
Inverse of the EdgeLength. |

`EdgeLength` |
The distance between the two nodes of an edge |

`edge_index` |
Index of the edge on the region |

`unitx` |
x component of the unit vector along an edge |

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

`unitz` |
z component of the unit vector along an edge (3D only) |

### 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. 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.

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 |
---|---|

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

`ElementNodeVolume` |
The node volume at either end of each element edge. |

### 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 5.

Model Type | Derivatives Required |
---|---|

Node Model | `model:variable` |

Edge Model | `model:variable@n0` , `model:variable@n1` |

Element Edge Model | `model:variable@en0` , `model:variable@en1` , `model:variable@en2` , `model:variable@en3` (3D) |

### 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.

### 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_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) | `nodemodel@r0` |
`nodemodel:variable@r0` |

Node Model (region 1) | `nodemodel@r1` |
`nodemodel:variable@r1` |

Interface Node Model | `inodemodel` |
`inodemodel:variable@r0` , `inodemodel:variable@r1` |

## Interface¶

### Interface models¶

Fig. 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 6.

```
devsim.interface_model(device="device", interface="interface",
name="continuousPotential", equation="Potential@r0-Potential@r1")
```

### 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")
```

### 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.

## Contact¶

### Contact models¶

Fig. 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.

### 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`

.

### 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, small-signal 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.

## 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 must expect to receive two arguments and return two lists. For example a procedure named `myassemble`

registered with

```
devsim.custom_equation(name="test1", procedure="myassemble")
```

must expect to receive two arguments

```
def myassemble(what, timemode):
.
.
.
return [rcv, rv]
```

where `what`

may be passed as one of

MATRIXONLYRHSMATRIXANDRHS

and `timemode`

may be passed as one of

DCTIME

When `timemode`

is `DC`

, the time-independent part of the equation is returned.
When `timemode`

is `TIME`

, the time-derivative part of the equation is returned. The simulator will scale the time-derivative terms with the proper frequency or time scale.

The return value from the procedure must return two lists 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]

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 `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. Row permutations, required for contact and interface boundary conditions, are automatically applied to the row numbers returned by the `Python`

procedure.

## 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")
```