%% ****** Start of file template.aps ****** %
%%
%%
%%   This file is part of the APS files in the REVTeX 4 distribution.
%%   Version 4.0 of REVTeX, August 2001
%%
%%
%%   Copyright (c) 2001 The American Physical Society.
%%
%%   See the REVTeX 4 README file for restrictions and more information.
%%
%
% This is a template for producing manuscripts for use with REVTEX 4.0
% Copy this file to another name and then work on that file.
% That way, you always have this original template file to use.
%
% Group addresses by affiliation; use superscriptaddress for long
% author lists, or if there are many overlapping affiliations.
% For Phys. Rev. appearance, change preprint to twocolumn.
% Choose pra, prb, prc, prd, pre, prl, prstab, or rmp for journal
%  Add 'draft' option to mark overfull boxes with black boxes
%  Add 'showpacs' option to make PACS codes appear
%\documentclass[aps,pra,twocolumn,superscriptaddress,floatfix]{revtex4-1}  % for review and submission
\documentclass[aip,apl,twocolumn,superscriptaddress,floatfix,reprint,numerical]{revtex4-1}
%\documentclass[aps,preprint,showpacs,superscriptaddress,groupedaddress]{revtex4}
%% for double-spaced preprint
% \widetext
\usepackage{graphicx}  % needed for figures
\usepackage[breaklinks]{hyperref}
\usepackage{amsmath}
% \usepackage{dcolumn}   % needed for some tables
% \usepackage{bm}        % for math
% \usepackage{amsmath}   % for math
% \usepackage{caption}
% \captionsetup{compatibility=false}
%\usepackage{stfloats}
% \usepackage{dblfnote}
%\usepackage[bottom]{footmisc}
% avoids incorrect hyphenation, added Nov/08 by SSR
\hyphenation{ALPGEN}
\hyphenation{EVTGEN}
\hyphenation{PYTHIA}


\begin{document}

\title{Random numbers from vacuum fluctuations}

\author{Yicheng Shi}
\affiliation{Department of Physics, National University of Singapore, 2
  Science Drive 3,  Singapore, 117542}
\affiliation{Center for Quantum Technologies, National University of
  Singapore, 3 Science Drive 2, Singapore, 117543}
\author{Brenda Chng}
\affiliation{Center for Quantum Technologies, National University of
  Singapore, 3 Science Drive 2, Singapore, 117543}
\author{Christian Kurtsiefer}
\affiliation{Department of Physics, National University of Singapore, 2
  Science Drive 3,  Singapore, 117542}
\affiliation{Center for Quantum Technologies, National University of
  Singapore, 3 Science Drive 2, Singapore, 117543}
\email[]{christian.kurtsiefer@gmail.com}
\date{\today}

\begin{abstract}
We implement a quantum random number generator based on a balanced
homodyne measurement of vacuum fluctuations of the electromagnetic field.
The digitized signal is directly processed with a fast randomness
extraction scheme based on a linear feedback shift register. The random bit
stream is continuously read in a computer at a rate of about 480 Mbit/s and
passes an extended test suite for random numbers.
\end{abstract}

% \pacs{}
\maketitle

\section{Introduction}
Various cryptographic schemes, classical or quantum, require high
quality and trusted random numbers for key generation and other aspects of the
protocols. 
In order to keep up with data rates in
modern communication schemes, these random numbers need to be generated at a
high rate~\cite{Gisin2002a}.
Equally, large amounts of random numbers are at the core of Monte Carlo
simulations~\cite{Metropolis:1987}.  
Algorithmically generated pseudo-random numbers are available at very high
rates, but are deterministic by definition and therefore
unsuitable for cryptographic purposes. For applications that require 
unpredictable random numbers, hardware random
number generators have been used in the past~\cite{GALTON1890} and more
recently~\cite{INTEL}. These involve measuring noisy physical processes, and
conversion of the outcome into random numbers. Since it is either practically
(e.g. for thermal noise sources) or fundamentally 
%(for certain quantum processes) 
impossible to predict the outcome of such measurements, these 
physically generated random numbers are considered ``truly'' random.

