17. Exchange coupling in a double quantum dot in FD-SOI—Part 1: Perturbation theory
17.1. Requirements
17.1.1. Software components
QTCAD
Gmsh
17.1.2. Geometry file
qtcad/examples/tutorials/meshes/dqdfdsoi.geo
17.1.3. Python script
qtcad/examples/tutorials/exchange_1.py
17.1.4. References
17.2. Briefing
In this tutorial, we will calculate the exchange interaction strength for a double quantum dot in the generic fully-depleted silicon-on-insulator (FD-SOI) geometry presented in A double quantum dot device in a fully-depleted silicon-on-insulator transistor. Here, we use the perturbation theory presented in Perturbation theory of exchange coupling in a double quantum dot to quickly evaluate the exchange interaction strength without requiring the costly exact diagonalization approach which will be presented in the next tutorial.
As explained in Perturbation theory of exchange coupling in a double quantum dot, in the symmetric configuration, the exchange interaction strength may be estimated within perturbation theory using only two quantities that are easily accessible through single-electron calculations: the tunnel coupling strength \(\Omega\) and the on-site Coulomb interaction energy \(U\). While the tunnel coupling strength was already considered in previous tutorials, here, we will evaluate the on-site Coulomb interaction energy. This will enable us to compute an estimate for the exchange interaction that we will later compare with an exact diagonalization calculation.
17.3. Header
We start by importing the necessary packages, classes, and functions.
import numpy as np
import pathlib
from qtcad.device import constants as ct
from qtcad.device import io
from qtcad.device.mesh3d import Mesh, SubMesh
from qtcad.device import SubDevice
from qtcad.device.poisson import SolverParams as PoissonSolverParams
from qtcad.device.poisson import Solver as PoissonSolver
from qtcad.device.schrodinger import SolverParams as SchrodingerSolverParams
from qtcad.device.schrodinger import Solver as SchrodingerSolver
from qtcad.device.many_body import Solver as ManyBodySolver
from qtcad.device.many_body import SolverParams as ManyBodySolverParams
from examples.tutorials.double_dot_fdsoi import get_double_dot_fdsoi
These packages, classes, and functions were all seen in previous tutorials.
17.4. Setting up the device
As in Tunnel coupling in a double quantum dot in FD-SOI—Part 1: Plunger gate tuning, we will set up the device using the
get_double_dot_fdsoi
function defined in A double quantum dot device in a fully-depleted silicon-on-insulator transistor.
We first set up paths for the mesh file and for output files.
# Set up paths
script_dir = pathlib.Path(__file__).parent.resolve()
# mesh path
path_mesh = script_dir / "meshes"
path_out = script_dir / "output"
We then load the strongest tunnel coupling that was calculated in Tunnel coupling in a double quantum dot in FD-SOI—Part 2: Tuning the barrier gate and saved into a text file.
# Load the tunnel coupling from strong-tunneling simulation
path_tunnel = path_out / "tunnel_coupling_low_barrier.txt"
tunnel_coupling = np.loadtxt(path_tunnel)
This tunnel coupling is one of the two parameters required to calculate exchange within the perturbation theory explained in Perturbation theory of exchange coupling in a double quantum dot. In the rest of this tutorial, we will thus focus on calculating the on-site Coulomb interaction \(U\).
The next step is to load the mesh. To expedite convergence of the non-linear Poisson solver, we load the mesh that was already obtained in Tunnel coupling in a double quantum dot in FD-SOI—Part 1: Plunger gate tuning thanks to adaptive refinement.
# Load the mesh
scaling = 1e-9
path_file = str(path_mesh / "refined_dqdfdsoi.msh")
mesh = Mesh(scaling, path_file)
We then set the gate biases. Here, instead of aiming at a symmetric configuration, we introduce a detuning between the plunger gates that is large enough to ensure that the first two single-electron eigenstates of the double dot are localized in the left and right dots, but still too small to distort the confinement potential in a way that would lead to Coulomb integrals that are not representative of the localized left and right states. After some trial and error, we find that a detuning of \(5~\textrm{meV}\) is adequate.
# Define the gate biases
detuning = 5e-3
back_gate_bias = -0.5
barrier_gate_1_bias = 0.5
plunger_gate_1_bias = 0.59
barrier_gate_2_bias = 0.57
plunger_gate_2_bias = 0.59 + detuning
barrier_gate_3_bias = 0.5
We then define the device and create a submesh object for the double quantum dot region.
# Define the device
dvc = get_double_dot_fdsoi(mesh, back_gate_bias, barrier_gate_1_bias,
plunger_gate_1_bias, barrier_gate_2_bias, plunger_gate_2_bias,
barrier_gate_3_bias)
# List of regions forming the double quantum dot region
dot_region_list = ["oxide_dot", "gate_oxide_dot", "buried_oxide_dot",
"channel_dot"]
# Create the submesh object for the dot region
submesh = SubMesh(dvc.mesh, dot_region_list)
17.5. Basic electrostatics and eigenstates in a localized configuration
Since the gate bias configuration is still close to the one employed in Tunnel coupling in a double quantum dot in FD-SOI—Part 1: Plunger gate tuning, we start by loading the final electric potential that was saved after solving non-linear Poisson in that previous run.
# Load the electric potential from the previous run to speed up Poisson
# calculation
phi = io.load(path_out / "tunnel_coupling_final_potential.hdf5",
var_name="phi")
dvc.set_potential(phi)
Using this potential as an initial guess will help speed up convergence of
the non-linear Poisson
Solver
.
We configure the non-linear Poisson and Schrödinger solver parameters with
# Configure the non-linear Poisson solver
params_poisson = PoissonSolverParams()
params_poisson.tol = 1e-3 # Convergence threshold (tolerance) for the error
# Configure the Schrödinger solver
params_schrod = SchrodingerSolverParams()
params_schrod.tol = 1e-6 # Tolerance on energies in eV
We then instantiate the Poissson Solver
and run its solve
method with
the initialize
argument set to False
to use the potential
that is already stored in the
Device
object as an initial guess.
# Instantiate Poisson solver
poisson_slv = PoissonSolver(dvc, solver_params=params_poisson)
# Solve Poisson's equation
poisson_slv.solve(initialize=False)
We then set the confinement potential from the solution of the non-linear
Poisson equation, create a SubDevice
object for the double quantum dot region, instantiate a Schrödinger
Solver
, run it, print the
resulting single-electron eigenenergies, and save the ground and first
excited state wavefunctions.
# Get the potential energy from the band edge for usage in the Schrodinger
# solver
dvc.set_V_from_phi()
# Create a subdevice for the dot region
subdvc = SubDevice(dvc, submesh)
# Create a Schrodinger solver, and run it
schrod_solver = SchrodingerSolver(subdvc, solver_params=params_schrod)
schrod_solver.solve()
subdvc.print_energies()
# Save the ground and first excited state wavefunctions into .vtu files
path_states = path_out / "exchange_localized_wavefunctions.vtu"
out_dict = {
"Ground state": subdvc.eigenfunctions[:,0],
"First excited state": subdvc.eigenfunctions[:,1]
}
io.save(path_states, out_dict, submesh)
Plotting the wavefunctions in ParaView, we find that the ground and first excited states are well localized in the right and left dots, respectively.
17.6. Coulomb integrals
Now that well-localized single-electron states have been found, we may readily calculate their Coulomb integrals to evaluate the on-site Coulomb interaction.
To do so, we instantiate a many-body
Solver
object.
# Instantiate a many-body solver
manybody_solver_params = ManyBodySolverParams()
manybody_solver_params.num_states = 2
manybody_solver_params.overlap = False
manybody_solver = ManyBodySolver(subdvc, solver_params=manybody_solver_params)
We set num_states
to 2
, because we are only interested in the Coulomb
integrals matrix in the basis spanned by the localized left and right states of
the double quantum dot. In addition, we set overlap
to False
, because
we are not interested in overlap integrals such as \(V_\mathrm{LRRL}\).
We then compute the Coulomb integrals matrix, print it, and save it in a text file.
# Compute the matrix of Coulomb integrals
coulomb_mat = manybody_solver.get_coulomb_matrix(overlap = False)
print('Coulomb interaction matrix in a nearly localized single-electron basis set')
print(coulomb_mat[0:2,0:2]/ct.e)
# Export the coulomb matrix
np.savetxt(path_out/"exchange_coulomb_localized.txt", coulomb_mat)
We obtain the following results.
Coulomb interaction matrix in a nearly localized single-electron basis set
[[0.01723048 0.00406259]
[0.00406259 0.0168071 ]]
Knowing from the above analysis that the ground state is localized in the right dot (corresponding to \(|\mathrm R\rangle\)), and that the first excited state is localized in the left dot (corresponding to \(|\mathrm L\rangle\)), we may associate this matrix with
Therefore, we have found \(V_L \approx 16.8~\textrm{meV}\), \(V_R \approx 17.2~\textrm{meV}\), and \(V_\mathrm{LR} = V_\mathrm{RL} \approx 4.1~\textrm{meV}\).
We remark that the inter-site Coulomb interaction, while significant, is a factor of four smaller than the on-site Coulomb interaction, and may thus be neglected in a first approximation, when interest mostly lies in the order of magnitude of the exchange interaction strength. In addition, the asymmetric bias used here to localize the eigenstates has introduced a small difference (\(\approx 0.5~\textrm{meV}\)) between \(V_\mathrm L\) and \(V_\mathrm R\) which we will neglect here.
17.7. Exchange interaction strength
Knowing both the tunnel coupling and Coulomb interaction strengths, we are ready to evaluate the exchange interaction strength using Eq. (9.8).
# Calculate the theoretical exchange coupling
exchange = (tunnel_coupling*ct.e)**2 / coulomb_mat[0,0]
print("Theoretical exchange interaction: {:.3f} MHz".format(exchange/ct.h/1e6))
# Export the theoretical exchange coupling
np.savetxt(path_out/"exchange_perturbation.txt", np.array([exchange]))
In addition to printing the value of the exchange interaction strength, we save it in a text file (in SI units) for further use in the next tutorial.
Expressing the exchange interaction strength in MHz, we find
Theoretical exchange interaction: 918.090 MHz
In the next tutorial, Exchange coupling in a double quantum dot in FD-SOI—Part 2: Exact diagonalization, we will compare this value with the output of an exact diagonalization approach.
17.8. Full code
__copyright__ = "Copyright 2024, Nanoacademic Technologies Inc."
import numpy as np
import pathlib
from qtcad.device import constants as ct
from qtcad.device import io
from qtcad.device.mesh3d import Mesh, SubMesh
from qtcad.device import SubDevice
from qtcad.device.poisson import SolverParams as PoissonSolverParams
from qtcad.device.poisson import Solver as PoissonSolver
from qtcad.device.schrodinger import SolverParams as SchrodingerSolverParams
from qtcad.device.schrodinger import Solver as SchrodingerSolver
from qtcad.device.many_body import Solver as ManyBodySolver
from qtcad.device.many_body import SolverParams as ManyBodySolverParams
from examples.tutorials.double_dot_fdsoi import get_double_dot_fdsoi
# Set up paths
script_dir = pathlib.Path(__file__).parent.resolve()
# mesh path
path_mesh = script_dir / "meshes"
path_out = script_dir / "output"
# Load the tunnel coupling from strong-tunneling simulation
path_tunnel = path_out / "tunnel_coupling_low_barrier.txt"
tunnel_coupling = np.loadtxt(path_tunnel)
# Load the mesh
scaling = 1e-9
path_file = str(path_mesh / "refined_dqdfdsoi.msh")
mesh = Mesh(scaling, path_file)
# Define the gate biases
detuning = 5e-3
back_gate_bias = -0.5
barrier_gate_1_bias = 0.5
plunger_gate_1_bias = 0.59
barrier_gate_2_bias = 0.57
plunger_gate_2_bias = 0.59 + detuning
barrier_gate_3_bias = 0.5
# Define the device
dvc = get_double_dot_fdsoi(mesh, back_gate_bias, barrier_gate_1_bias,
plunger_gate_1_bias, barrier_gate_2_bias, plunger_gate_2_bias,
barrier_gate_3_bias)
# List of regions forming the double quantum dot region
dot_region_list = ["oxide_dot", "gate_oxide_dot", "buried_oxide_dot",
"channel_dot"]
# Create the submesh object for the dot region
submesh = SubMesh(dvc.mesh, dot_region_list)
# Load the electric potential from the previous run to speed up Poisson
# calculation
phi = io.load(path_out / "tunnel_coupling_final_potential.hdf5",
var_name="phi")
dvc.set_potential(phi)
# Configure the non-linear Poisson solver
params_poisson = PoissonSolverParams()
params_poisson.tol = 1e-3 # Convergence threshold (tolerance) for the error
# Configure the Schrödinger solver
params_schrod = SchrodingerSolverParams()
params_schrod.tol = 1e-6 # Tolerance on energies in eV
# Instantiate Poisson solver
poisson_slv = PoissonSolver(dvc, solver_params=params_poisson)
# Solve Poisson's equation
poisson_slv.solve(initialize=False)
# Get the potential energy from the band edge for usage in the Schrodinger
# solver
dvc.set_V_from_phi()
# Create a subdevice for the dot region
subdvc = SubDevice(dvc, submesh)
# Create a Schrodinger solver, and run it
schrod_solver = SchrodingerSolver(subdvc, solver_params=params_schrod)
schrod_solver.solve()
subdvc.print_energies()
# Save the ground and first excited state wavefunctions into .vtu files
path_states = path_out / "exchange_localized_wavefunctions.vtu"
out_dict = {
"Ground state": subdvc.eigenfunctions[:,0],
"First excited state": subdvc.eigenfunctions[:,1]
}
io.save(path_states, out_dict, submesh)
# Instantiate a many-body solver
manybody_solver_params = ManyBodySolverParams()
manybody_solver_params.num_states = 2
manybody_solver_params.overlap = False
manybody_solver = ManyBodySolver(subdvc, solver_params=manybody_solver_params)
# Compute the matrix of Coulomb integrals
coulomb_mat = manybody_solver.get_coulomb_matrix(overlap = False)
print('Coulomb interaction matrix in a nearly localized single-electron basis set')
print(coulomb_mat[0:2,0:2]/ct.e)
# Export the coulomb matrix
np.savetxt(path_out/"exchange_coulomb_localized.txt", coulomb_mat)
# Calculate the theoretical exchange coupling
exchange = (tunnel_coupling*ct.e)**2 / coulomb_mat[0,0]
print("Theoretical exchange interaction: {:.3f} MHz".format(exchange/ct.h/1e6))
# Export the theoretical exchange coupling
np.savetxt(path_out/"exchange_perturbation.txt", np.array([exchange]))