L\'evy Distribution of Single Molecule Line Shape Cumulants in Low Temperature Glass
E. Barkai$^a$
R. Silbey$^a$ and G. Zumofen$^b$
$^a$ Department of Chemistry and Center for Material Science and Engineering
Cambridge
MA 02139.
$^b$ Physical Chemistry Laboratory
ETH-Z
CH-8092 Z$\ddot{u}$rich.
{Abstract} We investigate the distribution of single molecule line shape cumulants, $\kappa_1,\kappa_2,\cdots$, in low temperature glasses based on the sudden jump, standard tunneling model. We find that the cumulants are described by L\'evy stable laws, thus generalized central limit theorem is applicable for this problem.
PRL.tex
% use dvips -o name.ps name.dvi
%
%\documentstyle[preprint,eqsecnum,aps]{rrevtex}
%              ^^^^^^^^ this makes bigger letters and double-spacing
%\documentstyle[prl,aps,twocolumn]{revtex}
\documentstyle[prl,aps,epsf,twocolumn]{revtex}
%\documentstyle[12pt]{article}
\begin{document}
\title{
L\'evy Distribution of Single Molecule
Line Shape Cumulants in Low Temperature Glass
 }
\author{E. Barkai$^a$, R. Silbey$^a$ and G. Zumofen$^b$\\
$^a$ Department of Chemistry and Center for Material Science and Engineering,
Massachusetts Institute of Technology, Cambridge, MA 02139.\\
$^b$ Physical Chemistry Laboratory, Swiss Federal Institute of Technology, ETH-Z, CH-8092 Z$\ddot{u}$rich. \\
}
\date{\today}
\maketitle

\begin{abstract}

{\bf Abstract} 

We investigate the distribution of single molecule line shape
 cumulants, $\kappa_1,\kappa_2,\cdots$,
in low temperature glasses
based on the sudden jump, standard tunneling model. 
We find that the cumulants are described 
by L\'evy stable laws, thus generalized central limit theorem
is applicable for this problem.

\end{abstract}

Pacs: 05.40.-a, 05.40.Fb, 61.43.FS, 78.66.Jg \\

% 61.43.FS Glass

%\Huge{Draft Not For Distribution}

% Modern experimental techniques have made it possible to measure
%line shapes of single molecule (SM) embedded in solids \cite{MO}. 
%Line shapes of chemically  identical SMs in glassy materials 
%vary from one molecule to the other since each
%molecule is in a unique environment \cite{Tittel}.
%By investigating
%properties of line shapes of SMs in disordered condensed phase 
% one can learn
%about the dynamical and disorder 
%properties of the host medium 
%\cite{MO,Tittel,Reilly2,Zumofen,Geva,Pfluegl,Orrit,Barkai1,Barkai2}.  
% Since each line shape is a random function one must consider
%the statistical properties of the lines. 
%Geva and Skinner \cite{Geva} have modeled
%single molecule line shapes in low temperature glass
% based on the sudden jump approach of Anderson and Kubo 
%\cite{AK,Reilly2}.
%The model  assumes that low density dynamical defects 
%[e.g. spins,  two level systems (TLS)] are 
%distributed randomly and uniformally in the medium. These defects
%interact with the molecule via {\em long range interactions} (e.g., the dipole
%interaction) while the defect-defect interaction is neglected.
%
%The line of each SM  can
%be characterized by its cumulants $\kappa_j$, $(j=1,2,\cdots)$
%which will vary from one molecule to the other. In this letter we
%show
%that
%the probability densities $P(\kappa_1)$, $P(\kappa_2)$ are
%symmetrical and one sided L\'evy stable
%laws respectively, while
%higher order cumulants  
%are described 
%by
%L\'evy statistics only in the {\em slow modulation limit}.

 Recent experimental advances \cite{MO} have made it possible to measure
the spectral line shape of a single molecule 
(SM) embedded in a condensed phase. 
Because each molecule is in a unique static and dynamic environment,
the line shapes of chemically identical SMs vary from molecule 
to molecule \cite{Tittel}.
In this way, the dynamic properties of the host are encoded
in the distribution of single molecule spectral line shapes 
\cite{MO,Tittel,Reilly2,Zumofen,Geva,Pfluegl,Orrit,Barkai1,Barkai2}.
We examine the statistical properties of the line shapes
and show how these are related to the underlying microscopic  
dynamical events occurring in the condensed phase.

 We use the
Geva--Skinner \cite{Geva}
 model for the SM line shape in a low temperature
glass 
based on the sudden jump picture of 
Kubo and Anderson \cite{AK,Reilly}.
In this model, a random distribution of low--density (and
non--interacting) dynamical defects 
[e.g., spins or two level systems (TLS)] interacts with
the molecule via {\em long range 
interaction} 
(e.g., dipolar). 
We show that L\'evy statistics fully characterizes the
properties of the SM spectral line  both in the {\em fast}  and {\em slow
modulation limits}, while far from these limits L\'evy statistics
describes the mean and variance of the line shape. 
We then compare our
analytical results, derived in the slow modulation limit, with results 
obtained from numerical simulation.
The good agreement indicates that the slow modulation limit is correct
for the parameter set relevant to experiment.
%Specifically, the spectral line is characterized by its cumulants, 
%Below, we shall show that L\'evy statistics fully characterizes the
%properties of the SM spectral line. 
%Specifically, the spectral line is characterized by its cumulants, 
%$\kappa_j$
%$(j=1,2,\cdots)$ that vary from molecule to molecule.
%We show that the probability densities $P(\kappa_1)$
%and $P(\kappa_2)$ are symmetrical and one sided L\'evy stable laws,
%respectively, while higher order cumulants are described by L\'evy
%statistics only in the {\em slow modulation limit}.
%In the opposite limit of fast modulation we show that the distribution
%of line widths is a one sided L\'evy stable law.
%We then  compare our
%analytical results, derived in the slow modulation limit, with results 
%obtained from numerical simulation.
%The good agreement indicates that the slow modulation limit is correct,
%for the parameter set relevant to experiment.