Quantum random number generators (QRNG) belong to a class of hardware random
number generators where the source of randomness is the fundamentally unpredictable
outcome of quantum measurements. Early QRNGs were based on observing 
the decay statistics of radioactive nuclei \cite{gude1985,figotin2004random}. 
More recently, similar QRNGs based on Poisson statistics in optical photon
detection have been reported~\cite{Stipcevic:07,kwiat:09,Fuerst2010,Wayne:10,Wahl2011,Nie2014,Ren2011}. Different schemes use the 
randomness of a single photon scattered by a beam splitter into either of two
output ports~\cite{Jennewein2000,Stefanov2000}.
As the reflection/transmission of the photon is intrinsically 
random due to the quantum nature of the process, the unpredictability of the
generated numbers is ensured~\cite{Frauchiger2013}.
Other implementations of QRNGs measure the amplified spontaneous emission~\cite{Williams2010}, the
vacuum fluctuations of the electromagnetic field~\cite{Gabriel2010,Syed2011,Shen2010}, or
the intensity~\cite{Kanter2009,Sanguinetti2014} and phase noise of different light sources~\cite{Nie2015,Qi2010,Xu2012,Abellan2014,PhysRevLett.115.250403,Yuan2014,Zhou2015,Jofre2011}. 

In this paper we report on a QRNG based on
measuring vacuum fluctuations of a light filed as the source of
ramdomness~\cite{Gabriel2010,Syed2011,Shen2010}. Such measurements have a very
high bandwidth compared to schemes based on photon
counting~\cite{Stipcevic:07,kwiat:09,Fuerst2010,Wayne:10,Wahl2011,Nie2014,Ren2011}, and have a much simpler 
optical setup compared to phase noise
measurements~\cite{Nie2015,Qi2010,Xu2012,Abellan2014,PhysRevLett.115.250403,Yuan2014,Jofre2011,Zhou2015}. Coupled
with an efficient randomness extractor, we obtain an unbiased, uncorrelated
stream of random bits at a high rate.

\section{Implementation}
Figure~\ref{fig:schematic} schematically shows the setup of our QRNG. A
continuous wave laser
(wavelength 780\,nm) is used as the local oscillator (LO) for the vacuum
fluctuations entering the beam splitter at the
empty port. The output of the beam splitter is directed onto two 
photodiodes, and the photocurrent difference is processed further. 
\begin{figure}
\centerline{
  \includegraphics[width=0.7\columnwidth]{./plot/setup/setup1}
}
\caption{Schematic of the quantum random number generator. A polarizing beam
  splitter (PBS) distributes light from a 780\,nm laser equally onto two
  photodiodes, generating photocurrents $i_1$ and 
  $i_2$. The current difference $i_1-i_2$ is
  amplified, digitized, and processed to generate
  random numbers.}
\label{fig:schematic}
\end{figure}
This setup is known as a balanced homodyne
detector~\cite{jakeman1974,yuen1983} and maps the electrical field in the
second mode entering the beam splitter to the photocurrent difference
$i_1-i_2$. Here, the second input port is empty, so the homodyne measurement
is probing the vacuum state of the electromagnetic field. This field
fluctuates~\cite{Glauber:1963}, and is used as the source of randomness. As
the vacuum field is independent of external physical quantities, it can not be 
tampered with. Since the optical power impinging
on the two photodiodes is balanced, any power
fluctuation in the local oscillator will be simultaneously detected, and therefore cancel in the photocurrent
difference~\cite{yuen1983,schumaker1984}. In  
an alternative view, the laser beam can be seen as generating photocurrents
$i_1,i_2$ with a shot noise power proportional to the average optical power. The
shot noise currents from the diodes add up as they are
uncorrelated, while amplitude fluctuations in the laser intensity (referred to
as classical noise) do not affect the photocurrent difference.
%represented by the average current of the photodiodes

The power of the two output ports is balanced by rotating the laser diode
in front of a polarizing beam splitter (PBS). Light leaving 
the PBS is detected by a pair of 
reverse biased silicon pin photodiodes (Hamamatsu S5972) connected in
series to perform the current subtraction. The balancing of photocurrents
is monitored by observing the voltage drop across a resistor $R_{DC}$ providing a
DC path from the common node to ground. 
We achieve a 50\,dB rejection ratio of the classical noise 
from the laser intensity fluctuations by careful balancing. The fluctuations
$\Delta (i_1-i_2)$
above 20\,MHz are amplified by a transimpedance amplifier (Analog Devices
AD8015) followed by two wideband RF gain blocks (Mini Circuits MAR-6). 
The entire amplifier chain has a calculated effective transimpedance of $R_{\rm eff}\approx540\,\pm 118$\,k$\Omega$. 

