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.
108 \subsection{A given system can always be claimed as non-chaotic}
111 Let $\mathcal{X}$ be a set and $f: \mathcal{X} \to \X$.
112 If $\X$ is infinite, then $\left( \X_{\tau_\infty}, f\right)$ is not chaotic (for the Devaney's formulation), where $\tau_\infty$ is the discrete topology.
116 Let us prove it by contradiction, assuming that $\left(\X_{\tau_\infty}, f\right)$ is both transitive and regular.
118 Let $x \in \X$ and $\{x\}$ one of its neighborhood. This neighborhood must contain a periodic point for $f$, if we want that $\left(\X_{\tau_\infty}, f\right)$ is regular. Then $x$ must be a periodic point of $f$.
120 Let $I_x = \left\{ f^{(n)}(x), n \in \mathds{N}\right\}$. This set is finite because $x$ is periodic, and $\mathcal{X}$ is infinite, then $\exists y \in \mathcal{X}, y \notin I_x$.
122 As $\left(\X_{\tau_\infty}, f\right)$ must be transitive, for all open nonempty sets $A$ and $B$, an integer $n$ must satisfy $f^{(n)}(A) \cap B \neq \varnothing$. However $\{x\}$ and $\{y\}$ are open sets and $y \notin I_x \Rightarrow \forall n, f^{(n)}\left( \{x\} \right) \cap \{y\} = \varnothing$.
130 \section{Chaos on the order topology}
132 \subsection{The phase space is an interval of the real line}
134 \subsubsection{Toward a topological semiconjugacy}
136 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.
138 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}$.
142 The function $\varphi: \mathcal{S}_{10} \times\mathds{B}^{10} \rightarrow \big[ 0, 2^{10} \big[$ is defined by:
145 \varphi: & \mathcal{X}_{10} = \mathcal{S}_{10} \times\mathds{B}^{10}& \longrightarrow & \big[ 0, 2^{10} \big[ \\
146 & (S,E) = \left((S^0, S^1, \hdots ); (E_0, \hdots, E_9)\right) & \longmapsto & \varphi \left((S,E)\right)
149 \noindent where $\varphi\left((S,E)\right)$ is the real number:
151 \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$.
152 \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}.$
158 $\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:
163 Let $x \in \big[ 0, 2^{10} \big[$ and:
165 \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}$.
166 \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:
167 $\displaystyle{x = \lfloor x \rfloor + \sum_{k=0}^{+\infty} s^k 10^{-k-1}}$.
169 $e$ and $s$ are thus defined as follows:
172 e: & \big[ 0, 2^{10} \big[ & \longrightarrow & \mathds{B}^{10} \\
173 & x & \longmapsto & (e_0, \hdots, e_9)
179 s: & \big[ 0, 2^{10} \big[ & \longrightarrow & \llbracket 0, 9 \rrbracket^{\mathds{N}} \\
180 & x & \longmapsto & (s^k)_{k \in \mathds{N}}
185 We are now able to define the function $g$, whose goal is to translate the chaotic iterations $\Go$ on an interval of $\mathds{R}$.
188 $g:\big[ 0, 2^{10} \big[ \longrightarrow \big[ 0, 2^{10} \big[$ is defined by:
191 g: & \big[ 0, 2^{10} \big[ & \longrightarrow & \big[ 0, 2^{10} \big[ \\
193 & x & \longmapsto & g(x)
196 \noindent where g(x) is the real number of $\big[ 0, 2^{10} \big[$ defined bellow:
198 \item its integral part has a binary decomposition equal to $e_0', \hdots, e_9'$, with:
202 e(x)_i & \textrm{ if } i \neq s^0\\
203 e(x)_i + 1 \textrm{ (mod 2)} & \textrm{ if } i = s^0\\
207 \item whose decimal part is $s(x)^1, s(x)^2, \hdots$
214 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}}.$$
216 \subsubsection{Defining a metric on $\big[ 0, 2^{10} \big[$}
218 Numerous metrics can be defined on the set $\big[ 0, 2^{10} \big[$, the most usual one being the Euclidian distance recalled bellow:
221 \index{distance!euclidienne}
222 $\Delta$ is the Euclidian distance on $\big[ 0, 2^{10} \big[$, that is, $\Delta(x,y) = |y-x|^2$.
227 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:
232 Let $x,y \in \big[ 0, 2^{10} \big[$.
233 $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:
235 $\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}}$.
240 $D$ is a distance on $\big[ 0, 2^{10} \big[$.
244 The three axioms defining a distance must be checked.
246 \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$.
247 \item $D(x,y)=D(y,x)$.
248 \item Finally, the triangular inequality is obtained due to the fact that both $\delta$ and $\Delta(x,y)=|x-y|$ satisfy it.
253 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''.
254 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.
259 \subfigure[Function $x \to dist(x;1,234) $ on the interval $(0;5)$.]{\includegraphics[scale=.35]{DvsEuclidien.pdf}}\quad
260 \subfigure[Function $x \to dist(x;3) $ on the interval $(0;5)$.]{\includegraphics[scale=.35]{DvsEuclidien2.pdf}}
262 \caption{Comparison between $D$ (in blue) and the Euclidian distane (in green).}
263 \label{fig:comparaison de distances}
269 \subsubsection{The semiconjugacy}
271 It is now possible to define a topological semiconjugacy between $\mathcal{X}$ and an interval of $\mathds{R}$:
274 Chaotic iterations on the phase space $\mathcal{X}$ are simple iterations on $\mathds{R}$, which is illustrated by the semiconjugacy of the diagram bellow:
277 \left(~\mathcal{S}_{10} \times\mathds{B}^{10}, d~\right) @>G_{f_0}>> \left(~\mathcal{S}_{10} \times\mathds{B}^{10}, d~\right)\\
278 @V{\varphi}VV @VV{\varphi}V\\
279 \left( ~\big[ 0, 2^{10} \big[, D~\right) @>>g> \left(~\big[ 0, 2^{10} \big[, D~\right)
285 $\varphi$ has been constructed in order to be continuous and onto.
288 In other words, $\mathcal{X}$ is approximately equal to $\big[ 0, 2^\mathsf{N} \big[$.
295 \subsection{Study of the chaotic iterations described as a real function}
300 \subfigure[ICs on the interval $(0,9;1)$.]{\includegraphics[scale=.35]{ICs09a1.pdf}}\quad
301 \subfigure[ICs on the interval $(0,7;1)$.]{\includegraphics[scale=.35]{ICs07a95.pdf}}\\
302 \subfigure[ICs on the interval $(0,5;1)$.]{\includegraphics[scale=.35]{ICs05a1.pdf}}\quad
303 \subfigure[ICs on the interval $(0;1)$]{\includegraphics[scale=.35]{ICs0a1.pdf}}
305 \caption{Representation of the chaotic iterations.}
314 \subfigure[ICs on the interval $(510;514)$.]{\includegraphics[scale=.35]{ICs510a514.pdf}}\quad
315 \subfigure[ICs on the interval $(1000;1008)$]{\includegraphics[scale=.35]{ICs1000a1008.pdf}}
317 \caption{ICs on small intervals.}
323 \subfigure[ICs on the interval $(0;16)$.]{\includegraphics[scale=.3]{ICs0a16.pdf}}\quad
324 \subfigure[ICs on the interval $(40;70)$.]{\includegraphics[scale=.45]{ICs40a70.pdf}}\quad
326 \caption{General aspect of the chaotic iterations.}
331 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:
334 \label{Prop:derivabilite des ICs}
335 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\}$.
337 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$.
342 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$:
344 \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}$).
345 \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$.
347 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$.
351 Finally, chaotic iterations are elements of the large family of functions that are both chaotic and piecewise linear (like the tent map).
356 \subsection{Comparison of the two metrics on $\big[ 0, 2^\mathsf{N} \big[$}
358 The two propositions bellow allow to compare our two distances on $\big[ 0, 2^\mathsf{N} \big[$:
361 Id: $\left(~\big[ 0, 2^\mathsf{N} \big[,\Delta~\right) \to \left(~\big[ 0, 2^\mathsf{N} \big[, D~\right)$ is not continuous.
365 The sequence $x^n = 1,999\hdots 999$ constituted by $n$ 9 as decimal part, is such that:
367 \item $\Delta (x^n,2) \to 0.$
368 \item But $D(x^n,2) \geqslant 1$, then $D(x^n,2)$ does not converge to 0.
371 The sequential characterization of the continuity concludes the demonstration.
379 Id: $\left(~\big[ 0, 2^\mathsf{N} \big[,D~\right) \to \left(~\big[ 0, 2^\mathsf{N} \big[, \Delta ~\right)$ is a continuous fonction.
383 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$.
385 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.
388 The conclusion of these propositions is that the proposed metric is more precise than the Euclidian distance, that is:
391 $D$ is finer than the Euclidian distance $\Delta$.
394 This corollary can be reformulated as follows:
397 \item The topology produced by $\Delta$ is a subset of the topology produced by $D$.
398 \item $D$ has more open sets than $\Delta$.
399 \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$.
403 \subsection{Chaos of the chaotic iterations on $\mathds{R}$}
404 \label{chpt:Chaos des itérations chaotiques sur R}
408 \subsubsection{Chaos according to Devaney}
410 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:
412 \item $\left(\Go, \mathcal{X}_d\right)$ and $\left(g, \big[ 0, 2^{10} \big[_D\right)$ are semiconjugate by $\varphi$,
413 \item Then $\left(g, \big[ 0, 2^{10} \big[_D\right)$ is a system chaotic according to Devaney, because the semiconjugacy preserve this character.
414 \item But the topology generated by $D$ is finer than the topology generated by the Euclidian distance $\Delta$ -- which is the order topology.
415 \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}$.
418 This result can be formulated as follows.
421 \label{th:IC et topologie de l'ordre}
422 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.
425 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.
430 \section{Efficient prng based on chaotic iterations}
432 On parle du séquentiel avec des nombres 64 bits\\
434 Faire le lien avec le paragraphe précédent (je considère que la stratégie s'appelle $S^i$\\
436 In order to implement efficiently a PRNG based on chaotic iterations it is
437 possible to improve previous works [ref]. One solution consists in considering
438 that the strategy used $S^i$ contains all the bits for which the negation is
439 achieved out. Then instead of applying the negation on these bits we can simply
440 apply the xor operator between the current number and the strategy $S^i$. In
441 order to obtain the strategy we also use a classical PRNG.
446 \begin{minipage}{14cm}
447 unsigned int CIprng() \{\\
448 static unsigned int x = 123123123;\\
449 unsigned long t1 = xorshift();\\
450 unsigned long t2 = xor128();\\
451 unsigned long t3 = xorwow();\\
452 x = x\textasciicircum (unsigned int)t1;\\
453 x = x\textasciicircum (unsigned int)(t2$>>$32);\\
454 x = x\textasciicircum (unsigned int)(t3$>>$32);\\
455 x = x\textasciicircum (unsigned int)t2;\\
456 x = x\textasciicircum (unsigned int)(t1$>>$32);\\
457 x = x\textasciicircum (unsigned int)t3;\\
463 \caption{sequential Chaotic Iteration PRNG}
464 \label{algo:seqCIprng}
467 In Figure~\ref{algo:seqCIprng} a sequential version of our chaotic iterations
468 based PRNG is presented. This version uses three classical 64 bits PRNG: the
469 \texttt{xorshift}, the \texttt{xor128} and the \texttt{xorwow}. These three
470 PRNGs are presented in~\cite{Marsaglia2003}. As each PRNG used works with
471 64-bits and as our PRNG works with 32 bits, the use of \texttt{(unsigned int)}
472 selects the 32 least significant bits whereas \texttt{(unsigned int)(t3$>>$32)}
473 selects the 32 most significants bits of the variable \texttt{t}. This version
474 sucesses the BigCrush of the TestU01 battery [P. L’ecuyer and
477 \section{Efficient prng based on chaotic iterations on GPU}
479 On parle du passage du sequentiel au GPU
481 \section{Experiments}
483 On passe le BigCrush\\
484 On donne des temps de générations sur GPU/CPU\\
485 On donne des temps de générations de nombre sur GPU puis on rappatrie sur CPU / CPU ? bof bof, on verra
489 \bibliographystyle{plain}
490 \bibliography{mabase}