% ------------------------------------------------------------------------
\documentclass{llncs}
\input{prelude}
\begin{document} 
\title{Star Discrepancies for generalized Halton points}
% \titlerunning{}  % abbreviated title (for running head)
%                                     also used for the TOC unless
%                                     \toctitle is used

\author{Thomas Espitau\inst{1} \and Olivier Marty\inst{1}} 
%
% \authorrunning{} % abbreviated author list (for running head)
%
%%%% list of authors for the TOC (use if author list has to be modified)
% \tocauthor{}
%
\institute{
  $\mbox{}^1$ ENS Cachan \qquad
}

\maketitle

\makeatletter

\makeatother


\begin{abstract}
\end{abstract}


\section{Introduction}

\section{General architecture of the tool}

The testing tool is aimed to be modular: it is made of independents blocks that
are interfaced trough a scheduler. More precisely a master wrapper is written
in Python that calls a first layer which performs the chosen heuristic. This
layer is written in C++ for performances. The given discrepancy algorithm 
--- written in C --- is called when evaluations of a state is needed.
The wrapper dispatch the computations on the multi-core architecture of 
modern computers\footnote{for us, between 2 and 4 physical cores and 4 or 8 
virtual cores}. This basic architecture is described in figure~\ref{main_flow}.
Experiments were conducted on two machines: 
\begin{itemize}
  \item 2.4 GHz Intel Dual Core i5 hyper-threaded to 2.8GHz, 8 Go 1600 MHz DDR3.
  \item 2.8 GHz Intel Quad Core i7 hyper-threaded to 3.1GHz, 8 Go 1600 MHz DDR3.
\end{itemize}

\begin{figure}
\includegraphics[scale=0.6]{main_flow.pdf}
\caption{Tool overview}
  \label{main_flow}
\end{figure}

On these machines, some basic profiling has make clear that 
the main bottleneck of the computations is hiding in the \emph{computation
of the discrepancy}. The chosen algorithm and implantation of this 
cost function is the DEM-algorithm~\cite{Dobkin} of 
\emph{Magnus Wahlstr\o m}~\cite{Magnus}.\medskip

All the experiments has been conducted on dimension 2,3,4 
--- with a fixed Halton basis 7, 13, 29, 3 ---. Some minor tests have
been made in order to discuss the dependency of the discrepancy and 
efficiency of the heuristics with regards to the values chosen for the
prime base. The average results remains roughly identical when taking 
changing these primes and taking them in the range [2, 100]. For such
a reason we decided to pursue the full computations with a fixed 
basis.


\subsection{Algorithmic insights}

To perform an experiment we made up a 
loop above the main algorithm which calls the chosen heuristic multiple 
times in order to smooth up the results and obtain more exploitable datas.
Then an arithmetic mean of the result is performed on the values. In addition
extremal values are also given in order to construct error bands graphs.

A flowchart of the conduct of one experiment is described in the 
flowchart~\ref{insight_flow}. The number of iteration of the heuristic is 
I and the number of full restart is N. Th function Heuristic() correspond to
a single step of the chosen heuristic. We now present an in-depth view of
the implemented heuristics.

\begin{figure}
 \begin{mdframed}
\includegraphics[scale=0.4]{insight.pdf}
\caption{Flowchart of a single experiment}
\label{insight_flow}
\end{mdframed}
\end{figure}

Graph are presented not with the usual "mustache boxes" to show the 
error bounds, but in a more graphical way with error bands. The graph
of the mean result is included inside a band of the same color which
represents the incertitude with regards to the values obtained.

\section{Heuristics developed}

\subsection{Fully random search (Test case)}
 The first heuristic implemented is the random search. We generates
 random sets of Halton points and select the best set with regard to its
 discrepancy iteratively. The process is wrapped up in the 
 flowchart~\ref{random_flow}. In order to generate at each step a random 
 permutation, we transform it directly from the previous one.
  More precisely the permutation is a singleton object which have method 
  random, built on the Knuth Fisher Yates shuffle. This algorithm allows
  us to generate an uniformly chosen  permutation at each step. We recall 
  this fact and detail the algorithm in the following section.
