1 \documentclass{article}
2 \usepackage[utf8]{inputenc}
3 \usepackage[T1]{fontenc}
10 \usepackage[standard]{ntheorem}
12 % Pour mathds : les ensembles IR, IN, etc.
15 % Pour avoir des intervalles d'entiers
19 % Pour faire des sous-figures dans les figures
20 \usepackage{subfigure}
24 \newtheorem{notation}{Notation}
26 \newcommand{\X}{\mathcal{X}}
27 \newcommand{\Go}{G_{f_0}}
28 \newcommand{\B}{\mathds{B}}
29 \newcommand{\N}{\mathds{N}}
30 \newcommand{\BN}{\mathds{B}^\mathsf{N}}
33 \newcommand{\alert}[1]{\begin{color}{blue}\textit{#1}\end{color}}
35 \title{Efficient generation of pseudo random numbers based on chaotic iterations on GPU}
43 \section{Introduction}
45 Interet des itérations chaotiques pour générer des nombre alea\\
46 Interet de générer des nombres alea sur GPU
49 \section{Chaotic iterations}
51 Présentation des itérations chaotiques
56 \section{The relativity of disorder}
57 \label{sec:de la relativité du désordre}
59 \subsection{Impact of the topology's finenesse}
61 Let us firstly introduce the following notations.
64 $\mathcal{X}_\tau$ will denote the topological space $\left(\mathcal{X},\tau\right)$, whereas $\mathcal{V}_\tau (x)$ will be the set of all the neighborhoods of $x$ when considering the topology $\tau$ (or simply $\mathcal{V} (x)$, if there is no ambiguity).
70 \label{Th:chaos et finesse}
71 Let $\mathcal{X}$ a set and $\tau, \tau'$ two topologies on $\mathcal{X}$ s.t. $\tau'$ is finer than $\tau$. Let $f:\mathcal{X} \to \mathcal{X}$, continuous both for $\tau$ and $\tau'$.
73 If $(\mathcal{X}_{\tau'},f)$ is chaotic according to Devaney, then $(\mathcal{X}_\tau,f)$ is chaotic too.
77 Let us firstly establish the transitivity of $(\mathcal{X}_\tau,f)$.
79 Let $\omega_1, \omega_2$ two open sets of $\tau$. Then $\omega_1, \omega_2 \in \tau'$, becaus $\tau'$ is finer than $\tau$. As $f$ is $\tau'-$transitive, we can deduce that $\exists n \in \mathds{N}, \omega_1 \cap f^{(n)}(\omega_2) = \varnothing$. Consequently, $f$ is $\tau-$transitive.
81 Let us now consider the regularity of $(\mathcal{X}_\tau,f)$, \emph{i.e.}, for all $x \in \mathcal{X}$, and for all $\tau-$neighborhood $V$ of $x$, there is a periodic point for $f$ into $V$.
83 Let $x \in \mathcal{X}$ and $V \in \mathcal{V}_\tau (x)$ a $\tau-$neighborhood of $x$. By definition, $\exists \omega \in \tau, x \in \omega \subset V$.
85 But $\tau \subset \tau'$, so $\omega \in \tau'$, and then $V \in \mathcal{V}_{\tau'} (x)$. As $(\mathcal{X}_{\tau'},f)$ is regular, there is a periodic point for $f$ into $V$, and the regularity of $(\mathcal{X}_\tau,f)$ is proven.
88 \subsection{A given system can always be claimed as chaotic}
90 Let $f$ an iteration function on $\mathcal{X}$ having at least a fixed point. Then this function is chaotic (in a certain way):
93 Let $\mathcal{X}$ a nonempty set and $f: \mathcal{X} \to \X$ a function having at least a fixed point.
94 Then $f$ is $\tau_0-$chaotic, where $\tau_0$ is the trivial (indiscrete) topology on $\X$.
99 $f$ is transitive when $\forall \omega, \omega' \in \tau_0 \setminus \{\varnothing\}, \exists n \in \mathds{N}, f^{(n)}(\omega) \cap \omega' \neq \varnothing$.
100 As $\tau_0 = \left\{ \varnothing, \X \right\}$, this is equivalent to look for an integer $n$ s.t. $f^{(n)}\left( \X \right) \cap \X \neq \varnothing$. For instance, $n=0$ is appropriate.
102 Let us now consider $x \in \X$ and $V \in \mathcal{V}_{\tau_0} (x)$. Then $V = \mathcal{X}$, so $V$ has at least a fixed point for $f$. Consequently $f$ is regular, and the result is established.
111 \section{Chaos on the order topology}
113 \subsection{The phase space is an interval of the real line}
115 \subsubsection{Toward a topological semiconjugacy}
117 In what follows, our intention is to establish, by using a topological semiconjugacy, that chaotic iterations over $\mathcal{X}$ can be described as iterations on a real interval. To do so, we must firstly introduce some notations and terminologies.
119 Let $\mathcal{S}_\mathsf{N}$ be the set of sequences belonging into $\llbracket 1; \mathsf{N}\rrbracket$ and $\mathcal{X}_{\mathsf{N}} = \mathcal{S}_\mathsf{N} \times \B^\mathsf{N}$.
123 The function $\varphi: \mathcal{S}_{10} \times\mathds{B}^{10} \rightarrow \big[ 0, 2^{10} \big[$ is defined by:
126 \varphi: & \mathcal{X}_{10} = \mathcal{S}_{10} \times\mathds{B}^{10}& \longrightarrow & \big[ 0, 2^{10} \big[ \\
127 & (S,E) = \left((S^0, S^1, \hdots ); (E_0, \hdots, E_9)\right) & \longmapsto & \varphi \left((S,E)\right)
130 \noindent where $\varphi\left((S,E)\right)$ is the real number:
132 \item whose integral part $e$ is $\displaystyle{\sum_{k=0}^9 2^{9-k} E_k}$, that is, the binary digits of $e$ are $E_0 ~ E_1 ~ \hdots ~ E_9$.
133 \item whose decimal part $s$ is equal to $s = 0,S^0~ S^1~ S^2~ \hdots = \sum_{k=1}^{+\infty} 10^{-k} S^{k-1}.$
139 $\varphi$ realizes the association between a point of $\mathcal{X}_{10}$ and a real number into $\big[ 0, 2^{10} \big[$. We must now translate the chaotic iterations $\Go$ on this real interval. To do so, two intermediate functions over $\big[ 0, 2^{10} \big[$ must be introduced:
144 Let $x \in \big[ 0, 2^{10} \big[$ and:
146 \item $e_0, \hdots, e_9$ the binary digits of the integral part of $x$: $\displaystyle{\lfloor x \rfloor = \sum_{k=0}^{9} 2^{9-k} e_k}$.
147 \item $(s^k)_{k\in \mathds{N}}$ the digits of $x$, where the chosen decimal decomposition of $x$ is the one that does not have an infinite number of 9:
148 $\displaystyle{x = \lfloor x \rfloor + \sum_{k=0}^{+\infty} s^k 10^{-k-1}}$.
150 $e$ and $s$ are thus defined as follows:
153 e: & \big[ 0, 2^{10} \big[ & \longrightarrow & \mathds{B}^{10} \\
154 & x & \longmapsto & (e_0, \hdots, e_9)
160 s: & \big[ 0, 2^{10} \big[ & \longrightarrow & \llbracket 0, 9 \rrbracket^{\mathds{N}} \\
161 & x & \longmapsto & (s^k)_{k \in \mathds{N}}
166 We are now able to define the function $g$, whose goal is to translate the chaotic iterations $\Go$ on an interval of $\mathds{R}$.
169 $g:\big[ 0, 2^{10} \big[ \longrightarrow \big[ 0, 2^{10} \big[$ is defined by:
172 g: & \big[ 0, 2^{10} \big[ & \longrightarrow & \big[ 0, 2^{10} \big[ \\
174 & x & \longmapsto & g(x)
177 \noindent where g(x) is the real number of $\big[ 0, 2^{10} \big[$ defined bellow:
179 \item its integral part has a binary decomposition equal to $e_0', \hdots, e_9'$, with:
183 e(x)_i & \textrm{ if } i \neq s^0\\
184 e(x)_i + 1 \textrm{ (mod 2)} & \textrm{ if } i = s^0\\
188 \item whose decimal part is $s(x)^1, s(x)^2, \hdots$
195 In other words, if $x = \displaystyle{\sum_{k=0}^{9} 2^{9-k} e_k + \sum_{k=0}^{+\infty} s^{k} ~10^{-k-1}}$, then: $$g(x) = \displaystyle{\sum_{k=0}^{9} 2^{9-k} (e_k + \delta(k,s^0) \textrm{ (mod 2)}) + \sum_{k=0}^{+\infty} s^{k+1} 10^{-k-1}}.$$
197 \subsubsection{Defining a metric on $\big[ 0, 2^{10} \big[$}
199 Numerous metrics can be defined on the set $\big[ 0, 2^{10} \big[$, the most usual one being the Euclidian distance recalled bellow:
202 \index{distance!euclidienne}
203 $\Delta$ is the Euclidian distance on $\big[ 0, 2^{10} \big[$, that is, $\Delta(x,y) = |y-x|^2$.
208 This Euclidian distance does not reproduce exactly the notion of proximity induced by our first distance $d$ on $\X$. Indeed $d$ is finer than $\Delta$. This is the reason why we have to introduce the following metric:
213 Let $x,y \in \big[ 0, 2^{10} \big[$.
214 $D$ denotes the function from $\big[ 0, 2^{10} \big[^2$ to $\mathds{R}^+$ defined by: $D(x,y) = D_e\left(e(x),e(y)\right) + D_s\left(s(x),s(y)\right)$, where:
216 $\displaystyle{D_e(E,\check{E}) = \sum_{k=0}^\mathsf{9} \delta (E_k, \check{E}_k)}$, ~~and~ $\displaystyle{D_s(S,\check{S}) = \sum_{k = 1}^\infty \dfrac{|S^k-\check{S}^k|}{10^k}}$.
221 $D$ is a distance on $\big[ 0, 2^{10} \big[$.
225 The three axioms defining a distance must be checked.
227 \item $D \geqslant 0$, because everything is positive in its definition. If $D(x,y)=0$, then $D_e(x,y)=0$, so the integral parts of $x$ and $y$ are equal (they have the same binary decomposition). Additionally, $D_s(x,y) = 0$, then $\forall k \in \mathds{N}^*, s(x)^k = s(y)^k$. In other words, $x$ and $y$ have the same $k-$th decimal digit, $\forall k \in \mathds{N}^*$. And so $x=y$.
228 \item $D(x,y)=D(y,x)$.
229 \item Finally, the triangular inequality is obtained due to the fact that both $\delta$ and $\Delta(x,y)=|x-y|$ satisfy it.
234 The convergence of sequences according to $D$ is not the same than the usual convergence related to the Euclidian metric. For instance, if $x^n \to x$ according to $D$, then necessarily the integral part of each $x^n$ is equal to the integral part of $x$ (at least after a given threshold), and the decimal part of $x^n$ corresponds to the one of $x$ ``as far as required''.
235 To illustrate this fact, a comparison between $D$ and the Euclidian distance is given Figure \ref{fig:comparaison de distances}. These illustrations show that $D$ is richer and more refined than the Euclidian distance, and thus is more precise.
240 \subfigure[Function $x \to dist(x;1,234) $ on the interval $(0;5)$.]{\includegraphics[scale=.35]{DvsEuclidien.pdf}}\quad
241 \subfigure[Function $x \to dist(x;3) $ on the interval $(0;5)$.]{\includegraphics[scale=.35]{DvsEuclidien2.pdf}}
243 \caption{Comparison between $D$ (in blue) and the Euclidian distane (in green).}
244 \label{fig:comparaison de distances}
250 \subsubsection{The semiconjugacy}
252 It is now possible to define a topological semiconjugacy between $\mathcal{X}$ and an interval of $\mathds{R}$:
255 Chaotic iterations on the phase space $\mathcal{X}$ are simple iterations on $\mathds{R}$, which is illustrated by the semiconjugacy of the diagram bellow:
258 \left(~\mathcal{S}_{10} \times\mathds{B}^{10}, d~\right) @>G_{f_0}>> \left(~\mathcal{S}_{10} \times\mathds{B}^{10}, d~\right)\\
259 @V{\varphi}VV @VV{\varphi}V\\
260 \left( ~\big[ 0, 2^{10} \big[, D~\right) @>>g> \left(~\big[ 0, 2^{10} \big[, D~\right)
266 $\varphi$ has been constructed in order to be continuous and onto.
269 In other words, $\mathcal{X}$ is approximately equal to $\big[ 0, 2^\mathsf{N} \big[$.
276 \subsection{Study of the chaotic iterations described as a real function}
281 \subfigure[ICs on the interval $(0,9;1)$.]{\includegraphics[scale=.35]{ICs09a1.pdf}}\quad
282 \subfigure[ICs on the interval $(0,7;1)$.]{\includegraphics[scale=.35]{ICs07a95.pdf}}\\
283 \subfigure[ICs on the interval $(0,5;1)$.]{\includegraphics[scale=.35]{ICs05a1.pdf}}\quad
284 \subfigure[ICs on the interval $(0;1)$]{\includegraphics[scale=.35]{ICs0a1.pdf}}
286 \caption{Representation of the chaotic iterations.}
295 \subfigure[ICs on the interval $(510;514)$.]{\includegraphics[scale=.35]{ICs510a514.pdf}}\quad
296 \subfigure[ICs on the interval $(1000;1008)$]{\includegraphics[scale=.35]{ICs1000a1008.pdf}}
298 \caption{ICs on small intervals.}
304 \subfigure[ICs on the interval $(0;16)$.]{\includegraphics[scale=.3]{ICs0a16.pdf}}\quad
305 \subfigure[ICs on the interval $(40;70)$.]{\includegraphics[scale=.45]{ICs40a70.pdf}}\quad
307 \caption{General aspect of the chaotic iterations.}
312 We have written a Python program to represent the chaotic iterations with the vectorial negation on the real line $\mathds{R}$. Various representations of these CIs are given in Figures \ref{fig:ICs}, \ref{fig:ICs2} and \ref{fig:ICs3}. It can be remarked that the function $g$ is a piecewise linear function: it is linear on each interval having the form $\left[ \dfrac{n}{10}, \dfrac{n+1}{10}\right[$, $n \in \llbracket 0;2^{10}\times 10 \rrbracket$ and its slope is equal to 10. Let us justify these claims:
315 \label{Prop:derivabilite des ICs}
316 Chaotic iterations $g$ defined on $\mathds{R}$ have derivatives of all orders on $\big[ 0, 2^{10} \big[$, except on the 10241 points in $I$ defined by $\left\{ \dfrac{n}{10} ~\big/~ n \in \llbracket 0;2^{10}\times 10\rrbracket \right\}$.
318 Furthermore, on each interval of the form $\left[ \dfrac{n}{10}, \dfrac{n+1}{10}\right[$, with $n \in \llbracket 0;2^{10}\times 10 \rrbracket$, $g$ is a linear function, having a slope equal to 10: $\forall x \notin I, g'(x)=10$.
323 Let $I_n = \left[ \dfrac{n}{10}, \dfrac{n+1}{10}\right[$, with $n \in \llbracket 0;2^{10}\times 10 \rrbracket$. All the points of $I_n$ have the same integral prat $e$ and the same decimal part $s^0$: on the set $I_n$, functions $e(x)$ and $x \mapsto s(x)^0$ of Definition \ref{def:e et s} only depend on $n$. So all the images $g(x)$ of these points $x$:
325 \item Have the same integral part, which is $e$, except probably the bit number $s^0$. In other words, this integer has approximately the same binary decomposition than $e$, the sole exception being the digit $s^0$ (this number is then either $e+2^{10-s^0}$ or $e-2^{10-s^0}$, depending on the parity of $s^0$, \emph{i.e.}, it is equal to $e+(-1)^{s^0}\times 2^{10-s^0}$).
326 \item A shift to the left has been applied to the decimal part $y$, losing by doing so the common first digit $s^0$. In other words, $y$ has been mapped into $10\times y - s^0$.
328 To sum up, the action of $g$ on the points of $I$ is as follows: first, make a multiplication by 10, and second, add the same constant to each term, which is $\dfrac{1}{10}\left(e+(-1)^{s^0}\times 2^{10-s^0}\right)-s^0$.
332 Finally, chaotic iterations are elements of the large family of functions that are both chaotic and piecewise linear (like the tent map).
337 \subsection{Comparison of the two metrics on $\big[ 0, 2^\mathsf{N} \big[$}
339 The two propositions bellow allow to compare our two distances on $\big[ 0, 2^\mathsf{N} \big[$:
342 Id: $\left(~\big[ 0, 2^\mathsf{N} \big[,\Delta~\right) \to \left(~\big[ 0, 2^\mathsf{N} \big[, D~\right)$ is not continuous.
346 The sequence $x^n = 1,999\hdots 999$ constituted by $n$ 9 as decimal part, is such that:
348 \item $\Delta (x^n,2) \to 0.$
349 \item But $D(x^n,2) \geqslant 1$, then $D(x^n,2)$ does not converge to 0.
352 The sequential characterization of the continuity concludes the demonstration.
360 Id: $\left(~\big[ 0, 2^\mathsf{N} \big[,D~\right) \to \left(~\big[ 0, 2^\mathsf{N} \big[, \Delta ~\right)$ is a continuous fonction.
364 If $D(x^n,x) \to 0$, then $D_e(x^n,x) = 0$ at least for $n$ larger than a given threshold, because $D_e$ only returns integers. So, after this threshold, the integral parts of all the $x^n$ are equal to the integral part of $x$.
366 Additionally, $D_s(x^n, x) \to 0$, then $\forall k \in \mathds{N}^*, \exists N_k \in \mathds{N}, n \geqslant N_k \Rightarrow D_s(x^n,x) \leqslant 10^{-k}$. This means that for all $k$, an index $N_k$ can be found such that, $\forall n \geqslant N_k$, all the $x^n$ have the same $k$ firsts digits, which are the digits of $x$. We can deduce the convergence $\Delta(x^n,x) \to 0$, and thus the result.
369 The conclusion of these propositions is that the proposed metric is more precise than the Euclidian distance, that is:
372 $D$ is finer than the Euclidian distance $\Delta$.
375 This corollary can be reformulated as follows:
378 \item The topology produced by $\Delta$ is a subset of the topology produced by $D$.
379 \item $D$ has more open sets than $\Delta$.
380 \item It is harder to converge for the topology $\tau_D$ inherited by $D$, than to converge with the one inherited by $\Delta$, which is denoted here by $\tau_\Delta$.
384 \subsection{Chaos of the chaotic iterations on $\mathds{R}$}
385 \label{chpt:Chaos des itérations chaotiques sur R}
389 \subsubsection{Chaos according to Devaney}
391 We have recalled previously that the chaotic iterations $\left(\Go, \mathcal{X}_d\right)$ are chaotic according to the formulation of Devaney. We can deduce that they are chaotic on $\mathds{R}$ too, when considering the order topology, because:
393 \item $\left(\Go, \mathcal{X}_d\right)$ and $\left(g, \big[ 0, 2^{10} \big[_D\right)$ are semiconjugate by $\varphi$,
394 \item Then $\left(g, \big[ 0, 2^{10} \big[_D\right)$ is a system chaotic according to Devaney, because the semiconjugacy preserve this character.
395 \item But the topology generated by $D$ is finer than the topology generated by the Euclidian distance $\Delta$ -- which is the order topology.
396 \item According to Theorem \ref{Th:chaos et finesse}, we can deduce that the chaotic iterations $g$ are indeed chaotic, as defined by Devaney, for the order topology on $\mathds{R}$.
399 This result can be formulated as follows.
402 \label{th:IC et topologie de l'ordre}
403 The chaotic iterations $g$ on $\mathds{R}$ are chaotic according to the Devaney's formulation, when $\mathds{R}$ has his usual topology, which is the order topology.
406 Indeed this result is weaker than the theorem establishing the chaos for the finer topology $d$. However the Theorem \ref{th:IC et topologie de l'ordre} still remains important. Indeed, we have studied in our previous works a set different from the usual set of study ($\mathcal{X}$ instead of $\mathds{R}$), in order to be as close as possible from the computer: the properties of disorder proved theoretically will then be preserved when computing. However, we could wonder whether this change does not lead to a disorder of a lower quality. In other words, have we replaced a situation of a good disorder lost when computing, to another situation of a disorder preserved but of bad quality. Theorem \ref{th:IC et topologie de l'ordre} prove exactly the contrary.
411 \section{Efficient prng based on chaotic iterations}
413 On parle du séquentiel avec des nombres 64 bits\\
415 Faire le lien avec le paragraphe précédent (je considère que la stratégie s'appelle $S^i$\\
417 In order to implement efficiently a PRNG based on chaotic iterations it is
418 possible to improve previous works [ref]. One solution consists in considering
419 that the strategy used $S^i$ contains all the bits for which the negation is
420 achieved out. Then instead of applying the negation on these bits we can simply
421 apply the xor operator between the current number and the strategy $S^i$. In
422 order to obtain the strategy we also use a classical PRNG.
427 \begin{minipage}{14cm}
428 unsigned int CIprng() \{\\
429 static unsigned int x = 123123123;\\
430 unsigned long t1 = xorshift();\\
431 unsigned long t2 = xor128();\\
432 unsigned long t3 = xorwow();\\
433 x = x\textasciicircum (unsigned int)t1;\\
434 x = x\textasciicircum (unsigned int)(t2$>>$32);\\
435 x = x\textasciicircum (unsigned int)(t3$>>$32);\\
436 x = x\textasciicircum (unsigned int)t2;\\
437 x = x\textasciicircum (unsigned int)(t1$>>$32);\\
438 x = x\textasciicircum (unsigned int)t3;\\
444 \caption{sequential Chaotic Iteration PRNG}
445 \label{algo:seqCIprng}
448 In Figure~\ref{algo:seqCIprng} a sequential version of our chaotic iterations
449 based PRNG is presented. This version uses three classical 64 bits PRNG: the
450 \texttt{xorshift}, the \texttt{xor128} and the \texttt{xorwow}. These three
451 PRNGs are presented in~\cite{Marsaglia2003}. As each PRNG used works with
452 64-bits and as our PRNG works with 32 bits, the use of \texttt{(unsigned int)}
453 selects the 32 least significant bits whereas \texttt{(unsigned int)(t3$>>$32)}
454 selects the 32 most significants bits of the variable \texttt{t}. This version
455 sucesses the BigCrush of the TestU01 battery [P. L’ecuyer and
458 \section{Efficient prng based on chaotic iterations on GPU}
460 On parle du passage du sequentiel au GPU
462 \section{Experiments}
464 On passe le BigCrush\\
465 On donne des temps de générations sur GPU/CPU\\
466 On donne des temps de générations de nombre sur GPU puis on rappatrie sur CPU / CPU ? bof bof, on verra
470 \bibliographystyle{plain}
471 \bibliography{mabase}