%fast modulation, spectral lines are approximately
%Lorentzian, then the distribution of line widths can be used
%to characterize the statistical properties of the lines.
%For this case the distribution of line width is a one sided L\'evy
%stable law. 

 L\'evy stable distributions 
serve as
a natural generalization of the normal Gaussian distribution.
The importance of the Gaussian in statistical physics
stems from the central limit theorem.  L\'evy stable laws are used
when analyzing sums of the type $\sum x_i$, with 
$\{ x_i \}$ being independent identically distributed
 random variables
characterized by a diverging  variance.
In this case the ordinary Gaussian central limit theorem must be replaced
with the generalized central limit theorem.
With this generalization, 
L\'evy stable probability densities, $L_{\gamma,\eta}(x)$, replace
the  Gaussian of the standard central
limit theorem. 
Khintchine and L\'evy found
that stable characteristic functions,
$\hat{L}_{\gamma,\eta}(k)$,
 are of the form \cite{Feller}
%
\begin{equation}
\ln\left[ \hat{L}_{\gamma,\eta} \left( k \right) \right] = i \mu k 
- z_{\gamma}|k|^{\gamma} 
\left[ 1 - i \eta { k \over |k|} 
\tan\left( { \pi \gamma \over 2}\right)\right]
\label{eq1}
\end{equation}
%
for $0< \gamma\le 2$ (for the case $\eta \ne 0$, $\gamma=1$ see
\cite{Feller}). Four parameters are needed for a full description of 
a stable law. The constant $\gamma$ is called the characteristic exponent, 
the parameter $\mu$ is a location parameter which is unimportant in
the present case, $z_{\gamma} > 0$ is a scale parameter and $-1\le \eta \le 1$
is the index of symmetry. When $\eta=0$ $(\eta=\pm 1)$
the stable density $L_{\gamma,\eta}(x)$ is symmetrical (one sided).
%Today the term L\'evy behavior is often used in a wider sense than that
%implied by the functions in Eq. (\ref{eq1}) and is used to describe behavior
%with broad tails. Examples for  recent applications of 
%L\'evy statistics are sub recoil 
%laser cooling \cite{Cohen} and the characterization of geometrical properties
%of random fractals \cite{Hovi}. 
L\'evy statistics is known  to describe
several long range interaction systems in 
diverse fields such as  astronomy \cite{Feller}, turbulence \cite{Turbo}
and spin glass  \cite{Jap}. 
%All these phenomena result from long-range interactions with
%uniformally distributed defects: gravitational masses, vortexes,
%magnetic dipoles etc.
Stoneham's theory \cite{Stoneham}
of inhomogeneous line broadening in defected crystal, is based on long-range forces
and parts of it can be
interpreted in terms of L\'evy stable laws \cite{Barkai2}. 
Stoneham's approach \cite{Stoneham} is inherently static, while the SM
line shape model  
considers both dynamical and static contributions from
the defects.

% In \cite{Geva} the distribution of line widths was investigated
%based on numerical approach and it was shown that this model
%can reasonably reproduce single molecule experimental results,
%as well as hole burning and photon echo measurements 
%in low temperature glass.
%Based on simulations, we show that the L\'evy
%results are reasonable for this model and that indeed a slow modulation
%limit exists for a relevant choice of parameters.
%Our finding follows previous
%numerical work of 
%Pfluegl et al \cite{Pfluegl}  who
%showed  that the 
%distribution of line widths 
%is well approximated by the slow modulation limit.
%We note that line widths  and cumulants are not related
%by any known  formula for the case considered here
%\cite{remark1}.
%Here we use L\'evy statistics as a test for the validity
%of the slow modulation limit in low temperature glass.
%

