%% ****** 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
\documentclass[10pt]{article}
\usepackage[OE]{express}
\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{A compact random number generator based on optical homodyne detection}

\author{Yicheng Shi,\authormark{1} Brenda Chng,\authormark{1}, Christian Kurtsiefer,\authormark{1,2,*}}

\address{\authormark{1}Center for Quantum Technologies, National University of Singapore, 3 Science Drive 2, Singapore, 117543\\
\authormark{2}Department of Physics, National University of Singapore, 2 Science Drive 3, Singapore, 117542\\}

\email{\authormark{*}christian.kurtsiefer@gmail.com}
%\date{\today}
\bibliographystyle{unsrt}
\bibliography{references/reference}

\begin{abstract}
We implement a quantum random number generator based on a balanced homodyne measurement of vacuum fluctuations of the electromagnetic field. We used wave front splitting of the local oscillator instead of amplitude splitting in order to simplify the optical setup. The digitized noise signal is 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}


\section{Introduction}
Generating high quality and trusted random numbers is an essential task in various cryptographic schemes and many other fields such as Monte Carlo simulations~\cite{Metropolis:1987} and various randomized algorithms. Algorithmically generated pseudo-random numbers are available at very high rates and can be easily implemented in software, but they are deterministic in nature and therefore are not suitable for cryptographic purposes. As an alternative, hardware random number generators have been used~\cite{GALTON1890, INTEL}. They measure noisy physical processes and convert the outcome into random numbers. Since it is impossible to predict the outcome of such measurements, these physically generated random numbers are more trusted compared to pseudo-random numbers.

Quantum random number generators (QRNG) is a class of hardware random number generators whose source of randomness is the outcome of quantum measurements. Early implementations of QRNGs made use of decay statistics of radioactive nuclei~\cite{gude1985,figotin2004random}. A number of more recent implementations using quantum optical measurements have been reported. These include measuring photon number statistics~\cite{Stipcevic:07,kwiat:09,Fuerst2010,Wayne:10,Wahl2011,Nie2014,Ren2011}, scattering events of single photons by a beam splitter~\cite{Frauchiger2013}, amplified spontaneous emission of a fiber amplifier~\cite{Williams2010}. QRNGs based on measuring the intensity~\cite{Kanter2009,Sanguinetti2014} and phase noise~\cite{Nie2015,Qi2010,Xu2012,Abellan2014,PhysRevLett.115.250403,
Yuan2014,Zhou2015,Jofre2011} of different light sources have also been reported.

In this paper we report on a QRNG implementation based on measuring the vacuum fluctuations of the electromagnetic field, which has been reported in~\cite{Gabriel2010,Syed2011,Shen2010,Shi2016}. Such measurements are known for their high bandwidth and simple optical setup. In this paper, we simplify the optical setup of the homodyne detector down to only a laser diode and two photodiodes without using light splitting components. Combined with an efficient randomness extractor, we are able to generate unbiased and uncorrelated stream of random bits at a high rate, but now with a much simpler optical setup.

\section{Optical homodyne measurement by wavefront splitting}


The QRNG based on measuring vacuum fluctuations of the electromagnetic field uses an balanced homodyne detector. A conventional setup consists of a laser diode (LD) as a local oscillator, a 50:50 beam splitter (BS), and two photodiodes ($\rm{PD_1}$, $\rm{PD_2}$). A collimation lens and a mirror is also used to steer and guide the laser beam. 

Fig.~\ref{fig:splitting} (b) shows the schematic of a conventional balanced homodyne detector. The local oscillator (LO) in mode $a$ enters the BS and is directed onto two photodiodes and the photocurrent difference is measured. This setup maps fluctuations of electrical field in mode $b$ to the photocurrent difference $i_1 - i_2$. When probing the vacuum fluctuations of EM fields, the input mode $b$ of the BS is kept empty (i.e. mode b is at $|0\rangle$). 

