3. Simulating electrostatics using the linear Poisson solver

3.1. Requirements

3.1.1. Software components

• QTCAD

3.1.2. Python script

• qtcad/epractical_application/FDSOI/2-linear_poisson.py

3.1.3. References

• Mesh generation using QTCAD Builder

• Creating the FD-SOI device

3.2. Briefing

The goal of this part of the tutorial is to simulate the electrostatics over the fully-depleted silicon on insulator (FD-SOI) device created in the previous part of the tutorial Mesh generation using QTCAD Builder using the linear Poisson Solver. In particular, we will demonstrate how to leverage the adaptive meshing capabilities of this solver to automatically refine the mesh in regions that disproportionately contribute to the error in the potential phi. This mesh will therefore provide

• a strategically refined mesh over which we can perform subsequent simulations

• a solution to the linear Poisson equation which can be used as an initial guess for subsequent calculations.

3.3. Setup

3.3.1. Header

We start by importing the necessary Python modules

import pathlib

from qtcad.device import constants as ct
from qtcad.device.mesh3d import Mesh
from qtcad.device.poisson_linear import Solver as LinearPoissonSolver
from qtcad.device.poisson_linear import SolverParams as LinearPoissonSolverParams

from qtcad.device import io
from qtcad.device.analysis import plot_slice

from double_dot_fdsoi import get_double_dot_fdsoi, normal, origin

The main module relevant to this tutorial is the LinearPoissonSolver. Another important import is the get_double_dot_fdsoi function from double_dot_fdsoi.py which was introduced in Creating the FD-SOI device. This function creates and returns an instance of the Device class representing the FD-SOI device. Specifically:

• It creates the Device object over the Mesh generated in the previous part of the tutorial, Mesh generation using QTCAD Builder.

• It assigns the relevant materials to the different parts of the device.

• It applies the appropriate boundary conditions to the different surfaces of the device.

This function outputs the Device object representing the FD-SOI device along with two lists:

• List of region names corresponding to the single-electron transistor (SET).

• List of region names corresponding to the qubit quantum-dot region.

These lists will be used later when we solve the Schrödinger equation. For more details on the geometry of the device, see Mesh generation using QTCAD Builder and Creating the FD-SOI device

Note

In this tutorial, although the source and drain regions are included in the device geometry, we do not solve for the charges within them. Typically, source and drain regions act as reservoirs of semiclassical electrons and can be modeled by filling parabolic bands up to the Fermi energy according to the Fermi-–Dirac distribution. Explicit modeling of these regions can be carried out using the nonlinear PoissonSolver.

3.3.2. File paths

Next we define the paths to input and output files for the solver and its results:

# Files -----------------------------------------------------------------------
script_dir           = pathlib.Path(__file__).parent.resolve()

path_mesh            = script_dir / "meshes"
path_out             = script_dir / "output"

mesh_file            = str(path_mesh / f"dqdfdsoi.msh")
geo_file             = str(path_mesh / f"dqdfdsoi.xao")
ref_file             = str(path_mesh / f"refined_dqdfdsoi.msh")
pos_file             = str(path_mesh / f"dqdfdsoi.pos")

# output files
phi_out_file         = str(path_out / "phi" / f"phi.hdf5")
CB_out_file          = str(path_out / "phi" / f"CB.vtu")
CB_slice_file        = str(path_out / "phi" / f"CB_slice.png")

3.4. Creating the device

After the setup is complete, we can create the FD-SOI device using the get_double_dot_fdsoi function imported above.

# Create the device -----------------------------------------------------------
d, set_region, qd_region = get_double_dot_fdsoi(mesh_file_name="dqdfdsoi.msh")

The device is created over the mesh generated in Mesh generation using QTCAD Builder. The voltages applied to the gates controlling the qubit and the SET are set to the default values defined in double_dot_fdsoi.py.

3.5. Solving the linear Poisson equation

Now that the device has been created, we can proceed to set up and solve the linear Poisson equation using the LinearPoissonSolver. We start by defining the solver parameters using the SolverParams class.

# Solve Linear Poisson equation -----------------------------------------------

# Solver parameters
p_adapt_params                  = LinearPoissonSolverParams()

p_adapt_params.tol              = 1e-10
p_adapt_params.eta0             = 0.10
p_adapt_params.refined_region   = qd_region + set_region
p_adapt_params.h_refined        = 0.8
p_adapt_params.refined_mesh_filename = ref_file

The paramaters defined above are as follows:

• tol: The target tolerance (in \(\mathrm{V}\)) for the iterative solver used to invert the linear Poisson equation.

• eta0: The target relative accuracy of the potential. The mesh will continue to be adapted until the estimated relative accuracy of the potential reaches eta0 and the mesh is well-behaved in all regions.

• refined_region: A list of region names over which a minimal characteristic length will be enforced. This minimal characteristic length is specified by the h_refined parameter. Here the region names correspond to the qubit and SET regions defined when creating the device.

• h_refined: The minimal characteristic length in the refined regions specified via the refined_region parameter.