An important issue is the slow and fast modulation limits
\cite{Pfluegl,Reilly}. 
Briefly, the  fast (slow) modulation limit is valid if important
contributions to the line shape are from TLSs which satisfy
$\nu \ll K$ $(\nu \gg K)$, where 
$\nu$ is the frequency shift of the SM due to 
SM--TLS
interaction
and $K$ is the 
transition rate of the flipping TLS
(see details below).
In the fast modulation limit,
all (or most) lines are Lorentzian 
with a width 
that varies from one molecule to the other. For this case, the (L\'evy)
distribution
of {\em line widths} fully characterizes the statistical properties
of the lines. 
%Below we show that in this
%limit the distribution of line widths is a one sided
%L\'evy stable law. 
The second, more complicated, case corresponds
to the slow modulation
limit. Then the SM line is typically 
composed of several peaks
(splitting) and is not described well by a Lorentzian. 
If a SM shows splitting, one
can investigate the validity of the
standard tunneling model of glass \cite{AHV} in a direct way, 
since the splitting of
a line is directly associated with SM--TLS interaction \cite{Orrit}. 
%From an experimental point of view it is difficult, though not 
%impossible to determine whether two distinct peaks
%arise from two molecules
%or a split single molecule line. 
As mentioned, we demonstrate the existence
of a slow modulation limit in SM--glass system.

%Below, we compare our
%analytical results, derived in the slow modulation limit, with results 
%obtained from numerical simulation.
%The good agreement indicates that a slow modulation limit exists
%for the set of parameters relevant to experiment.
 
 Following \cite{Geva}
we assume a SM coupled to non identical independent TLSs
at distances ${\bf r}$ in dimension $d$.
Each TLS is characterized by its asymmetry variable $A$
and tunneling element $J$.
The energy of the TLS is $E=\sqrt{A^2 + J^2}$.
The TLSs are coupled to phonons
or other thermal excitations
such that the state of the TLS changes with time. The state of
the $n$ th TLS is described by an occupation parameter,
$\xi_n(t)$, that is equal $0$ or $1$ if the TLS is in its
ground or excited state respectively. The probability for finding the TLS
in its upper $\xi=1$ state, 
$p$, is given 
by the standard Boltzmann form
$p=1/\{ 1+ \exp[ E/(k_b T)]\}$. 
The transitions
between the ground and excited state are described by the
up and down transition rates $K_u,K_d$, which  are related to
each other by the standard detailed balance condition.

 The excitation of the $n$ th TLS shifts the SM's
transition frequency by $\nu_n$. Thus, the SM's
transition frequency is
%
\begin{equation}
\omega(t) = \omega_0 + \sum_{n=1}^{N_{\rm act}} \xi_n \left( t \right) \nu_n,
\label{eq2}
\end{equation}
%
where $N_{\rm act}$ is the number of active TLSs in the system
(see details below) and 
$\omega_0$ is the bare transition frequency that differs
from one molecule to the other depending on the
local static disorder.
We consider a wide class of frequency perturbations 
%
\begin{equation}
\nu = 2 \pi \alpha \Psi\left( \Omega\right) f\left( A, J \right)
{ 1 \over r^{\delta}},
\label{eq3}
\end{equation}
%
$\alpha$ is a coupling constant with units
$\left[ \mbox{Hz}\ \mbox{nm}^{\delta}\right]$,
$\Psi\left(\Omega \right)$ is a dimensionless
function of order unity, $\Omega$ is a vector of angles
determined by the orientations of the TLS and molecule
(in some simple cases $\Omega$ depends on polar angles
only),
$f(A,J)\ge 0$  is a dimensionless function of the internal
degrees of freedom of the fluctuating TLS,  
$\delta$
is the interaction exponent. 
 The line shape of the SM is given by
the complex Laplace transform of the 
relaxation function
%
\begin{equation}
I_{SM}\left( \omega \right)={1 \over \pi} \mbox{Re}
\left[  \int_0^\infty dt e^{i \omega t} 
%e^{ - {t \over 2  T_1}}
\Pi_{n=1}^{N_{\rm act}}\Phi_n\left( t \right) \right]
\label{eq4}
\end{equation}
%
provided that the natural life time of the SM excited state is long.
The relaxation function  of a single TLS was evaluated \cite{Reilly} based
on methods developed in \cite{AK}
%
\begin{equation}
\Phi\left( t \right) = e^{-\left( \Xi+ i p\nu\right) t}
\left[\cosh\left( \Omega t \right) + {\Xi\over \Omega}
\sinh\left( \Omega t \right) \right]
\label{eq5}
\end{equation}
%
with
%
%\begin{equation}
$\Omega = [K^2/ 4 - \nu^2 / 4 - i \left( p - 1/2 
 \right) \nu K]^{1/2}$,
%\label{eq6}
%end{equation}
%
%\begin{equation}
$\Xi = { K \over 2} - i \left( p - {1 \over 2} \right) \nu$
%\label{eq7}
%\end{equation}
%
and $K=K_u + K_d$.
For a bath of TLSs the line shape, Eq. (\ref{eq4}), 
is a formidable function of the random TLS
parameters $(r,\Omega,A,J)$ as well as the system
parameters $( \alpha , T, \mbox{etc})$.
In the fast modulation limit
$K\gg |\nu|$, one finds a simpler behavior: all lines are Lorentzian
with half width
%
\begin{equation}
\tilde{\Gamma}= \sum_n^{N_{act}} p_n(1-p_n)\nu_n^2/K_n
\label{eqGamma}
\end{equation}
%
which varies from one molecule to the other.
Eq. (\ref{eqGamma}) shows the well known
phenomena of motional narrowing.
In the slow modulation limit $K \ll |\nu|$
one finds  
$ \Phi(t)=1 - p + pe^{- i \nu t}$ implying that the line shape
of a molecule coupled to a single TLS 
is composed of two delta peaks, the line shape 
of a molecule coupled to
two TLSs is composed of four delta peaks, etc (splitting).

 The spectral line is characterized by its cumulants
$\kappa_j$ $(j=1,2,\cdots)$ that vary from one molecule to the
other, and we investigate the cumulant probability density
$P(\kappa_j)$.
 We have derived the cumulants of the SM line shape, 
and the first four cumulants are
presented in Table 1 \cite{remark}. 
%The exact calculation of the cumulants was 
%made possible using symbolic programming without which the
%calculation becomes cumbersome.
We observe that cumulants of order 
$j\le 2$ are real while generally cumulants of
order $j>2$ are complex,
implying that  the moments of the line shape
diverge when $j>2$. The summation, $\sum_n$, in Table $1$
is over the active TLSs, namely those TLSs which flip
on the time scale of observation $\tau$ (i.e., $K_n > 1/ \tau$).
 We consider the slow modulation limit, soon to be justified,
which means that we consider the case $K_n \ll \nu_n$.
To investigate this limit we set $K_n=0$ in Table $1$,
then all the cumulants
are real and are rewritten as $\kappa_j=\sum_n H_{j n} \nu_n^j$,
where $H_j$ are functions of $p$ only
and $H_1=p$, $H_2=p(1-p)$, $H_3=p(1-p)(2p-1)$ etc.
Note that for $\kappa_1$ and $\kappa_2$
no approximation 
is made.


\begin{tabular}{|l|l|}
\hline
%
j & $\kappa_j$  \\
%
$ $ & $$ \\
%
$1$ &$ \sum_n p_n \nu_n $  \\
%
$ $ & $$ \\
%
$2$ &$  \sum_n p_n\left(1 - p_n\right)  \nu_n^2  $   \\
%
$ $ & $$ \\
%
$3$ &$ \sum_n p_n\left(1- p_n\right)\left( 2 p_n - 1\right) \nu_n^3
+i \sum_n p_n(1-p_n)K_n \nu_n^2 $   \\
%
$ $ & $$ \\
%
$4$ &$ \sum_n p_n\left(-1+ p_n\right)
\left[ K_n^2 + \nu_n^2\left(-1 + 6p_n - 6 p_n^2\right) \right]\nu_n^2 -$\\
$ \ $ &$
2 i \sum_n  K_n \left( - 1 + p_n \right) p_n \left( - 1 + 2 p_n\right) \nu_n^3
$   \\
%
$ $ & $$ \\
%
\hline
\end{tabular}
$$\mbox{Table  1: Cumulants $\kappa_j$ of the SM line shape} $$
%

Let 
$\left\langle \cdot \right\rangle_{r \Omega A J}$  denote an averaging
over
the random TLS parameters.
The characteristic function
of the $j$ cumulant can be written in a form
%
$$ \left\langle \exp \left( i k \kappa_j\right)\right \rangle_{r \Omega A J} = $$
\begin{equation}
\exp\left[ - \rho_{\mbox{eff}} \left\langle \int d \Omega  \int_0^{\infty}
{d (r^d) \over d}\left( 1 - e^{i k B_j/ r^{\delta j }}\right)\right\rangle_{AJ}\right],
\label{eq8}
\end{equation} 
%
$ \rho_{\mbox{eff}}$ is the density of the active TLS and 
$B_j=(2 \pi \alpha)^j \Psi^j \left( \Omega \right) f^j \left(A,J\right)H_j$.
To derive Eq. (\ref{eq8}) we have
used the assumption of independent
TLSs uniformly distributed in the system.
 For odd $j$  cumulants we find
%
\begin{equation}
\left\langle \exp\left(i k \kappa_j \right) \right\rangle_{r \Omega A J}=
\hat{L}_{\gamma,0}( k ),
\label{eq9}
\end{equation}
%
with characteristic exponent $\gamma=d/(\delta j )$ and the scale parameter
%
\begin{equation}
z_{\gamma}= \rho_{\mbox{eff}} \left( 2 \pi \alpha\right)^{{d\over \delta}}
\left\langle f^{d/\delta}\left( A , J \right) |H_j|^{\gamma} \right\rangle_{AJ}
c_{\gamma} \int d \Omega |\Psi^j\left( \Omega\right)|^{\gamma}
\label{eq10}
\end{equation}
%
with $c_{\gamma} = \cos\left( \gamma \pi/2\right)\Gamma(1 - \gamma)$,
$c_1=\pi/2$.
Eq.  (\ref{eq9}) shows that
odd cumulants are described by symmetrical L\'evy stable density, i.e.,
$P(\kappa_j)=L_{\gamma,0}(\kappa_j)$. 
Two conditions must be satisfied for such
a behavior,  $0 <\gamma < 2$ and 
$\int d \Omega \sin\left[ \Psi^j(\Omega)\right]=0$.  The latter condition
gives the symmetry condition, $\eta=0$,  which means that negative
and positive contributions to $\kappa_j$ are equally probable.

 For even cumulants and $0< \gamma < 1$ we find
%
\begin{equation}
\left\langle \exp\left(i k \kappa_j \right) \right\rangle_{r \Omega A J}=
\hat{L}_{\gamma,\eta}(k)
\label{eq11}
\end{equation}
%
with a scale parameter Eq. (\ref{eq10})
and with L\'evy index of symmetry 
%
\begin{equation}
\eta= { \left\langle f^{j \gamma} \left( A, J \right) |H_j|^{\gamma} { H_j \over |H_j|} \right\rangle_{AJ} \over
\left\langle f^{j \gamma} \left( A, J \right) |H_j|^{\gamma} \right\rangle_{AJ}}.
\label{eq12}
\end{equation}
%
Eq. (\ref{eq11}) implies that even cumulants are distributed
according to $P(\kappa_j)=L_{\gamma,\eta}(\kappa_j)$.
We note that the asymmetrical L\'evy functions,
with $\eta\ne \pm1,0$, only rarely find
their applications in the literature.
The characteristic exponent $\gamma$ depends only on the general features
of the model (namely on $d$ and $\delta$). In contrast the
L\'evy index of symmetry  $\eta$ depends on the details
of the model and on system
parameters $(T,\mbox{etc})$. For $j=2$ we
have $H_j=|H_j|$ and then $\eta=1$ so
the L\'evy density is one sided, as is expected since
$\kappa_2 > 0$.

 As mentioned, in the fast modulation limit, 
the random line width in Eq. (\ref{eqGamma}) characterizes the
statistical properties of the spectral lines. 
Using the approach in Eqs. (\ref{eq8}-\ref{eq10})
one can show that $P(\tilde{\Gamma})=L_{d/(2 \delta),1}(\tilde{\Gamma})$
with the scale parameter $z_{d/(2 \delta)}$ given by Eq. (\ref{eq10})
with $j=2$ and $H_2=p(1-p)/K$.
 
 In what follows we exhibit our results and compare to 
simulations based on the standard
tunneling model of low temperature glass \cite{AHV}. We use system parameters
given by Geva and  Skinner \cite{Geva} to model
terrylene in polystryrene. 
The SM-TLS interaction is dipolar,
hence $\delta =3$, and we consider spatial dimension $d=3$.
The distribution
of the asymmetry parameter and tunneling element 
is $P(A)P(J)=N^{-1} J^{-1}$ for 
$2.8\times 10^{-7}{\rm K} < J <18 {\rm K}$ and
$ 0< A < 17{\rm K}$, $N$ denoting a normalization constant. 
We use $f(A,J)=A/E$ and define a  TLS to be active if
$K> 1/\tau$, $\tau=120$ sec is the time of experimental observation.
In this way the averaging $\left\langle \cdots \right\rangle_{AJ}$ becomes
$\tau$ independent. 
The rate of the TLS is given by 
$K = c J^2 E \coth\left({\beta E_n/2} \right)$ 
and $c=3.9 \times 10^{8} {\rm K}^{-3}$Hz 
is the TLS phonon coupling constant.
Additional system parameters are the coupling constant 
$\alpha=3.75\times 10^{11}$ $\mbox{nm}^{3}$Hz and the
TLS density $1.15\times 10^{-2}$ nm$^{-3}$.
According to
Eqs. (\ref{eq9})-(\ref{eq12}), only the scale parameter $z_{\gamma}$
depends on 
the orientation of the TLS and SM, through $\Psi(\Omega)$.
It is therefore reasonable to assume simple forms for $\Psi(\Omega)$.
We consider two examples,  
model $1$  ($\mbox{M}1$) for which  $\Psi(\Omega)$
is replaced with a two state variable (i.e., a spin model)
$\Psi=1$ or $\Psi=-1$ with equal probabilities of occurrence
and model $2$ ($\mbox{M}2$)
$\Psi(\Omega)=\cos(\theta)$, with $\theta$, the standard polar
coordinate, distributed uniformly.
With these definitions we calculate the symmetry index
$\eta$ and the scaling parameter $z_{\gamma}$ and compare between
the theory and numerical simulation.



%here prediction is distribution of kappa_1 is
%(1/Pi) z_1/(w^2 + z_1^2) gamma is calculated
%using Mathematica. We present four cases in a scaling
%form. 1. T=1.7 and epsilon = +- model       z_1=9.07357  GHz
%      2. T=1.7 and Psi(Omega) = cos(theta)  z_1=9.07/2
%      3. T=1.7 and epsilon = +- model
%      4. T=1 and Psi(Omega) = cos(theta)
%all curves collapse on a master curve (scaling)
%see more details fig2 below and Zu*/1math as well
%as Zu*/Pr*/cumu.f
% for the Graphics see Zu*/Pr*R/Graphics Fig1

\begin{figure}[htb]
\epsfxsize=20\baselineskip
\centerline{\vbox{
      \epsffile{diskap1.ps}  }}
\caption {%\protect\footnotesize
Scaled probability density of first cumulant
$P(\kappa_1)z_1$ versus  $\kappa_1/z_1$.
Symbols are the simulation results obtained from $4000$
molecules for different cases indicated in the figure.
%The star (dot)
%represents model 1 (model 2) for the case $T=1.7 {\bf K}$.
%while the plus  (diamond) represents model 1 (model 2)  
%and $T=1 {\bf K}$. 
The theory, plotted as a solid curve, 
predicts a Lorentzian density $P(\kappa_1)=L_{1,0}(\kappa_1)$
with a scaling parameter $z_1$ which varies from one set of data
to the other.
}
\label{fig1}
\end{figure}

 We consider the first two cumulants $\kappa_1$ and $\kappa_2$
(i.e, the line shape mean and variance). 
%We remind the reader
%that for this case our Eqs. (\ref{eq9}) and (\ref{eq11}) are exact
%within the framework of the model. 
Since $d=\delta$ we find $P(\kappa_1)=L_{1,0}(\kappa_1)$
which is the Lorentzian density,
and $P(\kappa_2) = L_{1/2,1}(\kappa_2)$ which is 
Smirnov's
density. We have considered two temperatures for the two models
M1 and M2. As shown in Fig. \ref{fig1} and \ref{fig2},
a scaling behavior is observed and all data collapse on the
L\'evy densities $L_{1,0}(\kappa_1)$ and $L_{1/2,1}(\kappa_2)$ respectively.
In Fig. \ref{fig1} and \ref{fig2} we have rejected TLSs within a sphere
of radius $r_{\rm min}=1$ nm,
demonstrating that our results are not sensitive to a short cutoff.
Also shown in the inset of Fig. \ref{fig2} is 
$P(\mbox{Re}[\kappa_3])$ which is distributed
according to $L_{1/3,0}(\mbox{Re}[\kappa_3])$ and
a scale parameter $z_{1/3}$ given in 
Eq. (\ref{eq10}).
The L\'evy behavior of $\kappa_1,\kappa_2$ and
$\mbox{Re}[\kappa_3]$
holds generally and is not limited to the slow modulation limit
since these random variables do not depend explicitly on the 
rates $K$. 

\begin{figure}[htb]
\epsfxsize=20\baselineskip
\centerline{\vbox{
      \epsffile{fig2.ps}  }}
\caption {%\protect\footnotesize
Same as Fig. (\protect{\ref{fig1}}) 
for the second cumulant. We show $P(\kappa_2)z_{1/2}^2$ versus
$\kappa_2/z_{1/2}^2$. The solid curve is Smirnov's density 
$L_{1/2,1}(\kappa_2)$. 
In the inset we show
the same as Fig. (\protect{\ref{fig1}}) for $\mbox{Re}[\kappa_3]]$,
the solid curve is
L\'evy density $L_{1/3,0}(\mbox{Re}[\kappa_{1/3}])$.
}
\label{fig2}
\end{figure}



 Consider the distribution of $\mbox{Re}[\kappa_4]$,
