Select Git revision
drive_gate_voltage.py 8.25 KiB
# Copyright 2022 Valentin Bruch <valentin.bruch@rwth-aachen.de>
# License: MIT
"""
Kondo FRTRG, variant of kondo model for oscillating coupling
Oscillations in the coupling between a quantum dot and its reservoirs can be
the result of an oscillating gate voltage in a quantum dot. This module
includes a basic implementation of this idea. The RG flow is not converging!
"""
from kondo import *
class KondoG(Kondo):
"""
Variant of Kondo where absolute value of J oscillates.
This may be the result of an oscillating gate voltage for a quantum dot in
the Kondo regime or of direct driving of the tunneling barrier.
"""
def __init__(self, j_driving_param, **options):
super().__init__(**options)
self.j_driving_param = j_driving_param
def initialize(self,
**solveopts : 'keyword arguments passed to solver',
):
"""
Arguments:
**solveopts: keyword arguments passed to the solver. Most relevant
are rtol and atol.
Get initial conditions for Γ, Z and G2 by numerically solving the
equilibrium RG equations from E=0 to E=iD and for all required Re(E).
Initialize G3, Iγ, δΓ, δΓγ.
"""
assert self.compact == 0
sqrtxx = np.sqrt(self.xL*(1-self.xL))
symmetry = 0 if settings.IGNORE_SYMMETRIES else 1
#### Initial conditions from exact results at T=V=0
# Get Γ, Z and J (G2) for T=V=0.
gamma0, z0, j0 = solveTV0_Utransformed(d=self.d, properties=self.global_properties, **solveopts)
# Write T=V=0 results to Floquet index n=0.
gammavalues = np.zeros(self.shape(), dtype=np.complex128)
zvalues = np.zeros(self.shape(), dtype=np.complex128)
jvalues = np.zeros(self.shape(), dtype=np.complex128)
# construct diagonal matrices
diag_idx = (..., *np.diag_indices(2*self.nmax+1))
if self.resonant_dc_shift:
gammavalues[diag_idx] = gamma0[...,self.resonant_dc_shift:-self.resonant_dc_shift]
zvalues[diag_idx] = z0[...,self.resonant_dc_shift:-self.resonant_dc_shift]
jvalues[diag_idx] = j0[...,self.resonant_dc_shift:-self.resonant_dc_shift]
else:
gammavalues[diag_idx] = gamma0
zvalues[diag_idx] = z0
jvalues[diag_idx] = j0
# Create Γ and Z with the just calculated initial values.
self.gamma = RGfunction(self.global_properties, gammavalues, symmetry=symmetry)
self.z = RGfunction(self.global_properties, zvalues, symmetry=symmetry)
# Create G2 from J: G2_{ij} = - 2 sqrt(x_i x_j) J
self.g2 = ReservoirMatrix(self.global_properties, symmetry=symmetry)
j_rgfunction = RGfunction(self.global_properties, jvalues, symmetry=symmetry)
self.g2[0,0] = -2*self.xL * j_rgfunction
self.g2[1,1] = -2*(1-self.xL) * j_rgfunction
init_matrix = gen_init_matrix(self.nmax, 0, resonant_dc_shift=self.resonant_dc_shift)
j_driving_scaled = self.j_driving_param * j0
init_matrix[np.arange(0, 2*self.nmax), np.arange(1, 2*self.nmax+1)] = j_driving_scaled[:-1]
init_matrix[np.arange(1, 2*self.nmax+1), np.arange(0, 2*self.nmax)] = j_driving_scaled[1:].conjugate()
j_LR = np.einsum('ij,...j->...ij', init_matrix, j0[...,2*self.resonant_dc_shift:])
j_RL = np.einsum('ji,...j->...ij', init_matrix.conjugate(), j0[...,:j0.shape[-1]-2*self.resonant_dc_shift])
j_LR = RGfunction(self.global_properties, j_LR)
j_RL = RGfunction(self.global_properties, j_RL)
self.g2[0,1] = -2*sqrtxx * j_LR
self.g2[1,0] = -2*sqrtxx * j_RL
## Initial conditions for G3
# G3 ~ Jtilde^2 with Jtilde = Z J
# Every entry of G3 will be of the following form (up to prefactors):
self.g3 = ReservoirMatrix(self.global_properties, symmetry=-symmetry)
g3_entry = np.zeros(self.shape(), dtype=np.complex128)
g3_entry[diag_idx] = 1j*np.pi * (jvalues[diag_idx]*zvalues[diag_idx])**2
g3_entry = RGfunction(self.global_properties, g3_entry, symmetry=-symmetry)
self.g3[0,0] = 2*self.xL * g3_entry
self.g3[1,1] = 2*(1-self.xL) * g3_entry
g30 = 1j*np.pi*(z0*j0)**2
g3_LR = np.einsum('ij,...j->...ij', init_matrix, g30[...,2*self.resonant_dc_shift:])
g3_RL = np.einsum('ji,...j->...ij', init_matrix.conjugate(), g30[...,:g30.shape[-1]-2*self.resonant_dc_shift])
g3_LR = RGfunction(self.global_properties, g3_LR)
g3_RL = RGfunction(self.global_properties, g3_RL)
self.g3[0,1] = 2*sqrtxx * g3_LR
self.g3[1,0] = 2*sqrtxx * g3_RL
## Initial conditions for current I^{γ=L} = J0 (1 - Jtilde)
if self.voltage_branches:
current_entry = np.diag( 2*sqrtxx * jvalues[self.voltage_branches][diag_idx] * (1 - jvalues[self.voltage_branches][diag_idx] * zvalues[self.voltage_branches][diag_idx] ) )
else:
current_entry = np.diag( 2*sqrtxx * jvalues[diag_idx] * (1 - jvalues[diag_idx] * zvalues[diag_idx] ) )
current_entry = RGfunction(self.global_properties, current_entry, symmetry=-symmetry)
self.current = ReservoirMatrix(self.global_properties, symmetry=-symmetry)
self.current[0,0] = 0*current_entry
self.current[0,0].symmetry = -symmetry
self.current[1,1] = self.current[0,0].copy()
if self.voltage_branches:
i0 = j0[self.voltage_branches] * (1 - j0[self.voltage_branches]*z0[self.voltage_branches])
else:
i0 = j0 * (1 - j0*z0)
i_LR = np.einsum('ij,...j->...ij', init_matrix, i0[2*self.resonant_dc_shift:])
i_RL = np.einsum('ji,...j->...ij', init_matrix.conjugate(), i0[:i0.size-2*self.resonant_dc_shift])
i_LR = RGfunction(self.global_properties, i_LR)
i_RL = RGfunction(self.global_properties, i_RL)
self.current[0,1] = 2*sqrtxx * i_LR
self.current[1,0] = -2*sqrtxx * i_RL
## Initial conditions for voltage-variation of Γ: δΓ
self.deltaGamma = RGfunction(
self.global_properties,
np.zeros((3,2*self.nmax+1,2*self.nmax+1) if self.voltage_branches else self.shape(), dtype=np.complex128),
symmetry = symmetry
)
## Initial conditions for voltage-variation of current-Γ: δΓ_L
self.deltaGammaL = RGfunction(
self.global_properties,
np.zeros((2*self.nmax+1, 2*self.nmax+1), dtype=np.complex128),
symmetry = symmetry
)
if self.resonant_dc_shift:
if self.voltage_branches:
self.deltaGammaL.values[diag_idx] = 3*np.pi*sqrtxx**2 * (j0[self.voltage_branches,self.resonant_dc_shift:-self.resonant_dc_shift]*z0[self.voltage_branches,self.resonant_dc_shift:-self.resonant_dc_shift])**2
else:
self.deltaGammaL.values[diag_idx] = 3*np.pi*sqrtxx**2 * (j0[self.resonant_dc_shift:-self.resonant_dc_shift]*z0[self.resonant_dc_shift:-self.resonant_dc_shift])**2
else:
if self.voltage_branches:
diag_values = 3*np.pi*sqrtxx**2 * (j0[self.voltage_branches]*z0[self.voltage_branches])**2
else:
diag_values = 3*np.pi*sqrtxx**2 * (j0*z0)**2
self.deltaGammaL.values[diag_idx] = diag_values
del diag_values
### Derivative of full current
self.yL = RGfunction(
self.global_properties,
np.zeros((2*self.nmax+1,2*self.nmax+1), dtype=np.complex128),
symmetry=-symmetry
)
### Full current, also includes AC current
self.gammaL = (self.vdc + self.omega*self.resonant_dc_shift) * self.deltaGammaL.reduced()
if self.voltage_branches:
gammaL_AC = 3*np.pi*sqrtxx**2 * self.vac/2 * (j0[self.voltage_branches]*z0[self.voltage_branches])**2
else:
gammaL_AC = 3*np.pi*sqrtxx**2 * self.vac/2 * (j0*z0)**2
if self.resonant_dc_shift:
gammaL_AC = gammaL_AC[...,self.resonant_dc_shift:-self.resonant_dc_shift]
idx = (np.arange(1, 2*self.nmax+1), np.arange(2*self.nmax))
self.gammaL.values[idx] = gammaL_AC[...,1:]
idx = (np.arange(0, 2*self.nmax), np.arange(1, 2*self.nmax+1))
self.gammaL.values[idx] = gammaL_AC[...,:-1]
self.global_properties.energy = 1j*self.d