469 lines
31 KiB
TeX
469 lines
31 KiB
TeX
\documentclass{ian}
|
|
|
|
\usepackage{graphicx}
|
|
\usepackage{array}
|
|
|
|
\begin{document}
|
|
|
|
\hbox{}
|
|
\hfil{\bf\LARGE
|
|
Exact solution of the Schr\"odinger equation\par
|
|
\hfil for photoemission from a metal\par
|
|
}
|
|
\vfill
|
|
|
|
\hfil{\bf\large Ovidiu Costin}\par
|
|
\hfil{\it Department of Mathematics, The Ohio State University}\par
|
|
\hfil{\tt\color{blue}\href{mailto:costin.9@osu.edu}{costin.9@osu.edu}}\par
|
|
\vskip20pt
|
|
|
|
\hfil{\bf\large Rodica Costin}\par
|
|
\hfil{\it Department of Mathematics, The Ohio State University}\par
|
|
\hfil{\tt\color{blue}\href{mailto:costin.10@osu.edu}{costin.10@osu.edu}}\par
|
|
\vskip20pt
|
|
|
|
\hfil{\bf\large Ian Jauslin}\par
|
|
\hfil{\it Department of Physics, Princeton University}\par
|
|
\hfil{\tt\color{blue}\href{mailto:ijauslin@princeton.edu}{ijauslin@princeton.edu}}\par
|
|
\vskip20pt
|
|
|
|
\hfil{\bf\large Joel L. Lebowitz}\par
|
|
\hfil{\it Department of Mathematics and Physics, Rutgers University}\par
|
|
\hfil{\tt\color{blue}\href{mailto:lebowitz@math.rutgers.edu}{lebowitz@math.rutgers.edu}}\par
|
|
\vskip20pt
|
|
|
|
\vfill
|
|
|
|
|
|
\hfil {\bf Abstract}\par
|
|
\medskip
|
|
{\small
|
|
We solve rigorously the time dependent Schr\"odinger equation describing electron emission from a metal surface by a laser field perpendicular to the surface.
|
|
We consider the system to be one-dimensional, with the half-line $x<0$ corresponding to the bulk of the metal and $x>0$ to the vacuum.
|
|
The laser field is modeled as a classical electric field oscillating with frequency $\omega$, acting only at $x>0$.
|
|
We consider an initial condition which is a stationary state of the system without a field, and, at time $t=0$, the field is switched on.
|
|
We prove the existence of a solution $\psi(x,t)$ of the Schr\"odinger equation for $t>0$, and compute the surface current.
|
|
The current exhibits a complex oscillatory behavior, which is not captured by the ``simple'' three step scenario.
|
|
As $t\to\infty$, $\psi(x,t)$ converges with a rate $t^{-\frac32}$ to a time periodic function with period $\frac{2\pi}{\omega}$ which coincides with that found by Faisal, Kami\'nski and Saczuk (Phys Rev A {\bf 72}, 023412, 2015).
|
|
However, for realistic values of the parameters, we have found that it can take quite a long time (over 50 laser periods) for the system to converge to its asymptote.
|
|
Of particular physical importance is the current averaged over a laser period $\frac{2\pi}\omega$, which exhibits a dramatic increase when $\hbar\omega$ becomes larger than the work function of the metal, which is consistent with the original photoelectric effect.
|
|
}
|
|
|
|
\vfill
|
|
|
|
\tableofcontents
|
|
|
|
\vfill
|
|
\eject
|
|
|
|
\setcounter{page}1
|
|
\pagestyle{plain}
|
|
|
|
|
|
\section{Introduction}
|
|
\indent In this note, we continue our investigation of the exact solution of the time-dependent Schr\"odinger equation describing the emission of electrons from a metal surface by an electric field \cite{CCe18}.
|
|
We use the Sommerfeld model of quasi-free electrons with a Fermi distribution of energies, confined by a step potential $U=\mathcal E_F+W$, where $\mathcal E_F$ is the Fermi energy and $W$ is the work function of the metal.
|
|
This setup is commonly used as a simple model for the process of emission, both for a constant field and an oscillating field \cite{FN28,YGR11,Fo16,Je17}.
|
|
In both cases one focuses attention on electrons, part of the Fermi sea, with energy $\frac12k^2$ (in atomic units), moving from the left towards the metal surface at $x=0$.
|
|
These are described by a wave function $e^{ikx}$, $k>0$, $x<0$.
|
|
When the electrons reach the metal surface they are either reflected or transmitted through the barrier.
|
|
\bigskip
|
|
|
|
\indent The time evolution of the wave function of an electron in such a beam is described by the one dimensional Schr\"odinger equation,
|
|
\begin{equation}
|
|
i\partial_t\psi(x,t)=H(x,t)\psi(x,t)=-\frac12\Delta\psi(x,t)+\Theta(x)(U-Ex\cos(\omega t))\psi(x,t)
|
|
\label{schrodinger}
|
|
\end{equation}
|
|
where $\Theta(x)=0$ for $x<0$ and $\Theta(x)=1$ for $x>0$, $E$ is the electric field perpendicular to the surface, $\frac\omega{2\pi}$ is the frequency and we are using atomic units $\hbar=m=1$.
|
|
Here, we ignore the carrier wave and Shottky effect.
|
|
We also only focus on the direction that is orthogonal to the metal surface, ignoring the other two directions \cite{KSH11}.
|
|
|
|
\indent In the absence of an external field, $E=0$, the Schr\"odinger equation (\ref{schrodinger}) has a ``stationary" solution $e^{-i\frac12k^2t}\varphi_0(x)$ in which there is, for $k^2<2U$, a reflected beam of the same energy and intensity as the incoming beam $e^{ikx}$ and an evanescent, exponentially decaying tail on the right,
|
|
\begin{equation} \label{initial}
|
|
\varphi_0(x)=
|
|
\left\{\begin{array}{ll}
|
|
e^{ikx}+R_0e^{-ikx}&\mathrm{for\ }x<0\\
|
|
T_0e^{- \alpha x}&\mathrm{for\ }x>0
|
|
\end{array}\right.
|
|
\end{equation}
|
|
Using the continuity of $\psi$ and its derivative at $x=0$ yields
|
|
\begin{equation}
|
|
\alpha:=\sqrt{2U-k^2}
|
|
,\quad
|
|
1+R_0=T_0=\frac{2ik}{ik-\sqrt{2U-k^2}}
|
|
.
|
|
\end{equation}
|
|
The current,
|
|
\begin{equation}
|
|
j(x,t):=\mathcal Im(\psi^*(x,t)\partial_x\psi(x,t))
|
|
\end{equation}
|
|
is zero and no electrons leave the metal.
|
|
Note that $\varphi_0$ is not square integrable, which will complicate matters when it comes to solving the Schr\"odinger equation.
|
|
|
|
\indent If we turn on a constant field, $E>0, \ \omega=0$, at $t=0$, the Schr\"odinger equation has a ``stationary'' solution ${\rm e}^{-\frac 12 k^2 t}\varphi_{E}(x)$ with $\varphi_E$ satisfying the equation
|
|
\begin{equation}
|
|
\frac12k^2\varphi_E=-\frac12\frac{\partial^2}{\partial x^2}\varphi_E+\Theta(x)(U-Ex)\varphi_E
|
|
\end{equation}
|
|
with incident beam $\Theta(-x)e^{ikx}$.
|
|
This equation was solved by Fowler and Nordheim \cite{FN28} in 1928 to describe the quantum tunneling of electrons through a triangular barrier.
|
|
Their solution, with some corrections, is commonly used to describe the stationary current produced by a field $E$ acting on the surface of a metal \cite{Je17}.
|
|
\bigskip
|
|
|
|
\indent In \cite{CCe18} we solved the time-dependent Schr\"odinger equation (\ref{schrodinger}) with initial condition $\varphi_0(x)$ for the constant field.
|
|
We showed that $\psi(x,t)$ converges, as $t\to\infty$, to the Fowler-Nordheim (FN) solution ${\rm e}^{-\frac 12 k^2 t}\varphi_{E}(x)$
|
|
The exact rate of convergence is $t^{-\frac32}$.
|
|
The ``effective'' time of approach to the FN solution was found to be, for realistic values of the parameters, of the order of femtoseconds.
|
|
|
|
\indent Here, we investigate solutions of (\ref{schrodinger}) with $E>0$, $\omega>0$, $\psi(x,0)=\varphi_0(x)$ as given in (\ref{initial}).
|
|
Since we are interested in describing situations in which the laser is turned on for only a short time, with pulses as short as a few periods, understanding the role of the initial state is important.
|
|
This problem turns out to be considerably more difficult than the constant field case, $\omega=0$.
|
|
\bigskip
|
|
|
|
\indent There are several studies in the literature about the solution of the time-dependent Schr\"odinger equation, see \cite{KLe18} and references therein.
|
|
The method most often used is the Crank-Nicolson \cite{CN47} numerical scheme.
|
|
In the case treated here, however, as we will discuss in more detail in the following, because of the discontinuity in the potential at $x=0$, the Crank-Nicolson algorithm is actually inaccurate for short times.
|
|
We will introduce an alternative scheme that does not suffer from the discontinuity of the potential.
|
|
|
|
\indent Of particular note is the work by Yalunin, Gulden and Roper \cite{YGR11}, where the initial condition and potential are the same as in (\ref{schrodinger}) and (\ref{initial}).
|
|
Our result is quite different from theirs for short time, but the values we obtain at longer times seem to agree qualitatively with those obtained in \cite{YGR11} for short times.
|
|
Furthermore, in \cite{YGR11}, there does not seem to be a distinction between long and short times, whereas we find that for some important physical quantities, like the current averaged over a laser period, one has to wait many periods to be close to the asymptotic value.
|
|
\bigskip
|
|
|
|
\indent We look for an exact solution of (\ref{schrodinger}) for $t>0$ such that $\psi(x,t)$ and its derivative $\partial_x\psi(x,t)$ are continuous at $x=0$ and be bounded for $|x|\to\infty$.
|
|
To do this, we first obtain the solutions for $x<0$ and for $x>0$ with general initial conditions and unknown time dependent boundary conditions $\psi(0,t)=\psi_0(t)$ and $\partial_x\psi(x,t)|_{x=0}=\psi_{x,0}(t)$.
|
|
A solution of (\ref{schrodinger}) with appropriate $x\to\pm\infty$ conditions will exist if and only if we can match these boundary conditions coming from the right and from the left.
|
|
This approach to solving such an equation as (\ref{schrodinger}) goes under the general name of initial-boundary-value (IBV) problem, developed in particular by Fokas for these types of problems \cite{Fo97,Ry18}.
|
|
|
|
\indent We prove the existence of a solution of (\ref{schrodinger}) by showing that a necessary and sufficient condition for this is that $\psi_0(t)$ satisfies an integral equation.
|
|
The existence of the solution $\psi(x,t)$ for all $x$ and $t$ then follows from the proof that this integral equation has a unique solution, which we provide.
|
|
Expressions for $\psi$ and for the current $j(x,t)$ are given by explicit integrals involving $\psi_0(t)$ \cite{CCe}.
|
|
These can be evaluated numerically.
|
|
To do so accurately for both short and moderate times, we developed a numerical scheme which goes beyond the Crank-Nicolson algorithm \cite{CN47}, see section\,\ref{Sec3}.
|
|
\bigskip
|
|
|
|
\indent $\psi(x,t)$ has the long time form $e^{-i\frac12k^2t}\mathcal U(x,t)$, where $\mathcal U$ is (up to a phase) the periodic in time solution of (\ref{schrodinger}) derived by Faisal et al \cite{FKS05,ZL16}.
|
|
We note that Faisal et al. did not prove the existence of a solution to their infinite set of equations.
|
|
Instead, they solved numerically a truncated set of these equations.
|
|
Our results prove the existence of such a solution.
|
|
|
|
\indent The approach to this asymptotic form $\mathcal U$ goes like $t^{-\frac32}$.
|
|
$\mathcal{U}(x,t)$ consists of resonance terms which can be associated, as was done in \cite{FKS05}, with the absorption (emission) of different number of photons, see Section\,\ref{longtime}.
|
|
These resonances show up dramatically when considering the dependence of the average current on $\omega$ in the vicinity of $\omega_c=W$, as seen in early experiments on the photoelectric effect where light with frequency $\omega>\omega_c$ was necessary to get electron emission.
|
|
|
|
|
|
\section{Existence of solutions of the Schr\"odinger equation}
|
|
\indent Here we give an outline of the proof derived in \cite{CCe}. To simplify the notations we rescale: let
|
|
\begin{equation}\label{rescale}
|
|
\bar t=\omega t,\ \ \bar x=\sqrt{2\omega}x,\ \ \bar U=U\omega^{-1},\ \ \bar E=E\omega^{-3/2}2^{-1/2}
|
|
\end{equation}
|
|
Equation (\ref{schrodinger}) becomes
|
|
\begin{equation}\label{eq12}
|
|
i\partial_{\bar t}\psi(\bar x,\bar t)=-\partial^2_{\bar x}\psi(\bar x,\bar t)+\Theta(\bar x)\left( \bar U-2\bar E \bar x\cos(\bar t)\right)\psi(\bar x,\bar t)
|
|
\end{equation}
|
|
As mentioned earlier, we will solve (\ref{eq12}) with the general initial condition
|
|
$\psi(\bar x,0):=f(\bar x)$ separately for $\bar x<0$, $\psi_-(\bar x,\bar t)$ and for $\bar x>0$, $\psi_+(\bar x,\bar t)$, and then match the values of $\psi_\pm(0,\bar t)$ and $\partial_{\bar x}\psi_\pm(0,\bar t)$.
|
|
\bigskip
|
|
|
|
\indent To obtain these solutions we first go from $\bar x$ to Fourier space.
|
|
We write
|
|
\begin{equation}
|
|
\hat\psi(\xi,\bar t)=\frac1{\sqrt{2\pi}}\int_{-\infty}^\infty d\bar x\ e^{-i\bar x\xi}\psi(\bar x,\bar t)
|
|
=\hat\psi_-(\xi,\bar t)+\hat\psi_+(\xi,\bar t)
|
|
\end{equation}
|
|
where $\hat\psi_\pm$ is the half-line Fourier transform along $\pm \bar x>0$.
|
|
\bigskip
|
|
|
|
\indent The equation for $\hat\psi_-(\xi,\bar t)$ has the form
|
|
\begin{equation}\label{eq5}
|
|
i\frac{\partial \hat{\psi}_-}{\partial\bar t}=\xi^2 \hat{\psi}_-- \frac1{\sqrt{2\pi}} \psi_{\bar x,0}(\bar t)-i\xi \frac1{\sqrt{2\pi}} \psi_0(\bar t)
|
|
\end{equation}
|
|
where
|
|
\begin{equation}
|
|
\psi_0(\bar t):=\psi(0,\bar t)
|
|
,\quad
|
|
\psi_{\bar x,0}(\bar t):=\partial_{\bar x}\psi(0,\bar t)
|
|
\end{equation}
|
|
are treated as unknown functions to be determined later by the matching conditions.
|
|
The solution of (\ref{eq5}) is given by
|
|
\begin{equation} \label{psim}
|
|
\hat{\psi}_-(\xi,\bar t)= {\rm e}^{-i\xi^2\bar t}\left\{ C_-(\xi) + \frac1{\sqrt{2\pi}} \int_0^{\bar t} {\rm e}^{i\xi^2\bar s}\ \left[ i\psi_{\bar x,0}(\bar s)-\xi \psi_0(\bar s) \right]\, d\bar s\right\}
|
|
\end{equation}
|
|
where
|
|
$C_-(\xi)=\int_{-\infty}^0{\rm e}^{-i\bar x\xi}f(\bar x)\, d\bar x$.
|
|
Inverting the Fourier transform we obtain, for all $\bar x<0$,
|
|
\begin{equation}\label{psiminus}
|
|
{\psi}_-(\bar x,\bar t)=h_-(\bar x,\bar t)+ \frac i{{2\pi}} \int_0^{\bar t} d\bar s\, \psi_{\bar x,0}(\bar s) H_0(\bar t-\bar s,\bar x)+ \frac i{{2\pi}} \int_0^{\bar t} d\bar s\, \psi_0(\bar s) \, \partial_{\bar x}H_0(\bar t-\bar s,\bar x)
|
|
\end{equation}
|
|
where
|
|
$$
|
|
h_-(\bar x,\bar t) =\frac 1{{2\pi}}\int_{-\infty}^0 f(\bar y) H_0(\bar t,\bar x-\bar y)\, d\bar y,
|
|
\quad
|
|
H_0(\bar t,\bar r)= \frac {\sqrt{\pi}}{\sqrt{i}{\sqrt{\bar t}}} {\rm e}^{\frac{i\bar r^2}{4\bar t}}
|
|
$$
|
|
\bigskip
|
|
|
|
\indent The equation for $\hat{\psi}_+(\xi,\bar t)$ is given by
|
|
$$i\frac{\partial \hat{\psi}_+}{\partial\bar t}=-i 2\bar E \cos(\bar t)\, \frac{\partial \hat{\psi}_+}{\partial \xi} +\left( \xi^2+\bar U\right) \hat{\psi}_+(\xi,\bar t) + \frac1{\sqrt{2\pi}} \psi_{\bar x,0}(\bar t)+i\xi \frac1{\sqrt{2\pi}} \psi_0(\bar t) $$
|
|
Solving for $\hat\psi_+(\xi,\bar t)$ and inverting the Fourier transforms, we obtain
|
|
\begin{equation}\begin{array}{>\displaystyle c}\label{calcpsip}
|
|
\psi_+(\bar x,\bar t)=h_+(\bar x,\bar t)+\\
|
|
\frac1{{2\pi}} {\rm e}^{2i\bar E \bar x\sin(\bar t)}\int_{-\infty}^{\infty} du\,{\rm e}^{i\bar xu-\Phi(u,\bar t)} \int_0^{\bar t} \, {\rm e}^{\Phi(u,\bar s)}\left[ -i \psi_{\bar x,0}(\bar s)+2\bar E\sin\bar s\, \psi_0(\bar s)\right]\, d\bar s\\
|
|
+ \frac1{{2\pi}} {\rm e}^{2i\bar E \bar x\sin(\bar t)}\int_{-\infty}^{\infty} du\,{\rm e}^{i\bar xu-\Phi(u,\bar t)} \int_0^{\bar t} \, {\rm e}^{\Phi(u,\bar s)}\,u\, \psi_0(\bar s)d\bar s
|
|
\end{array}\end{equation}
|
|
where
|
|
\begin{equation}\begin{array}{>\displaystyle c}
|
|
\Phi(u,\bar t)=i\left[ (u^2 +\bar U+2\bar E^2)\bar t - 4\bar E u \cos(\bar t) -\bar E^2 \sin(2\bar t) \right], \\
|
|
h_+(\bar x,\bar t)=\frac 1{2\sqrt{i \pi }}{\rm e }^{2i\bar x\bar E\sin(\bar t)-i(\bar U+2\bar E^2)\bar t+i\bar E^2\sin2\bar t}\int_0^{\infty }d\bar y\, f(\bar y)\, \frac 1{\sqrt{\bar t}}\,{\rm e}^{\frac{i(\bar x-\bar y-4\bar E+4\bar E\cos(\bar t))^2}{4\bar t}}
|
|
\end{array}\end{equation}
|
|
\bigskip
|
|
|
|
\indent Taking now the limit $\bar x\to0$ from left and right in (\ref{psiminus}) and in (\ref{calcpsip}) and setting $\psi_-(0,\bar t)=\psi_+(0,\bar t)$, and similarly for $\partial_{\bar x}\psi_\pm $ we obtain (in fact $\psi_-(0,\bar t)=\psi_0(\bar t)=\psi_+(0,\bar t)$ implies that $\partial_{\bar x}\psi_-(0,\bar t)=\partial_{\bar x}\psi_+(0,\bar t)$) one equation with the only unknown function $\psi_0(\bar t)$, which has to satisfy the integral equation
|
|
\begin{equation}\begin{array}{>\displaystyle c}\label{eqcontr}
|
|
\psi_0(\bar t)= h_+(0,\bar t)+ h_-(0,\bar t) -\frac{1}{\pi}\int_0^{\bar t}\left( \int_0^{\bar s} \frac{h_-(0,\bar s')}{\sqrt{\bar s-\bar s'}}\, d\bar s'\right) \, G_1(\bar s,\bar t)\,d\bar s\\
|
|
+\frac{1}{2\pi}\int_0^{\bar t} \left( \int_0^{\bar s} \frac{\psi_0(\bar s')}{\sqrt{\bar s-\bar s'}}\, d\bar s'\right) \, G_1(\bar s,\bar t)\, d\bar s \\
|
|
+\frac{\bar E}{\sqrt{i\pi}}\int_0^{\bar t} d\bar s \, \,\psi_0(\bar s)\, \frac 1{\sqrt{\bar t-\bar s}}\, \left( \sin\bar s\, + \frac{ \cos(\bar t)- \cos(\bar s)}{\bar t-\bar s}\right) \,{\rm e}^{iF_0(\bar s,\bar t)}
|
|
\end{array}\end{equation}
|
|
where
|
|
$$G_1(\bar s,\bar t)=\frac d{d\bar s}\frac{e^{iF_0(\bar s,\bar t)}-1}{\sqrt{\bar t-\bar s}}$$
|
|
with
|
|
\begin{equation}\label{Fzero}
|
|
F_0(s,\bar t)=\frac{4\bar E^2 (\cos(\bar t)- \cos(\bar s))^2}{\bar t-\bar s} -(\bar U_2+2\bar E^2)(\bar t-\bar s) +\bar E^2(\sin 2\bar t-\sin 2\bar s)
|
|
\end{equation}
|
|
We then prove that the integral equation (\ref{eqcontr}) has a unique solution, which implies that (\ref{eq12}) does as well.
|
|
\bigskip
|
|
|
|
\indent Setting the initial condition $f(\bar x)=\varphi_1(\bar x)$ (obtained from (\ref{initial}) after the rescaling (\ref{rescale}))
|
|
we obtain explicit functions in equation (\ref{eqcontr}) for $\psi_0(\bar t)$,
|
|
\begin{equation}\label{hminspec}
|
|
h_-(0,\bar t)= {{\rm e}^{-i{\bar k}^{2}\bar t}}\,\frac12\, \left[ 1+{\rm erf} \left( -i^{3/2}\bar k\sqrt {\bar t}\right) \right]
|
|
-R_0\, {{\rm e}^{-i{k}^{2}\bar t}}\,\frac12\, \left[ -1+{\rm erf} \left( -i^{3/2} \bar k\sqrt {\bar t}\right) \right]
|
|
\end{equation}
|
|
and
|
|
\begin{equation}\begin{array}{>\displaystyle c}\label{hpluspec}
|
|
h_+(0,\bar t)= T_0\,\exp[{-i(\bar U+2\bar E^2)\bar t+i\bar E^2\sin 2\bar t-4\,\bar E \alpha\,\cos\bar t +i{\alpha}^{2}\bar t+4\,\bar E \alpha}] \times \\
|
|
\frac 1{2} \left[1- {\rm erf} \left(2\,i^{3/2}\bar E\frac{\cos\bar t-1}{\sqrt {\bar t}}+ \alpha \sqrt{i\bar t} \right) \right]
|
|
\end{array}\end{equation}
|
|
where
|
|
$$\bar k=k/\sqrt{2\omega},\ \ \ \alpha=\sqrt{\bar U-\bar k^2},\ \ \ 1+R_0=T_0=\frac{2i\bar k}{i\bar k-\sqrt{\bar U-\bar k^2}}.$$
|
|
|
|
\section{Short time behavior: approach to the stationary state}\label{Sec3}
|
|
\indent We carried out numerical solution of (\ref{eqcontr}) using (\ref{hminspec}) and (\ref{hpluspec}) for a variety of values of $E$.
|
|
To compare with the work in \cite{YGR11}, we take $\frac {\hbar^2k^2}{2m}=\mathcal E_F=4.5$\,eV, $W=5.5$\,eV.
|
|
Unless otherwise specified, we also take $\omega=1.55\ \mathrm{eV}$.
|
|
The laser period $\tau=\frac{2\pi}\omega$ is then equal to $2.7\ \mathrm{fs}$.
|
|
\bigskip
|
|
|
|
\indent The details of the numerical algorithms, which take into account the singular behavior of $\psi(x,t)$ at $x=0$ and $t=0$, will be given in a separate publication.
|
|
\bigskip
|
|
|
|
\indent Figure \ref{fig:wavefunction} shows the real and imaginary parts of $e^{i\frac{k^2}2t}\psi_0(t)$, as well as $|\psi_0(t)|^2$ and the normalized current $\frac1kj(0,t)$.
|
|
While $|\psi_0|^2$ seems more in phase with the real part of$e^{i\frac{k^2}2t}\psi_0$, the current looks more in phase with the imaginary part.
|
|
Figure \ref{fig:current} shows the values of the current $j(0,t)$ passing through the origin as a function of time for different strengths of the field, all at $\omega=1.55$\, eV.
|
|
We see there a change of behavior as $E$ increases from $E=1$\,V/nm to $E=30$\,V/nm (the Keldysh parameter $\gamma=\frac{2\omega}E\sqrt W$ goes from 18.6 to 0.62).
|
|
For large values of $E$, fast oscillations appear, which become faster and larger as $E$ grows.
|
|
The oscillatory behavior within one period for strong fields is also observed in the approximate solution \cite{YGR11}, though the details vary.
|
|
We have no simple physical explanation for these fast oscillations.
|
|
They do not fit in with the ``simple three step'' scenario \cite{KSK92,Co93,KLe18}.
|
|
We note that the height of the first maximum is linear in $E$, while its location is almost independent of $E$.
|
|
\bigskip
|
|
|
|
\begin{figure}
|
|
\includegraphics[width=8cm]{real.pdf}
|
|
\hfill
|
|
\includegraphics[width=8cm]{imag.pdf}
|
|
\par\penalty10000\bigskip\penalty10000
|
|
\includegraphics[width=8cm]{density.pdf}
|
|
\hfill
|
|
\includegraphics[width=8cm]{current.pdf}
|
|
\caption{
|
|
The real (a) and imaginary (b) parts of $e^{i\frac{k^2}2t}\psi_0(t)$ the density (c) and the current (d) as a function of $\frac{t\omega}{2\pi}$ for 3 periods, with $E=15\ \mathrm V\cdot\mathrm{nm}^{-1}$, $\omega=1.55\ \mathrm{eV}$.
|
|
The dotted line represents the phase of the field: $\cos(\omega t)$.
|
|
}
|
|
\label{fig:wavefunction}
|
|
\end{figure}
|
|
|
|
\begin{figure}
|
|
\includegraphics[width=8cm]{current01.pdf}
|
|
\hfill
|
|
\includegraphics[width=8cm]{current05.pdf}
|
|
\par\penalty10000\bigskip\penalty10000
|
|
\includegraphics[width=8cm]{current18.pdf}
|
|
\hfill
|
|
\includegraphics[width=8cm]{current30.pdf}
|
|
\caption{
|
|
The normalized current $\frac jk$ as a function of $\frac{t\omega}{2\pi}$ for $\omega=1.55\ \mathrm{eV}$ and various values of the electric field: $E=1,\ 5,\ 18.6,\ 30\ \mathrm{V}\cdot\mathrm{nm}^{-1}$.
|
|
The Keldysh parameter $\gamma=\frac{2\omega}E\sqrt W$ for these fields is, respectively, $18.6$, $3.73$, $1.00$ and $0.621$.
|
|
The dotted line represents the phase of the electric field: $\cos(\omega t)$.
|
|
As the field increases, fast oscillations in the current appear.
|
|
These fast oscillations mostly occur while the field is increasing.
|
|
}
|
|
\label{fig:current}
|
|
\end{figure}
|
|
|
|
\indent In Figure \ref{fig:cn}, we show a comparison of our solution with the results obtained from a direct solution of (\ref{schrodinger}) via the Crank-Nicolson algorithm.
|
|
The agreement, especially for the location of the maxima and minima after very early times is very good.
|
|
At short times, however, the agreement is not so good.
|
|
This is expected since the Crank-Nicolson scheme is based on approximating derivatives by finite differences.
|
|
However, at short times, $\psi(0,t)\sim t^{\frac32}$, which has a singular second derivative at $0$, so $\partial_t\psi$ is poorly approximated at short times by finite differences.
|
|
Note, also, that the Crank-Nicolson method requires a cut-off in $x$, restricting $\psi(x,t)$ to $x\in[-a,a]$.
|
|
This causes distortions in $\psi$ due to reflections from these artificial boundaries, and so can only be used for short times (this can be avoided by using non-local boundary conditions, such as ``transparent boundary conditions'').
|
|
\bigskip
|
|
|
|
\begin{figure}
|
|
\includegraphics[width=8cm]{cn.pdf}
|
|
\hfill
|
|
\includegraphics[width=8cm]{cn_short.pdf}
|
|
\caption{
|
|
The current computed with our method (blue, color online), compared with the Crank-Nicolson algorithm (red), for $\omega=1.55\ \mathrm{eV}$ and $E=15\ \mathrm{V}\cdot\mathrm{nm}^{-1}$.
|
|
The maxima and minima seem to occur at the same time, and the agreement is somewhat good for $t>\frac\pi\omega$.
|
|
The graph on the right focuses on short times, for which the Crank-Nicolson algorithm produces a different, and quite surprising result: the current initially shoots down to negative values before rising back up, but the field is initially positive, which should drive the current up, not down.
|
|
The discrepancy is due to the fact that the wavefunction is singular at $t=0$ (it has a divergent second derivative), and that the Crank-Nicolson algorithm approximates derivatives using finite differences, which is good only for regular functions.
|
|
}
|
|
\label{fig:cn}
|
|
\end{figure}
|
|
|
|
\indent Also shown in the figures is the running average of the current
|
|
\begin{equation}
|
|
\left<j\right>_t:=\frac1\tau\int_{t-\tau}^t j(0,t)
|
|
,\quad
|
|
\tau:=\frac{2\pi}\omega
|
|
.
|
|
\end{equation}
|
|
It is seen there that while the current appears to approach its long time value rather rapidly it does not do so fully until much later.
|
|
As seen in Figure \ref{fig:avgcurrent}, where we plot $\left<j\right>$ for $\omega=6\ \mathrm{eV}$, $E=10\ \mathrm{V}\cdot\mathrm{nm}^{-1}$, $\gamma=9.6$, $\left<j\right>_t$ is still not all that close to its asymptotic value even when $t\approx48\tau$.
|
|
\bigskip
|
|
|
|
\indent The rate of the decay of the average current to its asymptotic value is evaluated in Figure \ref{fig:decay}.
|
|
In order to compute this rate without having to guess the asymptotic value, we proceed as follows.
|
|
At the end of every laser period $t_n=\frac{2\pi}\omega n$, we compute the minimal value $\mu_n$ and maximal value $M_n$ of the average current in the period $(\frac{2\pi}\omega(n-1),\frac{2\pi}\omega n]$.
|
|
The plot shoes $M_n-\mu_n$ as a function of $t$, and shows that $\left<j\right>_t\approx\left<j\right>_\infty+g(t)t^{-\frac32}$, where $g(t)$ is bounded and asymptotically constant.
|
|
\bigskip
|
|
|
|
\begin{figure}
|
|
\includegraphics[width=8cm]{avgcurrent.pdf}
|
|
\hfill
|
|
\includegraphics[width=8cm]{avgcurrent_zoom.pdf}
|
|
\caption{
|
|
The normalized average current $\frac1k\left<j\right>_t$ as a function if $\frac{t\omega}{2\pi}$ for $\omega=6\ \mathrm{eV}$ and $E=10\ \mathrm{V}\cdot\mathrm{nm}^{-1}$.
|
|
The two plots show the same data, but the first one focuses on times between the 12th and 48th period.
|
|
Here, $\omega$ is large enough that absorbing one photon suffices to overcome the work function.
|
|
Even after 48 periods, the average current has not converged to its asymptotic value.
|
|
}
|
|
\label{fig:avgcurrent}
|
|
\end{figure}
|
|
|
|
\begin{figure}
|
|
\hfil\includegraphics[width=8cm]{decay.pdf}
|
|
\caption{
|
|
The convergence rate to the asymptote of the average current as a function of $\frac{t\omega}{2\pi}$ on a log-log plot for $\omega=6$ and $E=10\ \mathrm{V}\cdot\mathrm{nm}^{-1}$.
|
|
The dots are computed at the end of each period $t_n=\frac{2\pi}\omega n$, and their value is the difference between the maximun $M_t$ and the minimum $\mu_t$ of the normalized average current $\frac1k\left<j\right>_t$ in the period immediately preceding $t_n$.
|
|
The red line is a plot of $0.0030\times(\frac{t\omega}{2\pi})^{-\frac32}$, which fits the data rather well.
|
|
}
|
|
\label{fig:decay}
|
|
\end{figure}
|
|
|
|
\indent In Figure \ref{fig:omega}, we show an estimate of the asymptotic average current as a function of $\omega$ in the vicinity of the one-photon threshold $\omega_c=W$ for $E=10\ \mathrm{V}\cdot\mathrm{nm}^{-1}$.
|
|
As we argued above, the average current converges slowly to its asymptotic value.
|
|
In order to estimate this asymptotic value without computing the current at very large times, we actually compute the average of the average current over a laser period:
|
|
\begin{equation}
|
|
\left<\left<j\right>\right>:=\frac1\tau\int_{T-\tau}^T dt\ \left<j\right>_t
|
|
.
|
|
\end{equation}
|
|
We see that there is a steep increase in $\left<\left<j\right>\right>$ as $\omega$ increases past $\omega_c$.
|
|
This is precisely what is observed in experiments on the photoelectric effect, where the emission of electrons from the metal surface has such a threshold \cite{Mi16}.
|
|
This became a key element in Einstein's ansatz of localized photons.
|
|
Here, in the semi-classical treatment of the laser field, this phenomenon appears as a resonance.
|
|
The ponderomotive term $\frac{E^2}{4\omega^2}$, which appears in the photon energy, see below, is negligible compared to the work function $W$.
|
|
\bigskip
|
|
|
|
\begin{figure}
|
|
\hfil\includegraphics[width=8cm]{omega.pdf}
|
|
\caption{
|
|
The average of the average of the normalized current $\frac1k\left<\left<j\right>_t\right>$ after 12 periods as a function of $\omega$, for $E=10\ \mathrm{V}\cdot\mathrm{nm}^{-1}$.
|
|
After 12 periods, the average current $\left<j\right>_t$ is still fluctuating quite a lot, so, in order to eliminate these fluctuations and get closer to the asymptotic value $\left<j\right>_\infty$, we average the average current over a laser period.
|
|
}
|
|
\label{fig:omega}
|
|
\end{figure}
|
|
|
|
\section{Long time behavior of $\psi(x,t)$}\label{longtime}
|
|
\indent The long time behavior of the system is given by the poles of the Laplace transform
|
|
$$\hat{\psi}(x,p)=\int_0^{\infty}{\rm e}^{-pt} \psi(x,t)\, dt$$
|
|
on the imaginary axis, as can be seen by taking the inverse Laplace transform
|
|
\begin{equation}
|
|
\psi(x,t)=\int_{-i\infty}^{i\infty}dp\ e^{pt}\hat\psi(x,p)
|
|
.
|
|
\end{equation}
|
|
Other poles, which necessarily have negative real parts, give rise to exponentially decaying terms while branch cuts generally contribute $t^{-\frac32}$ terms to the approach to the asymptotic state.
|
|
As mentioned earlier, the time asymptotic of $\psi$ is of the form
|
|
\begin{equation}\label{4p1}
|
|
\psi(x,t)\sim {\rm e}^{\frac 12 i k^2 t} \bar{\psi}(x,t)
|
|
\end{equation}
|
|
where $\bar{\psi}(x,t)$ is periodic in $t$ with period $2\pi/\omega$ and $\frac 12 k^2$ is the energy of the electrons in the incoming beam.
|
|
This corresponds to the poles of $\hat{\psi}(p,t)$ on the imaginary axis being at
|
|
\begin{equation}\label{4p2}
|
|
p=-\frac 12 i k^2+i\omega n,\ \ \ n\mathrm{\ \ integer}
|
|
\end{equation}
|
|
\bigskip
|
|
|
|
\indent We further show that $\bar{\psi}(x,t)$ coincides, after the change of gauge, with the wave function $\mathcal{U}(x,t)$ computed by Faisal et al \cite{FKS05}.
|
|
In that work it was assumed (\ref{schrodinger}) (in the magnetic gauge) has a solution of the form (\ref{4p1}) with an incoming beam $e^{ikx}$ for $x<0$, without considering any initial conditions.
|
|
Computing the residues at the poles on the imaginary axis we show that
|
|
\begin{equation}\label{4p3}
|
|
\bar\psi(x,t)\sim \left\{\begin{array}{>\displaystyle ll}
|
|
e^{ikx}
|
|
+\sum_{m\in\mathbb Z}e^{-im\omega t}e^{-ix\sqrt{k^2+2m\omega}}\mathcal R_m
|
|
&\mathrm{for\ }x<0
|
|
\\[0.5cm]
|
|
e^{i\frac E\omega x\sin\omega t}\sum_{n,m\in\mathbb Z}e^{-in\omega t}g_{n-m}^{(\kappa_m)}e^{-\kappa_mx}\mathcal T_m
|
|
&\mathrm{for\ }x>0
|
|
\end{array}\right.
|
|
\end{equation}
|
|
where
|
|
\begin{equation}
|
|
\kappa_m=\sqrt{2U-k^2+\frac{E^2}{2\omega^2}-2m\omega}
|
|
\label{kappa}
|
|
\end{equation}
|
|
and
|
|
\begin{equation}
|
|
g_{n-m}^{(\kappa_m)}=\frac\omega{2\pi}\int_0^{\frac{2\pi}\omega}dt\ e^{-i(n-m)\omega t}e^{\frac{i\frac{E^2}{4\omega^2}}\omega\sin(2\omega t)+\kappa_m\frac{2E}{\omega^2}\cos(\omega t)}
|
|
.
|
|
\end{equation}
|
|
This is exactly of the form considered in \cite{FKS05}.
|
|
$\mathcal R_n$ and $\mathcal T_m$ are computed by matching boundary values of $\psi(0,t)$ and $\psi_{x,0}(t)$ at $x=0$.
|
|
\bigskip
|
|
|
|
\indent A physical interpretation of (\ref{4p3}), see \cite{FKS05}, is that an electron in a beam coming from $-\infty$ and moving in the positive $x$-direction, ${\rm e}^{ikx},\ k>0$, absorbs or emits ``$m$ photons'' and is either reflected, transmitted or damped.
|
|
Transmission occurs when $m\omega>U+\frac {E^2}{4\omega^2}-\frac{k^2}2\equiv\omega_c$.
|
|
The $\frac{E^2}{4\omega^2}$ in (\ref{kappa}) corresponds to the ponderomotive energy of the electron in the laser field.
|
|
Damping occurs when the inequalities are in the opposite direction.
|
|
$\omega_c$ is the minimum frequency necessary to push the electron with kinetic energy $\frac12k^2$ (in the $x$-direction) over the potential barrier of height $U-\frac12k^2=W$.
|
|
For $x\gg0$, the current in the $m$-photon channel will have kinetic energy $m\omega-\omega_c$ and the current will be given by $\sqrt{m\omega-\omega_c}$.
|
|
This will also be equal to the average current at large $t$, which is independent of $x$.
|
|
This explains the picture in Figure 4 for $m=1$.
|
|
The larger $m$ values necessary for smaller $\omega$ are difficult to see.
|
|
\bigskip
|
|
|
|
\indent The asymptotic form (\ref{4p3}) is true for all initial conditions of the form $\Theta(-x)e^{ikx}+f_0(s)$ as long as $f_0(x)$ only contains terms which are square integrable.
|
|
The additional terms in $\psi(x,t)$ which come from $f_0(x)$ go to zero as $t\to\infty$.
|
|
This follows from the fact, proven in \cite{CCe}, that the Floquet operator associated to (\ref{schrodinger}) has no point spectrum.
|
|
\bigskip
|
|
|
|
|
|
|
|
\vfill
|
|
|
|
\delimtitle{\bf Acknowledgements}
|
|
We thank David Huse, Kevin Jensen and Donald Schiffler for very valuable discussions. This work was supported by AFOSR Grant No. FA9500-16-1-0037. O.C. was partially supported by the NSF-DMS (Grant No. 1515755). I.J. was partially supported by the NSF-DMS (Grants No. 31128155 and 1802170). J.L.L. thanks the Institute for Advanced Study for its hospitality.
|
|
\enddelim
|
|
|
|
\eject
|
|
|
|
\begin{thebibliography}{WWW99}
|
|
\small
|
|
\IfFileExists{bibliography/bibliography.tex}{\input bibliography/bibliography.tex}{}
|
|
\end{thebibliography}
|
|
|
|
|
|
\end{document}
|
|
|