which in the slow modulation limit is distributed according
to $L_{1/4,\eta}(\mbox{Re}[\kappa_4])$, 
Eq. (\ref{eq11}). The question
remains if such a slow modulation limit is valid for the standard
tunneling model parameters we are considering. The slow modulation limit
is expected to work when $K\ll |\nu|$. For large enough $r$
this inequality will fail; however,
depending on system parameters, we expect that
contributions from TLS situated far from the SM
are negligible. We also note
that according to the standard tunneling model the rates
$K$ are distributed over a broad range, albeit with
finite cutoffs that insure that the averaged rate 
is finite.
% Therefore the existence of a slow modulation limit
%cannot be generally true. 
To check if the slow modulation limit
is compatible with the standard tunneling model approach,
we compare 
our slow modulation results 
with those obtained by simulation 
in Fig. \ref{fig4}. 
We also show
simulation results in which all rates are set to zero ($K=0$).
For model
$\mbox{M}1$, we find that the deviation between simulation
and theory is small so the assumption of slow modulation
limit is justified. For model $\mbox{M} 2$, we see slightly
larger deviations between the theory and numerical results,
due to
the angular dependence of model $\mbox{M} 2$,
$\Psi(\Omega) = \cos(\theta)$, which reduces the typical frequency
shift $|\nu|$ compared to model $\mbox{M}1$. We conclude that
the present theory 
can be used as a criterion for the validity 
of the slow modulation limit. 