• refined_mesh_filename: The path to the file where the refined mesh will be saved. This refined Mesh will be reused in subsequent simulations to avoid repeating the adaptive meshing procedure.

Now, we solve the linear Poisson equation using

# Solve the linear Poisson equation
p_adapt_solver = LinearPoissonSolver(d=d, solver_params=p_adapt_params, geo_file=geo_file)
p_adapt_solver.solve()

Here, the LinearPoissonSolver is initialized using the Device object d, the solver parameters defined above, and the geometry file generated in the previous part of the tutorial Mesh generation using QTCAD Builder. Passing the geometry file triggers the adaptive meshing procedure while solving the linear Poisson equation using the solve method. This procedure will iteratively refine the mesh until the target relative accuracy specified by the eta0 parameter is reached and store the refined mesh in the file specified by the refined_mesh_filename parameter.

3.6. Saving the results

After solving the linear Poisson equation, we can save the results in various formats. Typically, the .hdf5 format is used to save the potential for later use. We will later load this potential and use it as an initial guess when solving the Schrödinger–Poisson equation in Simulating the Coulomb peak positions of an SET. We also use the .vtu format to visualize the conduction-band edge in Paraview, see Visualizing QTCAD quantities with ParaView.

# Save results ----------------------------------------------------------------

# Save the Poisson equation solution in .hdf5 format
arrays_dict = {"phi": d.phi}
io.save(phi_out_file, arrays_dict)

# Save the conduction band edge in .vtu format
CB = d.cond_band_edge() / ct.e
io.save(CB_out_file, {"EC [eV]": CB}, d.mesh)
# Save a slice of the conduction band edge in .png format
label_CB = "Conduction-band edge [eV]"
plot_slice(d.mesh, CB, normal=normal, origin=origin, cb_axis_label=label_CB,
    path=CB_slice_file, show_figure=False)

For quick visualization, we also save (in .png format) a slice of the conduction-band edge perpendicular to the growth (z) direction, 1 nm below the gate oxides using the plot_slice. function. The generated slice is shown below:

Fig. 3.6.1 Slice of the conduction-band edge taken perpendicular to the \(z\) direction, \(1\,\mathrm{nm}\) below the gate oxides. A discontinuity appears in the conduction-band edge at the silicon–silicon-oxide interface, which leads to visible artifacts in the plot due to the abrupt change in material properties and the interpolation used in the visualization.

The conduction band edge shows minima under the plunger gates, with a deeper minimum under the SET plunger gate where we expect the accumulation of many electrons.

3.7. Full code

__copyright__ = "Copyright 2022-2025, Nanoacademic Technologies Inc."

import pathlib
from qtcad.device import constants as ct
from qtcad.device.mesh3d import Mesh
from qtcad.device.poisson_linear import Solver as LinearPoissonSolver
from qtcad.device.poisson_linear import SolverParams as LinearPoissonSolverParams
from qtcad.device import io
from qtcad.device.analysis import plot_slice
from double_dot_fdsoi import get_double_dot_fdsoi, normal, origin

# Files -----------------------------------------------------------------------
script_dir           = pathlib.Path(__file__).parent.resolve()

path_mesh            = script_dir / "meshes"
path_out             = script_dir / "output"

mesh_file            = str(path_mesh / f"dqdfdsoi.msh")
geo_file             = str(path_mesh / f"dqdfdsoi.xao")
ref_file             = str(path_mesh / f"refined_dqdfdsoi.msh")
pos_file             = str(path_mesh / f"dqdfdsoi.pos")

# output files
phi_out_file         = str(path_out / "phi" / f"phi.hdf5")
CB_out_file          = str(path_out / "phi" / f"CB.vtu")
CB_slice_file        = str(path_out / "phi" / f"CB_slice.png")

# Create the device -----------------------------------------------------------
d, set_region, qd_region = get_double_dot_fdsoi(mesh_file_name="dqdfdsoi.msh")

# Solve Linear Poisson equation -----------------------------------------------

# Solver parameters 
p_adapt_params                  = LinearPoissonSolverParams()

p_adapt_params.tol              = 1e-10
p_adapt_params.eta0             = 0.10
p_adapt_params.refined_region   = qd_region + set_region
p_adapt_params.h_refined        = 0.8
p_adapt_params.refined_mesh_filename = ref_file

# Solve the linear Poisson equation
p_adapt_solver = LinearPoissonSolver(d=d, solver_params=p_adapt_params, geo_file=geo_file)
p_adapt_solver.solve()

# Save results ----------------------------------------------------------------

# Save the Poisson equation solution in .hdf5 format
arrays_dict = {"phi": d.phi}
io.save(phi_out_file, arrays_dict)

# Save the conduction band edge in .vtu format
CB = d.cond_band_edge() / ct.e
io.save(CB_out_file, {"EC [eV]": CB}, d.mesh)
# Save a slice of the conduction band edge in .png format
label_CB = "Conduction-band edge [eV]"
plot_slice(d.mesh, CB, normal=normal, origin=origin, cb_axis_label=label_CB,
    path=CB_slice_file, show_figure=False)