\begin{figure}
\includegraphics[width=0.85\columnwidth]{./plot/spectrum/spectrum2}
\caption{Amplified noise levels measured into a resolution bandwidth
  $B=1\,$kHz. The total noise is measured from the
  photocurrent difference $i_1-i_2$ with equal optical power impinging 
  on both photodiodes and approaches the theoretical shot noise level 
  of -52\,dBm (dashed trace) given by (\ref{eq:shot_noise}). The 
  current $i_1$ of a single photodiode reveals colored classical 
  noise. The electronic noise is measured without any optical 
  input.}
\label{fig:spectrum}
\end{figure}

To ensure that the fluctuations at the amplifier output are
dominated by quantum noise, the spectral power density is measured (see
Fig.~\ref{fig:spectrum}). With an optical power of 3.1\,mW 
received by each photodiode corresponding to an average photocurrent
$I=1.7$\,mA, we observe a noise power of $P=-53.5\,$dBm (at 75\,MHz) in
a bandwidth of $B=1\,$kHz. This is about 1.5\,dB lower
than the theoretically expected shot noise value (dashed trace)
\begin{equation}\label{eq:shot_noise}
P=\frac{4eIB{R_{\rm eff}}^2}{Z}\approx-52 \, \rm{dBm}\,,
\end{equation}
where $e$ is the electron charge and $Z=50\,\Omega$ the load impedance.\cite{Schottky1918a}
The difference is compatible with uncertainties in determining the
transimpedance of the amplifier.
The measured total noise has a relatively flat power density in the
range of 20 to 120\,MHz, with high pass filters in the circuit
suppressing low frequency fluctuations. The high end of the
pass band is defined by the cutoff frequency of the amplifier. To illustrate
the effectiveness of removing classical noise in the photocurrents, the
spectral power density of the photocurrent generated from a single diode is also
shown. Strong spectral peaks at various radio frequencies appear to 
enter the system probably via the laser diode current. For completeness, the 
spectral power density of the electronic noise is recorded without any 
optical input and found to be at least 10\,dB below the total noise 
level, i.e., the total noise is dominated by quantum fluctuations.

The amplified total noise is digitized into signed 16 bit words
$x_i$ at a sampling rate of 60\,MHz with an analog to digital converter
(ADC, Analog Devices AD9269-65). The sampling rate is lower than 
the cut-off frequency of the noise signal to avoid temporal
correlation between samples. As shown in
Fig.~\ref{fig:autocorrelation}, the normalized autocorrelation
\begin{equation}
A(d)=\langle x_i\,x_{i+d}\rangle_n/\langle x_i^2\rangle_n
\end{equation}
evaluated over $n=10^9$ samples shows that the absolute value of the 
autocorrelation $|A(d)|$ for non-zero delay
($d>0$) is below $1.2\times10^{-3}$, which is slightly smaller than what has been observed in other experiments.\cite{Mitchell2015,Ma2012,Nie2015}
The residual correlation above the $2\sigma$ confidence level for 
$d\lesssim60$ is a consequence of the finite bandwidth of the signal, as stated by
the Wiener-Khinchin theorem.
\begin{figure}
  \includegraphics[width=0.85\columnwidth]{./plot/autocorrelation_1g/autocorrelation_1g}
  \caption{Autocorrelation of the total noise signal sampled at 60\,MHz,
    computed over $10^9$ samples (solid line), compared with the $2\sigma$
    confidence level (dashed line).} 
\label{fig:autocorrelation}
\end{figure}

\section{Entropy estimation}
The total noise measured before the ADC contains both
quantum and electronic noise. To determine how much
randomness from non-classical origin can be safely extracted, it is necessary
to estimate the entropy $H(X_q)$ contributed by the quantum process.