%\begin{figure}[htb]
%\epsfxsize=20\baselineskip
%\centerline{\vbox{
%      \epsffile{disrealkap3.ps}  }}
%\caption {%\protect\footnotesize
%The same as Fig. (\protect{\ref{fig1}}) 
%for $\mbox{Re}[\kappa_3]$. The probability density
%function $P(\kappa_3)$ is described by the L\'evy stable 
%density $L_{1/3,0}(\mbox{Re}[\kappa_3])$. 
%$8000$ SM were considered and $r_{min}=0$.
%}
%\label{fig3}
%\end{figure}
%
\begin{figure}[htb]
\epsfxsize=20\baselineskip
\centerline{\vbox{
      \epsffile{diskap4.ps}  }}
\caption {%\protect\footnotesize
Same as Fig. 
(\protect{\ref{fig1}}) 
for $Re[\kappa_4]$ and for temperature $T=1.7 {\rm K}$.
The symbols are the simulation results obtained for four different
cases as indicated in the figure and for $12000$
molecules.
The solid curve is $L_{1/4,\eta}[\mbox{Re}(\kappa_4)]$ with
an index of asymmetry $\eta=0.6104$.
In the inset we show the power law
tail of the scaled probability density on a log-log plot.
}
\label{fig4}
\end{figure}

 Depending on system parameters, L\'evy statistics