\begin{figure}
 \begin{mdframed}
\includegraphics[scale=0.4]{flow_rand.pdf}
\caption{Flowchart of the random search}
  \label{random_flow}
\end{mdframed}
\end{figure}


  \subsubsection{The Knuth-Fisher-Yates shuffle}

The Fisher–Yates shuffle is an algorithm for generating a random permutation 
of a finite sets. The Fisher–Yates shuffle is unbiased, so that every 
permutation is equally likely. We present here the Durstenfeld variant of 
the algorithm, presented by Knuth in \emph{The Art of Computer programming}
vol. 2~\cite{Knuth}.
The algorithm's time complexity is here $O(n)$, compared to $O(n^2)$ of 
the naive implementation.

\begin{algorithm}[H]
  \SetAlgoLined
  \SetKwFunction{Rand}{Rand}
  \SetKwFunction{Swap}{Swap}
  \KwData{A table T[1..n]}
  \KwResult{Same table T, shuffled}
  \For{$i\leftarrow 1$ \KwTo $n-1$}{
     $j \leftarrow$ \Rand{$[1,n-i]$}\;
     \Swap{$T[i], T[i+j]$}\;
    }
  \caption{KFY algorithm}
\end{algorithm}


\begin{lemma}
  The resulting permutation of KFY is unbiased.
\end{lemma}
\begin{proof}
  Let consider the set $[1,\ldots n]$ as the vertices of a random graph 
  constructed as the trace of the execution of the algorithm: 
  an edge $(i,j)$ exists in the graph if and only if the swap of $T[i]$ and
  $T[j]$ had been executed. This graph encodes the permutation represented by
  $T$. To be able to encode any permutation the considered graph must be 
  connected --- in order to allow any pairs of points to be swapped ---.
  Since by construction every points is reached by an edge, and that there 
  exists exactly $n-1$ edges, we can conclude directly that any permutation can
  be reached by the algorithm. Since the probability of getting a fixed graph 
  of $n-1$ edges with every edges of degree at least one is $n!^{-1}$, the 
  algorithm is thus unbiased.

\end{proof}


\subsubsection{Results and stability}
We first want to analyze the dependence of the results on the number of 
iterations of the heuristic, in order to discuss its stability. 
The results are compiled in the figures~\ref{rand_iter2},~\ref{rand_iter3},
restricted to a number of points between 80 and 180.
We emphasize on the fact the lots of datas appears on the graphs, 
and the error bands representation make them a bit messy. These graphs
were made for extensive internal experiments and parameters researches.
The final wrap up graphs are much more lighter and only presents the best 
results obtained.
As expected from a fully random search, the error bands are very large for 
low number of iterations ($15\%$ of the value for 400 iterations) and tends
to shrink with a bigger number of iterations (around $5\%$ for 1600 iterations).
This shrinkage is a direct consequence of well known concentrations bounds
(Chernoff and Asuma-Hoeffding).
The average results are quite stable, they decrease progressively with 
the growing number of iterations, but seems to get to a limits after 1000 
iterations. This value acts as a threshold for the interesting number of iterations.
As such interesting results can be conducted with \emph{only} 1000 iterations, 
without altering too much the quality of the set with regards to its
discrepancy and this heuristic.

\begin{figure}
\includegraphics[scale=0.3]{Results/random_iter.png}
\caption{Dependence on iterations, dimension 2}
\label{rand_iter2}
\end{figure}
\begin{figure}
\includegraphics[scale=0.3]{Results/random_iter_3.png}
\caption{Dependence on iterations, dimension 3}
\label{rand_iter3}
\end{figure}

\subsection{Evolutionary heuristic: Simulated annealing and local search}
The second heuristic implemented is a randomized local search with 
simulated annealing. This heuristic is inspired by the physical 
process of annealing in metallurgy.
Simulated annealing interprets the physical slow cooling as a 
slow decrease in the probability of accepting worse solutions as it 
explores the solution space. 
More precisely the neighbours are here the permutations which can be obtained
by application of exactly one transposition of the current permutation.
The selection phase is dependant on the current temperature:
after applying a random transposition on the current permutation, either
the discrepancy of the corresponding Halton set is decreased and the 
evolution is kept, either it does not but is still kept with 
a probability $e^{\frac{\delta}{T}}$ where $\delta$ is the difference
between the old and new discrepancy, and $T$ the current temperature.
The all algorithm is described in the flowchart~\ref{flow_rec}.