A core requirement for a balanced homodyne detection is the mixing of mode $a$ and $b$ at the beam splitter, which is governed by the following matrix
\begin{equation} \label{eq:beam_splitter_matrix}
\begin{pmatrix}
\tilde{E}_c \\
\tilde{E}_d
\end{pmatrix}
= M
\begin{pmatrix}
\tilde{E}_a \\
\tilde{E}_b
\end{pmatrix}\quad
\rm{where}\quad
M = 
\begin{pmatrix}
1 & -1 \\
1 & 1
\end{pmatrix}
\end{equation}
Where $\tilde{E}_{a,\ b,\ c,\ d}$ represents the oscillating electrical fields at modes $a,\ b,\ c,\ d$. This matrix relation is ensured by the boundary conditions of the electromagnetic fields between the dielectric media of the beam splitter. The mixing of mode $a$ and $b$ also requires good overlap of their spatial profile as well as their frequency. Although measuring the vacuum state does not need to go through this non-trivial process, careful placement and tuning of optical components are still needed.

\begin{figure}[!ht]
\centering\includegraphics[width=0.8\columnwidth]{./plot/wavefront_splitting/wavefront_splitting}
\caption{Splitting mechanisms of the two different implementations. A conventional balanced homodyne detection scheme (a) relies on the beam splitter matrix relation between the input modes $a$, $b$ and the output modes $c$, $d$. In the wavefront-splitting implementation (b), this is replaced by spatially splitting the elliptical transverse mode of a laser beam.}
\label{fig:splitting}
\end{figure}

We can simplify the setup of the balanced homodyne detector by using a different mode decompositon of the local oscillator, thus replacing the beam splitter matrix in (\ref{eq:beam_splitter_matrix}). Fig.\ref{fig:splitting} (b) shows the simplified setup that we propose. A pair of square shaped photodiodes are placed adjacent to each other are directly exposed to the laser beam, each receiving approximately half of the laser beam spot. The optical mode from the laser diode is typically of elliptical transverse profile and for simplicity, the electrical field amplitude at the photodiode surface can be approximated by

\begin{equation}\label{eq:elliptical}
\tilde{E}_{l}(x,y,t) = \hat{E}_{0} \cdot g(x,y,t) = \hat{\epsilon}E_0( e^{-i\omega t} \cdot e^{-{x^2}/{{w_x}^2}} \cdot e^{-{y^2}/{{w_y}^2}})
\end{equation}

where $\epsilon$ is the polarization vector, and $E_0$ is the global field amplitude. The term $e^{-{x^2}/{{w_x}^2}} \cdot e^{-{y^2}/{{w_y}^2}}$ describes a transverse beam profile of a 2D gaussian distribution with different parameters $\omega_x$, $\omega_y$ along $x, y$ directions. The function $g(x,y,t) = e^{-i\omega t} \cdot e^{-{x^2}/{{w_x}^2}} \cdot e^{-{y^2}/{{w_y}^2}}$ describes the optical mode of the local oscillator. 

We now introduce a somewhat similar mode function $h(x,y,t)$ which is manually defined as