may become sensitive to the finite cutoff  $r_{\rm min}$,
Physically, the cutoff can be important since the
power low interaction is not supposed to work well
for short distances \cite{Pfluegl}.
Our results were derived for $r_{\rm min}=0$,
while for small though finite $r_{\rm min}$ 
one can find intermittency behavior, i.e.,
the ratio $\langle \kappa_2^2 \rangle / \langle \kappa_2\rangle^2$
(as well as similar dimensionless ratios) is very
large. 
When $r_{\rm min}$ is large one finds a Gaussian behavior.
The phenomena of intermittency in the context of a reaction
of a SM in a random environment was investigated
in \cite{Wang}.
Generally high order cumulants are more sensitive to finite
cutoff and for results in Fig (\ref{fig4}) $r_{\rm min}=0$
was chosen to see the proper decay laws in the
wing.

% The L\'evy behavior found in the present study can
%be partially understood using a
%simple argument. 
%Consider the cases for which the cumulants  can be written as 
%$\kappa_j=\sum_n x_n y_n$ where $x_n=1/r_n^{\delta j}$
%and $y_n$ are real random variables defined in Table 1.
%It is easy to show,
%based on the assumption of uniform distribution of TLSs,
%that the random variable $x_n$ is distributed 
%according $P(x)\sim x^{-( 1 + \gamma)}$. 
%This is why we found a L\'evy behavior 
%with characteristic exponent $\gamma$
%which depends only on the dimensionality and the
%interaction exponent.
%The random variables $\{ y_n \}$
%determine the 
%non universal features of the model, i.e.,  the scale parameter
%$z_{\gamma}$ and the index of symmetry $\eta$ 
%(provided that the distribution of $y$ is not too broad).

% For $j=1$ our theory is a (slight) generalization
%of Stoneham's \cite{Stoneham} theory of inhomogeneous 
%line broadening. Briefly,  this theory \cite{Stoneham}
%treats an ensemble of molecules as though each member of
%the ensemble is described by a line of zero width.
%This corresponds to the assumption $\kappa_j=0$ for $j>1$.
%On the other hand our theory shows that distribution of $\kappa_j$
%is wide and in fact even the average $\langle \kappa_j \rangle$ 
%does not generally exist. It is interesting to note that
%some parts of the
%theory of inhomogeneous line
%broadening can also be interpreted in terms of
%L\'evy statistics. For example, in \cite{Stoneham} it was shown
%that broadening due to
%Electric fields, (i.e., $\delta=2$,$d=3$,$j=1$) is described
%by a function which we identify as $L_{3/2,0}$, also known as
%Holtsmark distribution.