\begin{figure}
 \begin{mdframed}
\includegraphics[scale=0.4]{flow_recuit.pdf}
\caption{Flowchart of the simulated annealing local search heuristic}
\label{flow_rec}
\end{mdframed}
\end{figure}

\subsubsection{Dependence on the temperature}
First experiments were made to select the best initial temperature.
Results are compiled in graphs~\ref{temp_2},~\ref{temp3},\ref{temp3_z}.
Graphs~\ref{temp_2},~\ref{temp3} represents the results obtained respectively
in dimension 2 and 3 between 10 and 500 points. The curve obtained is 
characteristic of the average evolution of the discrepancy optimization 
algorithms for Halton points sets: a very fast decrease for low number of 
points --- roughly up to 80 points --- and then a very slow one 
after~\cite{Doerr}.
The most interesting part of these results are concentrated between 80 and 160
points were the different curves splits. The graph~\ref{temp3_z} is a zoom 
of~\ref{temp3} in this window. We remark on that graph that the lower the 
temperature is, the best the results are.

\begin{figure}
\includegraphics[scale=0.3]{Results/resu_2_temp.png}
\caption{Dependence on initial temperature: D=2}
  \label{temp_2}
\end{figure}

\begin{figure}
\includegraphics[scale=0.3]{Results/resu_temp3.png}
\caption{Dependence on initial temperature: D=3}
  \label{temp3}
\end{figure}

\begin{figure}
\includegraphics[scale=0.3]{Results/resu_temp3_zoom.png}
\caption{Dependence on initial temperature (zoom): D=3}
  \label{temp3_z}
\end{figure}


\subsubsection{Stability with regards to the number of iterations}

As for the fully random search heuristic we investigated the stability
of the algorithm with regards to the number of iterations. We present here
the result in dimension 3 in the graph~\ref{iter_sa}. Once again we
restricted the window between 80 and 180 points were curves are split.
An interesting phenomena can be observed: the error rates are somehow 
invariant w.r.t.\ the number of iteration and once again the 1000 iterations
threshold seems to appear --- point 145 is a light split between iteration 
1600 and the others, but excepted for that point, getting more than 1000
iterations tends be be a waste of time. The error rate is for 80 points the
biggest and is about $15\%$ of the value, which is similar to the error
rates for fully random search with 400 iterations.

\begin{figure}
\includegraphics[scale=0.3]{Results/sa_iter.png}
\caption{Dependence on iterations number for simulated annealing : D=3}
  \label{iter_sa}
\end{figure}

\subsection{Genetic (5+5) search}


\begin{figure}
 \begin{mdframed}
  \label{rand_flow}
\includegraphics[scale=0.4]{crossover_flow.pdf}
\caption{Flowchart of the crossover algorithm.}
\end{mdframed}
\end{figure}

\begin{figure}
  \label{rand_flow}
\includegraphics[scale=0.3]{Results/res_gen_2.png}
\caption{Dependence on iterations number: D=3}
\end{figure}

\begin{figure}
  \label{rand_flow}
\includegraphics[scale=0.3]{Results/res_gen_2_zoom.png}
\caption{Dependence on iterations number: D=3}
\end{figure}
\begin{figure}
  \label{rand_flow}
\includegraphics[scale=0.3]{Results/res_gen_3_zoom.png}
\caption{Dependence on iterations number: D=3}
\end{figure}

\begin{figure}
  \label{rand_flow}
\includegraphics[scale=0.3]{Results/res_gen_4_zoom.png}
\caption{Dependence on iterations number: D=3}
\end{figure}

\section{Results}


\section{Conclusion}
\bibliographystyle{alpha}
\bibliography{bi}
\end{document}