Therefore, we assume that the
measured total noise signal $X_t=X_q+X_e$ is the sum of independent random
variables $X_q$ for the quantum noise, and $X_e$ for the electronic
noise which includes the photodetector, amplifier and digitizer noise\cite{Ma2012, Sanguinetti2014}.
All three variables $X_q$, $X_e$ and $X_t$ are assumed to have discrete values
between $-2^{15}$ and $2^{15}-1$. 
We take the
worst case scenario that an adversary has full knowledge of the electronic
noise, i.e., is able to predict the exact outcome of $X_e$
at any moment. In this case, the amount of quantum-based randomness in the
acquired total noise signal is quantified by the conditional entropy
$H(X_t|X_e)$, i.e., the entropy in the total signal, given full knowledge of
the electronic noise $X_e$. As the variables are assumed to be additive and
independent, the conditional entropy is
$H(X_t|X_e)=H(X_q+X_e|X_e)=H(X_q|X_e)=H(X_q)$.  

\begin{figure}
\includegraphics[width=0.85\columnwidth]{./plot/distribution/distribution}
\caption{Probability distribution of the measured total noise with
  variance ${{\sigma}_t}^2$ (a), electronic noise with variance
  ${{\sigma}_e}^2$ (b), and the estimated quantum noise with variance
  ${{\sigma}_q}^2$ (c). The filled areas in (a), (b) show the actual
  measurements over $10^9$ samples, the solid lines fit to Gaussian 
  distributions.}