% The Holtsmark distribution finds its applications in diverse fields
%like turbulence \cite{Turbo} and astronomy \cite{Feller}.
%Furukawa et al \cite{Jap}
%have used M$\ddot{o}$ssbauer spectroscopy
%and found that distribution
%of internal magnetic fields
%in Au(Fe) spin glass alloys
%is L\'evy stable.
%All these effects result from the long range interaction
%of uniformally distributed
% defects, gravitational masses, vortexes, magnetic dipoles etc. 
%Here we went one step further and considered dynamical
%situation in which a function (i.e., the line shape) is the random
%quantity of interest. As shown here 
%% the interplay between dynamics (i.e., the rates
%$K$)  and static  shifts (i.e., $\nu$) is a possible cause for 
%deviations from L\'evy statistics in a fast modulation limit.
%We believe that approaches similar to ours might find their
%application in other fields.

% Finally, we note that L\'evy statistics can be used to analyze
%other statistical properties of SMs
%in disordered media.
%We have shown that, in the slow modulation limit,
%the averaged line shape is 
%a L\'evy stable function. 
%Other results based on models similar to ours
%\cite{Fleury,Pfluegl,Geva1} can be explained based upon 
%generalized central limit theorem. 
%We shall address these issues in future
%publication.

 To conclude, we showed that the generalized central limit theorem
can be used to analyze distribution of cumulants of SM
line shapes in glass. 
We note that besides cumulants,
L\'evy statistics can be used to analyze other
statistical properties of SMs in disordered media \cite{Barkai2}. 

{\bf Acknowledgment} EB thanks the ETH and Prof. Wild for
their hospitality. This work was supported by the NSF.

\begin{thebibliography}{99}

\bibitem{MO}  W. E. Moerner and M. Orrit,
{\em Science} {\bf 283} 5408 (1999).
T. Plakhotnik, E. Donley and U. P. Wild,
{\em Annu. Rev. Phys. Chem.} {\bf 48}, 181 (1997)
and references therein.


\bibitem{Tittel} L. Fleury, A. Zumbusch, M. Orrit, R. Brown
and J. Bernard, {\em J. Lumin.} {\bf 56} 15
(1993).
J. Tittel, R. Kettner, Th. Basche, C. Brauchle,
H. Quante, K. Mullen,
{\em J. Lumin.} {\bf 64} 1
(1995).
M. Vacha, Y. Liv, H. Nakatsuka, T. Tani,
{\em  J. Chem. Phys.}  {\bf 106} 8324 (1997).
B. Kozankiewicz, J. Bernard and  M. Orrit,  {\em J. Chem. Phys,
} {\bf 101} 9377 (1994).

\bibitem{Reilly2} P. D. Reilly and J. L. Skinner,
{\em Phys. Rev. Lett.}  {\bf 71} 4257 (1993).

\bibitem{Zumofen}
G. Zumofen and J. Klafter,
{\em Chem. Phys. Lett.}  {\bf 219} 303 (1994).

\bibitem{Geva}  E. Geva and J. L. Skinner,
{\em J.  Phys.  Chem. B}  {\bf 101} 8920 (1997).
E. Geva and J. L. Skinner,
{\em Chem. Phys. Lett.}  {\bf 287} 125 (1998).

\bibitem{Pfluegl} W. Pfluegl, F. L. H. Brown and R. J. Silbey,
{\em J.  Chem. Phys.} {\bf 108} 6876 (1998).

\bibitem{Orrit} A. M. Boiron, P. Tamarat, B. Lounis, R. Brown and M. Orrit,
{\em J. Chem. Phys.} {\bf 247} 119 (1999).

\bibitem{Barkai1}
E. Barkai and R. Silbey,
Chem. Phys. Lett., {\bf 310} 287 (1999).

\bibitem{Barkai2}
E. Barkai and R. Silbey, (submitted).

\bibitem{AK} P. W. Anderson, {\em J. Phys. Soc. Jpn.} {\bf 9}, 316 (1954).
R. Kubo, {\em ibid} {\bf 6} 935 (1954).
% R. Kubo Adv. Chem. Phys. {\bf 15} 101 (1969)

% Spectroscopy of a chromophore coupled
% to a lattice of dynamic two level systems 1. Absorption
% line shape
% 959 964
\bibitem{Reilly} P. D. Reilly and J. L. Skinner,
{\em J. Chem. Phys.}  {\bf 101} (2) 959 (1994).

\bibitem{Feller} W. Feller,{\em  An introduction to probability Theory and
Its Applications} Vol. 2 (John Wiley and Sons 1970).

\bibitem{Turbo} S. G. Llewellyn Smith and S. T. Gille {\em Phys. Rev. Lett.}
{\bf 81} 5249 (1998) 
and references therein.


\bibitem{Jap}
 Y. Furukawa, Y. Nakai and N. Kunitomi,
{\em J. of the Phys. Soc. Japan} {\bf 62} No.1. 306 (1993).

% Theory of single molecule optical line shape
% distributions in low temperature glasses
% p 8920 8932

%\bibitem{Brown}  F. L. H. Brown and R. J. Silbey
%{\em J.  Chem. Phys.} {\bf 108}  7434 (1998)


%\bibitem{remark1}
%In fact the line width is not well defined mathematically,
%and there is no general formula for line width of a line shape
%given in Eqs. (\ref{eq4},\ref{eq5})

%\bibitem{Cohen} Cohen Tan {\em Phys. Rev. Lett.} {\bf 83} 3696 (1999)

%\bibitem{Hovi} J. P. Hovi, A. Aharony, D. Stauffer and B. B. Mandelbrot, {
%\em Phys. Rev. Lett.}, {\bf 77}, 877, (1996).

\bibitem{Stoneham} A. M. Stoneham, {\em Rev. Mod. Phys.} {\bf 41}
82 (1969).

\bibitem{AHV} P. W. Anderson, B. I. Halperin, C. M.  Varma,  
{\em Philos. Mag.} {\bf 25}, 1 (1971). W. A. Philips,
J. Low. Temp. Phys. {\bf 7}, 351 (1972). 