\begin{equation}\label{eq:half_elliptical}
h(x,y,t)=g(x,y,t) \cdot
\left\{\begin{array}{lr}
+1 &for \quad y>0 \\
-1 &for \quad y<0
\end{array}
\right.
\end{equation}

It is clear that $h(x,y,t)$ is orthogonal to $g(x,y,t)$ since $\int g(x,y,t) \cdot h(x,y,t)  \mathrm{d} x^3 = 0$, thus the two modes can be considered as independent harmonic oscillators. Consider a second field in mode $h$, $\tilde{E}_{v} = \hat{\epsilon}E_1 \cdot h(x,y,t)$ is being mixed with the local oscillator at mode $g$. The electrical field on the "upper half" ($y>0$) and "down half" ($y<0$) can be written as

\begin{align}\label{eq:mixing}
\tilde{E}_{up} &= \hat{\epsilon}(E_0 + E_1) \cdot \left\{\begin{array}{lr}
g(x,y,t) &for \quad y>0 \\
0 &for \quad y<0
\end{array}
\right. \\ \nonumber
\tilde{E}_{down} &= \hat{\epsilon}(E_0 - E_1) \cdot \left\{\begin{array}{lr}
0 &for \quad y>0 \\
g(x,y,t) &for \quad y<0
\end{array}
\right.
\end{align}

It it easy to see that Equation (\ref{eq:mixing}) reproduce the beam splitter matrix in ~(\ref{eq:beam_splitter_matrix}). The photocurrent difference $i_1 - i_2$ shown in Fig.\ref{fig:splitting} (b) now maps to the vacuum fluctuations of an electromagnetic field with a mode function $h(x,y,t)$, which is kept at vacuum state. The mixing of the two modes is automatically ensured since the vacuum field is present everywhere.

In contrast to a conventional homodyne detector based on amplitude splitting of fields through beam splitter, we use the splitting of wave front of the laser beam and thus eliminated the usage of beam splitting components. The setup is now much simplified compared to those in ~\cite{Gabriel2010,Syed2011,Shen2010,Shi2016} and the non-trivial alignment procedures can be avoided.



\section{Implementation}

\begin{figure}[!ht]
\centering\includegraphics[width=0.5\columnwidth]{./plot/schematic/schematic}
\caption{Amplified noise levels measured into a resolution bandwidth $B = 60$ kHz. The red trace is the amplified photocurrent difference $i_1 - i_2$, with equal optical power impinging on both photodiodes. The blue trace corresponds to the electronic noise which is measured without any optical input.}
\label{fig:schematic}
\end{figure}

Figure~\ref{fig:schematic} schematically shows the balanced homodyne detection setup of our QRNG. A continuous wave laser (wavelength 780\,nm) is used as the local oscillator for the vacuum fluctuations. The beam from the laser diode impinges directly onto a pair of photodiodes (OSRAM SFH2701). The sensitive area of the two photodiodes are two $0.6 \times 0.6$ mm squares and are placed next to each other with a 1 mm gap in between.
\begin{figure}[!ht]
\centering\includegraphics[width=0.7\columnwidth]{./plot/spectrum/spectrum}
\caption{Amplified noise levels measured into a resolution bandwidth $B = 100$ kHz. The red trace is the amplified photocurrent difference $i_1 - i_2$, with equal optical power impinging on both photodiodes. The blue trace corresponds to the electronic noise which is measured without any optical input.}
\label{fig:spectrum}
\end{figure}

The laser diode casts an elliptical beam spot of about 2.5 mm long and 0.7 mm wide, which covers the two photodiodes. A fraction of the optical power is received by the two photodiodes. By carefully adjusting the position of the beam spot, we are able to balance the optical power received by the two diodes. The fluctuations of photocurrent difference $\Delta (i_1-i_2)$ is amplified by a transimpedance amplifier (Analog Devices AD8015) followed by two wideband RF differential amplifiers (Analog Devices AD8351). The entire amplifier chain has a calculated effective transimpedance of $R_{\rm eff}\approx1\,$M$\Omega$. 

Fig~\ref{fig:spectrum} shows the measured total noise output (red trace) which has a relative flat power density range from about 10 MHz to 150 MHz. The lower end of the band is set by the AC coupling capacitors in the gain block while the high end is determined by the cut-off frequency of the amplifiers. As a comparison, the electronic noise (blue trace) is measured with the laser diode switched off. The signal to noise ration is found to be over 10  dB at lower frequencies (0-50 MHz) and around 5 dB at higher frequencies (50-100 MHz) and we conclude that the total noise is dominated by quantum fluctuations. 

\begin{figure}[!ht]
\centering\includegraphics[width=0.7\columnwidth]{./plot/autocorrelation/autocorrelation}
\caption{Autocorrelation of the total noise signal sampled at 200\,MHz, computed over $10^7$ samples (solid line), compared with the $2\sigma$ confidence level (dashed line).}
\label{fig:autocorrelation}
\end{figure}

The amplified noise signal is digitized into signed 12 bit words at a sampling rate of 200 MHz with an analog to digital converter (ADC, Analog Devices AD9634). The normalized autocorrelation evaluated over $10^7$ samples is shown in Fig~\ref{fig:autocorrelation}. The autocorrelation measured up to a delay of $d=100$ is on the order of $10^{-4}$ which is below the $2\sigma$ confidence level. Residual correlation is observed for $d<10$ and is a consequence of the finite bandwidth of the noise signal (150 MHz) and the high sampling rate of the ADC (200 MHz).

\section{Entropy estimation and randomness extraction}
We estimate the entropy in the raw data to help determine the amount of extractable randomness from our QRNG. Two different definitions of entropy are used here. An upper bound of randomness is given by the Shannon entropy, and the min-entropy is computed to set an lower bound.  

For this setup, we followed the same assumptions made in our previous work\cite{Shi2016} which assumes 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. The three variables $X_q$, $X_e$ and $X_t$ are assumed to follow Gaussian distributions and take discrete values between $-2^{11}$ and $2^{11}-1$. Considering the worst case scenario that an adversary has full knowledge of the electronic noise, the conditional Shannon entropy in this case is 
\begin{equation}\label{eq:conditional}
H(X_t|X_e)=H(X_q+X_e|X_e)=H(X_q|X_e)=H(X_q)
\end{equation}
A variance of $\sigma_q^2=\sigma_t^2-\sigma_e^2\approx 531.6^2$ is calculated for $X_q$. For such a Gaussian distribution with $\sigma_q\gg1$, the Shannon entropy can be computed as
\begin{align}\label{eq:shannon}
H_S(X_q)&=\sum_{x=-2^{11}}^{2^{11}-1} -p_{q}(x)\log_2 p_{q}(x)\\ \nonumber
&\approx \int\limits_{-\infty}^{+\infty}\!\!-{f(x)}\log_2{f(x)}\,\mathrm{d}x=\log_2(\sqrt{2{\pi}e}\,\sigma_q)\\ \nonumber
&\approx 11.1 \mathrm{bits}
\end{align}

\begin{figure}[!ht]
\centering\includegraphics[width=0.7\columnwidth]{./plot/histogram/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^7$ samples, the solid lines fit to Gaussian distributions.}
\label{fig:distribution}
\end{figure}

The min-entropy of the quantum noise $X_q$ is computed as
\begin{align}\label{eq:min_entropy}
H_{\infty}(X_q)&=-\log_2 (\mathrm{max}[p_{q}(x)])\\ \nonumber
&\approx\log_2 ({\sqrt{2\pi}\sigma_{q}})\\ \nonumber
&\approx 10.38 \mathrm{bits}
\end{align}

The Shannon entropy $H_S(X_q)$ and min-entropy $H_{\infty}(X_q)$ set up the upper and lower bound of extractable randomness. We use a randomness extractor based on a Linear Feedback Shift Register (LFSR) which has been reported in our previous work~\cite{Shi2016}. 
\iffalse
The randomness extractor utilizes a maximum length LFSR with 63 cells and a two-element feedback path. It's state 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. The addition of row vector $R_{t}$ describes the injection of one raw random bit $r_t$ into $S_t$ at time $t$. 
\fi
The extractor is equivalent to multiplying an input stream of 63 bits to a $63 \times 63$ Toeplitz matrix generated from a LFSR and is shown to be a valid hashing function~\cite{Krawczyk1994}. The low complexity of this extractor allows it to be easily implemented either in high speed or low power technology. This scheme can be par-
allelized using 126 register cells, capable of receiving up
to 63 injected raw bits per clock cycle to even further speed up the process.


\section{Performance}
We apply the statistical test suite from NIST~\cite{NISTtestsuite}, and the ``Die-harder'' randomness test battery~\cite{dieharder} to evaluate the quality of the extracted random numbers. Our RNG output 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 200\,MHz and randomness extraction ratio of 66\%; this is limited by the speed of the data transmission protocol (USB2.0). Although significantly higher generation rates have been reported recently\cite{Nie2015,Williams2010,Xu2012}, our design is by far the most compact and with moderate
effort, our random number generation rate can be greatly
increased by extending the bandwidth of the photodi-
odes, amplifiers, and digitizer devices.

\section{Conclusion}

In summary, we demonstrated a random number generation scheme by measuring the vacuum fluctuations of
the electromagnetic field. By using wave front splitting instead of amplitude splitting, we eliminated the usage of any beam splitting optics in our setup. 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.

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