\label{fig:distribution}
\end{figure}
The variance of the total noise, $\sigma_t^2$, is given by the sum of the
variances $\sigma_q^2$ for the quantum noise, and $\sigma_e^2$ for the
electronic noise.  Over $10^9$ samples, we find 
$\sigma_t=4504.41$ and $\sigma_e=1481.8$, measured with the laser switched off
(see Fig.~\ref{fig:distribution}).
Note that for the total noise, the observed distribution is slightly skewed
compared to a Gaussian distribution [solid line in
Fig.~\ref{fig:distribution}(a)], possibly due to a distortion in the
digitizer. Assuming the quantum noise $X_q$ has a Gaussian
distribution~\cite{Glauber:1963}, we would assign a variance
$\sigma_q^2=\sigma_t^2-\sigma_e^2\approx 4253.7^2$.
To estimate the entropy for a Gaussian distribution, we use the Shannon 
entropy
\begin{equation}\label{eq:shannon}
H_S(X_q)=\sum_{x=-2^{15}}^{2^{15}-1} -p_{q}(x)\log_2 p_{q}(x)\,,
%\quad{\rm with}\quad
%P_q(x)=\frac{1}{\sqrt{2\pi\sigma_q^2}}e^{-\frac{x^2}{2\sigma_q^2}}\,,
\end{equation}
where $p_q(x)$is the probability distribution of the quantum noise $X_q$. Since $\sigma_q\gg1$, $H_S(X_q)$ can be well approximated by 
\begin{equation} \label{eq:differential_entropy}
\int\limits_{-\infty}^{+\infty}\!\!-{f(x)}\log_2{f(x)}\,\mathrm{d}x=\log_2(\sqrt{2{\pi}e}\,\sigma_q)\,,
\end{equation}
where $f(x)$ is a Gaussian probability density function with variance
$\sigma_q^2$, and $e$ the base of the natural logarithm~\footnote{One can show
  that
  $|H(X_q)-H'(X_q)|<\log_2{(\sqrt{2\pi}\sigma_q)}/(\sqrt{2\pi}\sigma_q)\approx 0.0013$\,bit  for $\sigma_q=4108$.}.
This yields 14.1 bits of entropy per 16 bit sample. We also evaluate
the min-entropy of this distribution,
\begin{equation}\label{eq:min_entropy}
H_{\infty}(X_q)=-\log_2 (\mathrm{max}[p_{q}(x)])\approx\log_2 ({\sqrt{2\pi}\sigma_{q}})\,, 
\end{equation}
where $\mathrm{max}[p_{q}(x)]$is the maximum value of the probability distribution of
$X_q$. This yields a min-entropy of 13.4 bits per 16 bit sample.

The Shannon entropy $H_S(X_q)$ serves as an upper
bound of extractable randomness, %, i.e. the maximum entropy one can
%extract from the source of randomness under the assumption of a Gaussian
%distribution of the independent random variables $X_q$ and $X_e$. 
while the min-entropy sets a lower bound, i.e.
the least amount of randomness possessed by each sample. An
alternative estimation of the entropy in $X_q$ assumes that electronic noise is
not only known to a third party, but also could be tampered
with\cite{Syed2011, Haw2014, Mitchell2015}.

\section{Randomness extraction}
In many applications, random numbers are required to be not only
unpredictable, but also uniformly distributed. As such, the raw ADC output 
%data at the amplifier output 
cannot be directly used. %since they are non-uniformly distributed. 
Randomness extractors %is the essential process required to
convert non-uniformly distributed raw data into a uniformly distributed binary
stream without correlations~\cite{miklos:86}. Although there is no
deterministic universal randomness extractor\cite{miklos:86},  
various practical implementations have been reported.
Examples are Trevisan's extractor, a Toeplitz hashing extractor\cite{Ma2012},
random matrix multiplications\cite{Sanguinetti2014,Frauchiger2013a}, or 
a family of secure hashing algorithms (SHA)~\cite{Gabriel2010}. 

\begin{figure}
\includegraphics[width=0.9\columnwidth]{./plot/new_extractor/extractor1}
\caption{Schematic of a LFSR-based randomness extractor. Eight bits (from 
the shaded positions) are extracted for every 16 bits of input.}
\label{fig:extractor}
\end{figure}
In this work, we use a randomness extractor based on a Linear
Feedback Shift Register (LFSR) as shown in Fig.~\ref{fig:extractor},
equivalent to a cyclic redundancy check (CRC)~\cite{CRC}. The LFSRs are well
known for generating long pseudo-random streams with little 
computational resources, and are in widespread use in communication 
applications for spectrum 
whitening~\cite{Krawczyk1994,Barkan2007,raey,wells2004random,4271377}.

We use a maximum length LFSR with 63 cells and a 
two-element feedback path. Its state %The state of LFSR stored in the shift 
%register 
at any time $t$ is a row vector 
$S_t$ of 63 bits, with a recursion relation
\begin{align}
S_{t+1} &= 	S_{t} M + R_{t}\label{eq:lfsrRecursion0}\\
		&= 	\begin{pmatrix}
			s_0, & s_1, & s_2, & \cdots & s_{62}
			\end{pmatrix}
			\begin{pmatrix}
			0 & 0 & \cdots & 0 & 1 \\
			1 & 0 & \cdots & 0 & 1 \\
			0 & 1 & \cdots & 0 & 0 \\
			\vdots  & \vdots  & \ddots & \vdots & \vdots \\
			0 & 0 & \cdots & 1 & 0
			\end{pmatrix}\nonumber\\
			& \quad +
			\begin{pmatrix}
			0, & 0, & 0, & \cdots & 0,& r_t
			\end{pmatrix}\nonumber\\
		&=	\begin{pmatrix}
			s_1, & s_2, & s_3, & \cdots & s_{62},& s_{0} + s_{1} + r_{t}
			\end{pmatrix}\,,\nonumber
\end{align}
where an elementary addition represents a binary xor, and a multiplication a
binary and operation.
%where
%\begin{align}
%s_{63} = s_{0} + s_{1} + r_{t}\,.
%\end{align}

The $63 \times 63$ matrix $M$ represents the shift and feedback 
operation on the LFSR state. The addition of row vector $R_{t}$ describes
the injection of one raw random bit $r_t$ into $S_t$.
After $n$ cycles, the LFSR state becomes

\begin{align}
S_{t+n} = S_{t}M^{n} + \underbrace{R_{t}M^{n-1} + \cdots + 
R_{t+n-1}}_{A}	\label{eq:lfsrRecursion1}\,.
\end{align}
%\textbf{
Row vector $A$ %in (\ref{eq:lfsrRecursion1}) can be alternatively 
can be expressed as a matrix product
%}
\begin{align}
A &= \begin{pmatrix}
	r_t, & r_{t+1}, & \cdots, & r_{t+n-1}
	\end{pmatrix}
	T\,,
\end{align}
with 
\begin{align}
T  =	  
	  \begin{pmatrix}
	S'M^{n-1}\\
	S'M^{n-2}\\
	\vdots\\
	S'I
	\end{pmatrix}, \quad
	S' = \begin{pmatrix}
	  0, & 0, & \cdots & 0,& 1
	  \end{pmatrix}\,.\label{eq:lfsrRecursion2}
\end{align}
%\textbf{
Matrix $T$ in (\ref{eq:lfsrRecursion2}) is a $63 \times 63$ Toeplitz matrix
with rows generated from a LFSR sequence (\ref{eq:lfsrRecursion0}) with
$R_t=0$ and initial state $S_t=S'$. 
It was shown that multiplying an input stream by such a Toeplitz
matrix can be used as a hashing function that generates an almost-uniform
output~\cite{Krawczyk1994}.
%}

%\textbf{
In our setup, we serially inject the 16 bits from each ADC output
word into the LFSR, but extract only 8 bits 
$s_{i}$ provided by stream $S_t$ (at positions $62, 60, \cdots 48$ after the
injection) in a parallelized topology. This is equivalent to a privacy
amplification
process~\cite{Bennett1995}, and ensures that no residual correlations due to
the non-uniform input distribution or any classical noise that may be known to
an adversary are present in the output stream, because
the extraction ratio of $50\%$ is lower than the $13.4/16\approx84\%$ allowed  by the min entropy
(\ref{eq:min_entropy}).% as a lower bound of extractable entropy.
%as an upper bound for extractable entropy ($14.1/16\approx88\%$) and 
%$13.4/16\approx84\%$ given by the Shannon
%entropy (\ref{eq:differential_entropy}) and min-entropy 
%(\ref{eq:min_entropy}), respectively. 
%This ensures that there are no residual correlations left in the output stream.
%}

\begin{figure}
\includegraphics[width=0.8\columnwidth]{./plot/histogram/multiplot}
\caption{Distribution of random data before (blue) and after (red) the randomness extractor, shown in time domain (left) and histogram (right).}
\label{fig:histogram}
\end{figure}

A merit of this extractor is its low complexity. Unlike many other secure
hashing algorithms, it can be easily implemented either in high speed or low
power technology. Therefore, the extraction process does not limit the random
number generation rate. This scheme can be parallelized using 126
regsiter cells, capable of receiving up to 63 injected raw bits
per clock cycle  while still following the extractor 
equation (\ref{eq:lfsrRecursion0}). With a CPLD operating at a 
clock frequency of 400 MHz, this algorithm would be able to process up to
$25\times10^9$ raw input bits per second.

\section{Performance}

To evaluate the quality of the extracted random numbers, we apply two suites
of randomness tests: the statistical test suite from
NIST~\cite{NISTtestsuite}, and the ``Die-harder'' randomness test
battery~\cite{dieharder}. The output of our RNG passed both tests consistently
when evaluated over a sample of 400\,Gigabit in the sense that occasional weak
outcomes of some tests do not repeat.

Our implementation has an output rate of about 480\,Mbit/s of uniformly
distributed random bits, with the digitizer unit sampling at 60\,MHz and
randomness extraction ratio of 50\%; this is limited by the speed of the
data transmission protocol (USB2.0). While significantly higher generation
rates have been reported recently\cite{Nie2015,Williams2010,Xu2012}, our 
design in comparison is simpler both in hard- and software implementation.
With moderate effort, our random number generation rate can be greatly increased by extending
the bandwidth of the photodiodes, amplifiers, 
and digitizer devices, while maintaining the relatively
simple randomness extraction mechanism. Practically, the
resolution-bandwidth product of the ADC limits the random bit
generation rate.
 
\section{Conclusion}
In summary, we demonstrated a random number generation scheme by measuring
the vacuum fluctuations of the electromagnetic field. By estimating the amount
of usable entropy from quantum noise and using an efficient randomness
extractor based on a linear feedback shift register, we can generate
uniformly distributed random numbers at a high rate from a fundamentally
unpredictable quantum measurement.

We acknowledge the support of this work by the National Research Foundation
(partly under grant No. NRF-CRP12-2013-03) \& Ministry of Education in
Singapore, partly through the Academic Research Fund MOE2012-T3-1-009.

%\bibliographystyle{apsrev4-1}
\bibliography{references/reference}

\end{document}
%
% ****** End of file template.aps ******