\bibitem{remark} In Table 1 we use $\omega_0=0$, for $\omega_0\ne 0$
we have $\kappa_1=\sum_n p_n \nu_n + \omega_0$ while the higher order 
cumulants are $\omega_0$ independent. The distribution of $\omega_0$ depends
on the static disorder \cite{Stoneham} while we are considering
the dynamic disorder due to the flipping TLSs. 

\bibitem{Wang} J. Wang and P.  Wolynes, {\em Phys. Rev. Lett.} {\bf 74} 
4317
(1995).
 

\end{thebibliography}

\end{document}

cond-mat0002045.extracted_bib
\begin{thebibliography}{99}

\bibitem{MO}  W. E. Moerner and M. Orrit,
{\em Science} {\bf 283} 5408 (1999).
T. Plakhotnik, E. Donley and U. P. Wild,
{\em Annu. Rev. Phys. Chem.} {\bf 48}, 181 (1997)
and references therein.


\bibitem{Tittel} L. Fleury, A. Zumbusch, M. Orrit, R. Brown
and J. Bernard, {\em J. Lumin.} {\bf 56} 15
(1993).
J. Tittel, R. Kettner, Th. Basche, C. Brauchle,
H. Quante, K. Mullen,
{\em J. Lumin.} {\bf 64} 1
(1995).
M. Vacha, Y. Liv, H. Nakatsuka, T. Tani,
{\em  J. Chem. Phys.}  {\bf 106} 8324 (1997).
B. Kozankiewicz, J. Bernard and  M. Orrit,  {\em J. Chem. Phys,
} {\bf 101} 9377 (1994).

\bibitem{Reilly2} P. D. Reilly and J. L. Skinner,
{\em Phys. Rev. Lett.}  {\bf 71} 4257 (1993).

\bibitem{Zumofen}
G. Zumofen and J. Klafter,
{\em Chem. Phys. Lett.}  {\bf 219} 303 (1994).

\bibitem{Geva}  E. Geva and J. L. Skinner,
{\em J.  Phys.  Chem. B}  {\bf 101} 8920 (1997).
E. Geva and J. L. Skinner,
{\em Chem. Phys. Lett.}  {\bf 287} 125 (1998).

\bibitem{Pfluegl} W. Pfluegl, F. L. H. Brown and R. J. Silbey,
{\em J.  Chem. Phys.} {\bf 108} 6876 (1998).

\bibitem{Orrit} A. M. Boiron, P. Tamarat, B. Lounis, R. Brown and M. Orrit,
{\em J. Chem. Phys.} {\bf 247} 119 (1999).

\bibitem{Barkai1}
E. Barkai and R. Silbey,
Chem. Phys. Lett., {\bf 310} 287 (1999).

\bibitem{Barkai2}
E. Barkai and R. Silbey, (submitted).

\bibitem{AK} P. W. Anderson, {\em J. Phys. Soc. Jpn.} {\bf 9}, 316 (1954).
R. Kubo, {\em ibid} {\bf 6} 935 (1954).
% R. Kubo Adv. Chem. Phys. {\bf 15} 101 (1969)

% Spectroscopy of a chromophore coupled
% to a lattice of dynamic two level systems 1. Absorption
% line shape
% 959 964
\bibitem{Reilly} P. D. Reilly and J. L. Skinner,
{\em J. Chem. Phys.}  {\bf 101} (2) 959 (1994).

\bibitem{Feller} W. Feller,{\em  An introduction to probability Theory and
Its Applications} Vol. 2 (John Wiley and Sons 1970).

\bibitem{Turbo} S. G. Llewellyn Smith and S. T. Gille {\em Phys. Rev. Lett.}
{\bf 81} 5249 (1998) 
and references therein.


\bibitem{Jap}
 Y. Furukawa, Y. Nakai and N. Kunitomi,
{\em J. of the Phys. Soc. Japan} {\bf 62} No.1. 306 (1993).

% Theory of single molecule optical line shape
% distributions in low temperature glasses
% p 8920 8932

%\bibitem{Brown}  F. L. H. Brown and R. J. Silbey
%{\em J.  Chem. Phys.} {\bf 108}  7434 (1998)


%\bibitem{remark1}
%In fact the line width is not well defined mathematically,
%and there is no general formula for line width of a line shape
%given in Eqs. (\ref{eq4},\ref{eq5})

%\bibitem{Cohen} Cohen Tan {\em Phys. Rev. Lett.} {\bf 83} 3696 (1999)

%\bibitem{Hovi} J. P. Hovi, A. Aharony, D. Stauffer and B. B. Mandelbrot, {
%\em Phys. Rev. Lett.}, {\bf 77}, 877, (1996).

\bibitem{Stoneham} A. M. Stoneham, {\em Rev. Mod. Phys.} {\bf 41}
82 (1969).

\bibitem{AHV} P. W. Anderson, B. I. Halperin, C. M.  Varma,  
{\em Philos. Mag.} {\bf 25}, 1 (1971). W. A. Philips,
J. Low. Temp. Phys. {\bf 7}, 351 (1972). 

\bibitem{remark} In Table 1 we use $\omega_0=0$, for $\omega_0\ne 0$
we have $\kappa_1=\sum_n p_n \nu_n + \omega_0$ while the higher order 
cumulants are $\omega_0$ independent. The distribution of $\omega_0$ depends
on the static disorder \cite{Stoneham} while we are considering
the dynamic disorder due to the flipping TLSs. 

\bibitem{Wang} J. Wang and P.  Wolynes, {\em Phys. Rev. Lett.} {\bf 74} 
4317
(1995).
 

\end{thebibliography}