diff --git a/docs/src/assets/four_terminal_T.png b/docs/src/assets/four_terminal_T.png new file mode 100644 index 00000000..5993b90e Binary files /dev/null and b/docs/src/assets/four_terminal_T.png differ diff --git a/docs/src/assets/four_terminal_g_big.png b/docs/src/assets/four_terminal_g_big.png new file mode 100644 index 00000000..31212b1f Binary files /dev/null and b/docs/src/assets/four_terminal_g_big.png differ diff --git a/docs/src/tutorial/greenfunctions.md b/docs/src/tutorial/greenfunctions.md index 021ae29f..d208e084 100644 --- a/docs/src/tutorial/greenfunctions.md +++ b/docs/src/tutorial/greenfunctions.md @@ -139,9 +139,9 @@ Let us define the classical example of a multiterminal mesoscopic junction. We c We first define a single lead `Greenfunction` and the central Hamiltonian ```julia -julia> glead = LP.square() |> hopping(1) |> supercell((1, 0), region = r -> abs(r[2]) <= 5/2) |> greenfunction(GS.Schur(boundary = 0)); +julia> glead = LP.square() |> onsite(4) - hopping(1) |> supercell((1, 0), region = r -> abs(r[2]) <= 5/2) |> greenfunction(GS.Schur(boundary = 0)); -julia> hcentral = LP.square() |> hopping(1) |> supercell(region = RP.circle(10) | RP.rectangle((22, 5)) | RP.rectangle((5, 22))) +julia> hcentral = LP.square() |> onsite(4) - hopping(1) |> supercell(region = RP.circle(10) | RP.rectangle((22, 5)) | RP.rectangle((5, 22))) ``` The two rectangles overlayed on the circle above create the stubs where the leads will be attached: ```@raw html @@ -250,7 +250,7 @@ This particular example made use of the `GS.Bands` solver (the default for `L>1` julia> qplot(bands(g)) ``` ```@raw html -Bands used by the GS.Bands solver +Bands used by the GS.Bands solver ``` These bands can be adjusted by passing arguments such as ``` @@ -265,12 +265,12 @@ julia> h = LP.square() |> onsite(4) - hopping(1) |> supercell(region = r -> norm julia> g = h |> attach(@onsite(ω -> -im), region = r -> r[1] ≈ 47) |> greenfunction; -julia> g_1_to_all = sum(abs2, g(0.1)[siteselector(), 1], dims = 2); +julia> gx1 = sum(abs2, g(0.1)[siteselector(), 1], dims = 2); -julia> qplot(h, hide = :hops, sitecolor = (i, r) -> gω[i], siteradius = (i, r) -> gω[i], minmaxsiteradius = (0, 2), sitecolormap = :balance) +julia> qplot(h, hide = :hops, sitecolor = (i, r) -> gx1[i], siteradius = (i, r) -> gx1[i], minmaxsiteradius = (0, 2), sitecolormap = :balance) ``` ```@raw html -Green function from a contact on the right +Green function from a contact on the right ``` !!! warning "Caveat for multiorbital systems" diff --git a/docs/src/tutorial/observables.md b/docs/src/tutorial/observables.md index f3aa80af..23be712b 100644 --- a/docs/src/tutorial/observables.md +++ b/docs/src/tutorial/observables.md @@ -5,7 +5,7 @@ We are almost at our destination now. We have defined a `Lattice`, a `Model` for Currently, we have the following observables built into Quantica.jl (with more to come in the future) - `ldos`: computes the local density of states at specific energy and sites - `current`: computes the local current density along specific directions, and at specific energy and sites - - `transmission`: computes the transmission probability between contacts + - `transmission`: computes the total transmission between contacts - `conductance`: computes the differential conductance `dIᵢ/dVⱼ` between contacts `i` and `j` - `josephson`: computes the supercurrent and the current-phase relation through a given contact in a superconducting system @@ -26,7 +26,7 @@ LocalSpectralDensitySolution{Float64} : local density of states at fixed energy julia> qplot(h, hide = :hops, sitecolor = ρ, siteradius = ρ, minmaxsiteradius = (0, 2), sitecolormap = :balance) ``` ```@raw html -LDOS +LDOS ``` Note that `ρ[sites...]` produces a vector with the LDOS at sites defined by `siteselector(; sites...)` (`ρ[]` is the ldos over all sites). We can also define a `kernel` to be traced over orbitals to obtain the spectral density of site-local observables (see `diagonal` slicing in the preceding section). @@ -47,13 +47,56 @@ CurrentDensitySolution{Float64} : current density at a fixed energy and arbitrar julia> qplot(h, siteradius = 0.08, sitecolor = :black, siteoutline = 0, hopradius = J, hopcolor = J, minmaxhopradius = (0, 2), hopcolormap = :balance, hopdarken = 0) ``` ```@raw html -Current density with magnetic flux +Current density with magnetic flux ``` + !!! note "Remember to construct supercell before applying position-dependent fields" Note that we built the supercell before applying the model with the magnetic flux. Not doing so would make the gauge field be repeated in each unit cell when expanding the supercell. This was mentioned in the section on Hamiltonians, and is a common mistake when modeling systems with position dependent fields. ## Transmission +The transmission `Tᵢⱼ` from contact `j` to contact `i` can be computed using `transmission`. This function accepts a `GreenSlice` between the contact. Let us recover the four-terminal setup of the preceding section, but let's make it bigger this time +```julia +julia> hcentral = LP.square() |> onsite(4) - hopping(1) |> supercell(region = RP.circle(100) | RP.rectangle((202, 50)) | RP.rectangle((50, 202)))^C + +julia> glead = LP.square() |> onsite(4) - hopping(1) |> supercell((1, 0), region = r -> abs(r[2]) <= 50/2) |> greenfunction(GS.Schur(boundary = 0)); + +julia> Rot = r -> SA[0 -1; 1 0] * r; # 90º rotation function + +julia> g = hcentral |> + attach(glead, region = r -> r[1] == 101) |> + attach(glead, region = r -> r[1] == -101, reverse = true) |> + attach(glead, region = r -> r[2] == 101, transform = Rot) |> + attach(glead, region = r -> r[2] == -101, reverse = true, transform = Rot) |> + greenfunction; + +julia> gx1 = sum(abs2, g(0.04)[siteselector(), 1], dims = 2); + +julia> qplot(hcentral, hide = :hops, siteoutline = 1, sitecolor = (i, r) -> gx1[i], siteradius = (i, r) -> gx1[i], minmaxsiteradius = (0, 2), sitecolormap = :balance) +``` +```@raw html +Green function from right lead +``` + +It's apparent from the plot that the transmission from right to left (`T₂₁` here) at this energy of `0.04` is larger than from right to top (`T₃₁`). Is this true in general? Let us compute the two transmissions as a function of energy. To show the progress of the calculation we can use a monitor package, such as `ProgressMeter` +```julia +julia> using ProgressMeter + +julia> T₂₁ = transmission(g[2,1]); T₃₁ = transmission(g[3,1]); ωs = subdiv(0, 4, 200); + +julia> T₂₁ω = @showprogress [T₂₁(ω) for ω in ωs]; T₃₁ω = @showprogress [T₃₁(ω) for ω in ωs]; +Progress: 100%|████████████████████████████████████████████████████████████████████████████████████| Time: 0:01:05 +Progress: 100%|████████████████████████████████████████████████████████████████████████████████████| Time: 0:01:04 + +julia> f = Figure(); a = Axis(f[1,1], xlabel = "ω/t", ylabel = "T(ω)"); lines!(a, ωs, T₂₁ω, label = L"T_{2,1}"); lines!(a, ωs, T₃₁ω, label = L"T_{3,1}"); axislegend("Transmission", position = :lt); f +``` +```@raw html +Total transmission from right contact +``` + +!!! note "Total transmission vs transmission probability" + Note that `transmission` gives the total transmission, which is the sum of the transmission probability from each orbital in the source contact to any other orbital in the drain contact. As such it is not normalized to 1, but to the number of source orbitals. It also gives the local conductance from a given contact in units of $$e^2/h$$ according to the Landauer formula, $$G\_j = e^2/h \sum_i T_{ij}(eV)$$. + ## Conductance ## Josephson diff --git a/src/show.jl b/src/show.jl index e728addb..710eac7f 100644 --- a/src/show.jl +++ b/src/show.jl @@ -355,7 +355,7 @@ $i To contact : $(currentcontact(parent(T)))") end Base.summary(::Transmission) = - "Transmission: total transmission probability between two different contacts" + "Transmission: total transmission between two different contacts" #endregion @@ -432,4 +432,4 @@ Base.summary(::CurrentDensitySolution{T}) where {T} = Base.summary(::CurrentDensitySlice{T}) where {T} = "CurrentDensitySlice{$T} : current density at a fixed location and arbitrary energy" -#endregion \ No newline at end of file +#endregion