6
$\begingroup$

I have a single-rectangular potential barrier in graphene. The transmission $T(E,\theta)$ is given by $$T(E,\theta) = \left[1+\left(\frac{V}{\hbar v_Fk_x}\right)^2\tan^2\theta\sin^2(k_xD)\right]^{-1}$$ where $k_x=\sqrt{(E-V)^2-E^2\sin^2\theta}/\hbar v_F$. I would like to plot conductance $G(E)$ versus energy $E$, with $E_{min}=-1 \, \mathrm{eV}, E_{max}=1.5 \, \mathrm{eV}$, using the Landauer-Buttiker approach: $$G(E)=\int_{-\pi/2}^{\pi/2}T(E,\theta)\cos\theta d\theta.$$ The constants are $V=200\,\mathrm{meV}$, $D=30\,\mathrm{nm}$ and $v_F = 10^8\,\mathrm{cm/s}$. More details about the model can be found in Physica E, 61 (2014) 118.

$\endgroup$
3
  • $\begingroup$ What are the units of E? The units of $k_x$ must also be meV, so how does one create a dimensionless quantity from 1 meV * 30 nm? $\endgroup$
    – Carl Woll
    Sep 26, 2018 at 20:24
  • $\begingroup$ In the original article, the formula contains the Planck constant and the characteristic velocity $v_F=10^8 cm/s$. $\endgroup$ Sep 27, 2018 at 5:18
  • $\begingroup$ Exactly, I considered $\hbar=v_F=1$ for simplicity. $\endgroup$ Sep 27, 2018 at 15:47

2 Answers 2

12
$\begingroup$

We reconstruct the original formula from the paper with dimensional units, as a result, an additional dimensionless coefficient $k$ appears in the formula, and the units of measurement of energy $\mathrm{eV}$ and the length of $\mathrm{nm}$.

ℏ = 1.054572*^-27; 
vf = 1*^8;
d0 = 1*^-7;
eV = 1.6021766208*^-12;
k = eV*d0/ℏ/vf;

T[x_, E0_, V_, d_] := ((E0 - V)^2 - E0^2 Sin[x]^2)/((E0 - V)^2 - 
  E0^2 Sin[x]^2 + V^2 Sin[k*d Sqrt[(E0 - V)^2 - E0^2 Sin[x]^2]]^2 Tan[x]^2)

DensityPlot[
 T[x, E0, .2, 30], {E0, -.1, .5}, {x, -π/2 + 10^-2, π/2 - 10^-2}, 
 PlotPoints -> 100, PerformanceGoal -> "Quality", PlotRange -> All, 
 FrameLabel -> {"E, eV", "θ"}, PlotLabel -> T]

G[E0_, V_, d_] := 
NIntegrate[Cos[x]*T[x, E0, V, d], {x, -Pi/2, Pi/2}, 
  Method -> "DoubleExponential", PrecisionGoal -> 8, 
  MinRecursion -> 4]   


lst = Table[{E0, G[E0, .2, 30]}, {E0, -1, 1.5, .001}];

ListLinePlot[lst, AxesLabel -> {"E, eV", "G"}]

fig1

$\endgroup$
4
  • 2
    $\begingroup$ This is pretty quick: With[{V = 1/5, d = 30}, Plot[NIntegrate[Cos[x] T[x, E0, V, d], {x, -Pi/2, Pi/2}, AccuracyGoal -> 8, Method -> "Trapezoidal"], {E0, -1, 3/2}]]. $\endgroup$ Sep 27, 2018 at 7:32
  • $\begingroup$ @J. M. Thank you, very good solution of the problem. $\endgroup$ Sep 27, 2018 at 8:53
  • $\begingroup$ Thank you so much, @AlexTrounev, You helped me beyond what I had asked! $\endgroup$ Sep 27, 2018 at 15:53
  • $\begingroup$ @Alex Nascimento You're welcome. I'm glad that this code turned out to be useful for you. $\endgroup$ Sep 27, 2018 at 16:05
7
$\begingroup$

The following is basically @Alex's answer, but making use of the units framework, and using a style more similar to the question.

First, define:

V = .2;
d = 30;

Next, determine the conversion factor using Ctrl-= to input the units (I used an image to see how it looks in a notebook):

enter image description here

(Alternatively, type in the FullForm equivalent):

conv = (Quantity["Electronvolts"]Quantity["Nanometers"])/(Quantity["ReducedPlanckConstant"]Quantity[10^8,"Centimeters"/"Seconds"])

1.5192675

(another shorter version that relies on Quantity interpretations)

Quantity["eV"] Quantity["nm"]/(Quantity["ℏ"] Quantity[10^8,"cm/s"])

1.519267

Then, define k and T:

k[e_, θ_] := Sqrt[(e-V)^2 - e^2 Sin[θ]^2]
T[e_,θ_] := (1 + V^2/k[e,θ]^2 Tan[θ]^2 Sin[k[e,θ] d conv]^2)^-1

Next, define G but restrict to only numeric arguments:

G[e_?NumericQ] := NIntegrate[T[e,θ] Cos[θ], {θ, -Pi/2, Pi/2}]

Finally, plot:

Plot[G[e], {e, -1, 1.5}, AxesOrigin->{0,0}]

enter image description here

$\endgroup$
2
  • $\begingroup$ thanks for the great help, but I'd like to know where the unit conversion factor comes in the code. $\endgroup$ Sep 27, 2018 at 16:24
  • $\begingroup$ @AlexNascimento See update. $\endgroup$
    – Carl Woll
    Sep 27, 2018 at 16:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.