X-Git-Url: https://bilbo.iut-bm.univ-fcomte.fr/and/gitweb/prng_gpu.git/blobdiff_plain/77973c61177ccd9f68459a2095c29d48d11af26f..4bdc449d7d1ab2508a4f72db158e0130f61fe2ef:/prng_gpu.tex diff --git a/prng_gpu.tex b/prng_gpu.tex index e177e24..73bc958 100644 --- a/prng_gpu.tex +++ b/prng_gpu.tex @@ -4,8 +4,33 @@ \usepackage{fullpage} \usepackage{fancybox} \usepackage{amsmath} +\usepackage{amscd} \usepackage{moreverb} \usepackage{commath} +\usepackage[standard]{ntheorem} + +% Pour mathds : les ensembles IR, IN, etc. +\usepackage{dsfont} + +% Pour avoir des intervalles d'entiers +\usepackage{stmaryrd} + +\usepackage{graphicx} +% Pour faire des sous-figures dans les figures +\usepackage{subfigure} + +\usepackage{color} + +\newtheorem{notation}{Notation} + +\newcommand{\X}{\mathcal{X}} +\newcommand{\Go}{G_{f_0}} +\newcommand{\B}{\mathds{B}} +\newcommand{\N}{\mathds{N}} +\newcommand{\BN}{\mathds{B}^\mathsf{N}} +\let\sur=\overline + +\newcommand{\alert}[1]{\begin{color}{blue}\textit{#1}\end{color}} \title{Efficient generation of pseudo random numbers based on chaotic iterations on GPU} \begin{document} @@ -25,6 +50,306 @@ Interet de générer des nombres alea sur GPU Présentation des itérations chaotiques +\section{blabla} + +\subsection{The phase space is an interval of the real line} + +\subsubsection{Toward a topological semiconjugacy} + +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. + +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}$. + + +\begin{definition} +The function $\varphi: \mathcal{S}_{10} \times\mathds{B}^{10} \rightarrow \big[ 0, 2^{10} \big[$ is defined by: +$$ +\begin{array}{cccl} +\varphi: & \mathcal{X}_{10} = \mathcal{S}_{10} \times\mathds{B}^{10}& \longrightarrow & \big[ 0, 2^{10} \big[ \\ + & (S,E) = \left((S^0, S^1, \hdots ); (E_0, \hdots, E_9)\right) & \longmapsto & \varphi \left((S,E)\right) +\end{array} +$$ +\noindent where $\varphi\left((S,E)\right)$ is the real number: +\begin{itemize} +\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$. +\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}.$ +\end{itemize} +\end{definition} + + + +$\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: + + +\begin{definition} +\label{def:e et s} +Let $x \in \big[ 0, 2^{10} \big[$ and: +\begin{itemize} +\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}$. +\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: +$\displaystyle{x = \lfloor x \rfloor + \sum_{k=0}^{+\infty} s^k 10^{-k-1}}$. +\end{itemize} +$e$ and $s$ are thus defined as follows: +$$ +\begin{array}{cccl} +e: & \big[ 0, 2^{10} \big[ & \longrightarrow & \mathds{B}^{10} \\ + & x & \longmapsto & (e_0, \hdots, e_9) +\end{array} +$$ +\noindent and +$$ +\begin{array}{cccl} +s: & \big[ 0, 2^{10} \big[ & \longrightarrow & \llbracket 0, 9 \rrbracket^{\mathds{N}} \\ + & x & \longmapsto & (s^k)_{k \in \mathds{N}} +\end{array} +$$ +\end{definition} + +We are now able to define the function $g$, whose goal is to translate the chaotic iterations $\Go$ on an interval of $\mathds{R}$. + +\begin{definition} +$g:\big[ 0, 2^{10} \big[ \longrightarrow \big[ 0, 2^{10} \big[$ is defined by: +$$ +\begin{array}{cccl} +g: & \big[ 0, 2^{10} \big[ & \longrightarrow & \big[ 0, 2^{10} \big[ \\ +& \\ + & x & \longmapsto & g(x) +\end{array} +$$ +\noindent where g(x) is the real number of $\big[ 0, 2^{10} \big[$ defined bellow: +\begin{itemize} +\item its integral part has a binary decomposition equal to $e_0', \hdots, e_9'$, with: +$$ +e_i' = \left\{ +\begin{array}{ll} +e(x)_i & \textrm{ if } i \neq s^0\\ +e(x)_i + 1 \textrm{ (mod 2)} & \textrm{ if } i = s^0\\ +\end{array} +\right. +$$ +\item whose decimal part is $s(x)^1, s(x)^2, \hdots$ +\end{itemize} +\end{definition} + +\bigskip + + +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}}.$$ + +\subsubsection{Defining a metric on $\big[ 0, 2^{10} \big[$} + +Numerous metrics can be defined on the set $\big[ 0, 2^{10} \big[$, the most usual one being the Euclidian distance recalled bellow: + +\begin{notation} +\index{distance!euclidienne} +$\Delta$ is the Euclidian distance on $\big[ 0, 2^{10} \big[$, that is, $\Delta(x,y) = |y-x|^2$. +\end{notation} + +\medskip + +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: + + + +\begin{definition} +Let $x,y \in \big[ 0, 2^{10} \big[$. +$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: +\begin{center} +$\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}}$. +\end{center} +\end{definition} + +\begin{proposition} +$D$ is a distance on $\big[ 0, 2^{10} \big[$. +\end{proposition} + +\begin{proof} +The three axioms defining a distance must be checked. +\begin{itemize} +\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$. +\item $D(x,y)=D(y,x)$. +\item Finally, the triangular inequality is obtained due to the fact that both $\delta$ and $\Delta(x,y)=|x-y|$ satisfy it. +\end{itemize} +\end{proof} + + +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''. +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. + + +\begin{figure}[t] +\begin{center} + \subfigure[Function $x \to dist(x;1,234) $ on the interval $(0;5)$.]{\includegraphics[scale=.35]{DvsEuclidien.pdf}}\quad + \subfigure[Function $x \to dist(x;3) $ on the interval $(0;5)$.]{\includegraphics[scale=.35]{DvsEuclidien2.pdf}} +\end{center} +\caption{Comparison between $D$ (in blue) and the Euclidian distane (in green).} +\label{fig:comparaison de distances} +\end{figure} + + + + +\subsubsection{The semiconjugacy} + +It is now possible to define a topological semiconjugacy between $\mathcal{X}$ and an interval of $\mathds{R}$: + +\begin{theorem} +Chaotic iterations on the phase space $\mathcal{X}$ are simple iterations on $\mathds{R}$, which is illustrated by the semiconjugacy of the diagram bellow: +\begin{equation*} +\begin{CD} +\left(~\mathcal{S}_{10} \times\mathds{B}^{10}, d~\right) @>G_{f_0}>> \left(~\mathcal{S}_{10} \times\mathds{B}^{10}, d~\right)\\ + @V{\varphi}VV @VV{\varphi}V\\ +\left( ~\big[ 0, 2^{10} \big[, D~\right) @>>g> \left(~\big[ 0, 2^{10} \big[, D~\right) +\end{CD} +\end{equation*} +\end{theorem} + +\begin{proof} +$\varphi$ has been constructed in order to be continuous and onto. +\end{proof} + +In other words, $\mathcal{X}$ is approximately equal to $\big[ 0, 2^\mathsf{N} \big[$. + + + + + + +\subsection{Study of the chaotic iterations described as a real function} + + +\begin{figure}[t] +\begin{center} + \subfigure[ICs on the interval $(0,9;1)$.]{\includegraphics[scale=.35]{ICs09a1.pdf}}\quad + \subfigure[ICs on the interval $(0,7;1)$.]{\includegraphics[scale=.35]{ICs07a95.pdf}}\\ + \subfigure[ICs on the interval $(0,5;1)$.]{\includegraphics[scale=.35]{ICs05a1.pdf}}\quad + \subfigure[ICs on the interval $(0;1)$]{\includegraphics[scale=.35]{ICs0a1.pdf}} +\end{center} +\caption{Representation of the chaotic iterations.} +\label{fig:ICs} +\end{figure} + + + + +\begin{figure}[t] +\begin{center} + \subfigure[ICs on the interval $(510;514)$.]{\includegraphics[scale=.35]{ICs510a514.pdf}}\quad + \subfigure[ICs on the interval $(1000;1008)$]{\includegraphics[scale=.35]{ICs1000a1008.pdf}} +\end{center} +\caption{ICs on small intervals.} +\label{fig:ICs2} +\end{figure} + +\begin{figure}[t] +\begin{center} + \subfigure[ICs on the interval $(0;16)$.]{\includegraphics[scale=.3]{ICs0a16.pdf}}\quad + \subfigure[ICs on the interval $(40;70)$.]{\includegraphics[scale=.45]{ICs40a70.pdf}}\quad +\end{center} +\caption{General aspect of the chaotic iterations.} +\label{fig:ICs3} +\end{figure} + + +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: + +\begin{proposition} +\label{Prop:derivabilite des ICs} +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\}$. + +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$. +\end{proposition} + + +\begin{proof} +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$: +\begin{itemize} +\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}$). +\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$. +\end{itemize} +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$. +\end{proof} + +\begin{remark} +Finally, chaotic iterations are elements of the large family of functions that are both chaotic and piecewise linear (like the tent map). +\end{remark} + + + +\subsection{Comparison of the two metrics on $\big[ 0, 2^\mathsf{N} \big[$} + +The two propositions bellow allow to compare our two distances on $\big[ 0, 2^\mathsf{N} \big[$: + +\begin{proposition} +Id: $\left(~\big[ 0, 2^\mathsf{N} \big[,\Delta~\right) \to \left(~\big[ 0, 2^\mathsf{N} \big[, D~\right)$ is not continuous. +\end{proposition} + +\begin{proof} +The sequence $x^n = 1,999\hdots 999$ constituted by $n$ 9 as decimal part, is such that: +\begin{itemize} +\item $\Delta (x^n,2) \to 0.$ +\item But $D(x^n,2) \geqslant 1$, then $D(x^n,2)$ does not converge to 0. +\end{itemize} + +The sequential characterization of the continuity concludes the demonstration. +\end{proof} + + + +A contrario: + +\begin{proposition} +Id: $\left(~\big[ 0, 2^\mathsf{N} \big[,D~\right) \to \left(~\big[ 0, 2^\mathsf{N} \big[, \Delta ~\right)$ is a continuous fonction. +\end{proposition} + +\begin{proof} +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$. + +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. +\end{proof} + +The conclusion of these propositions is that the proposed metric is more precise than the Euclidian distance, that is: + +\begin{corollary} +$D$ is finer than the Euclidian distance $\Delta$. +\end{corollary} + +This corollary can be reformulated as follows: + +\begin{itemize} +\item The topology produced by $\Delta$ is a subset of the topology produced by $D$. +\item $D$ has more open sets than $\Delta$. +\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$. +\end{itemize} + + +\subsection{Chaos of the chaotic iterations on $\mathds{R}$} +\label{chpt:Chaos des itérations chaotiques sur R} + + + +\subsubsection{Chaos according to Devaney} + +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: +\begin{itemize} +\item $\left(\Go, \mathcal{X}_d\right)$ and $\left(g, \big[ 0, 2^{10} \big[_D\right)$ are semiconjugate by $\varphi$, +\item Then $\left(g, \big[ 0, 2^{10} \big[_D\right)$ is a system chaotic according to Devaney, because the semiconjugacy preserve this character. +\item But the topology generated by $D$ is finer than the topology generated by the Euclidian distance $\Delta$ -- which is the order topology. +\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}$. +\end{itemize} + +This result can be formulated as follows. + +\begin{theorem} +\label{th:IC et topologie de l'ordre} +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. +\end{theorem} + +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. + + + + \section{Efficient prng based on chaotic iterations} On parle du séquentiel avec des nombres 64 bits\\ @@ -47,12 +372,12 @@ unsigned int CIprng() \{\\ unsigned long t1 = xorshift();\\ unsigned long t2 = xor128();\\ unsigned long t3 = xorwow();\\ - x = x\^\ (unsigned int)t1;\\ - x = x\^\ (unsigned int)(t2$>>$32);\\ - x = x\^\ (unsigned int)(t3$>>$32);\\ - x = x\^\ (unsigned int)t2;\\ - x = x\^\ (unsigned int)(t1$>>$32);\\ - x = x\^\ (unsigned int)t3;\\ + x = x\textasciicircum (unsigned int)t1;\\ + x = x\textasciicircum (unsigned int)(t2$>>$32);\\ + x = x\textasciicircum (unsigned int)(t3$>>$32);\\ + x = x\textasciicircum (unsigned int)t2;\\ + x = x\textasciicircum (unsigned int)(t1$>>$32);\\ + x = x\textasciicircum (unsigned int)t3;\\ return x;\\ \} \end{minipage} @@ -63,12 +388,14 @@ unsigned int CIprng() \{\\ \end{figure} In Figure~\ref{algo:seqCIprng} a sequential version of our chaotic iterations -based PRNG is presented. This version uses three classical 64-bits PRNG: the +based PRNG is presented. This version uses three classical 64 bits PRNG: the \texttt{xorshift}, the \texttt{xor128} and the \texttt{xorwow}. These three -PRNGs are presented in~\cite{Marsaglia2003}. As each PRNG used works with -64-bits and as our PRNG works with 32-bits, the use of \texttt{(unsigned int)} +PRNGs are presented in~\cite{Marsaglia2003}. As each PRNG used works with +64-bits and as our PRNG works with 32 bits, the use of \texttt{(unsigned int)} selects the 32 least significant bits whereas \texttt{(unsigned int)(t3$>>$32)} -selects the 32 bits most significants bits of the variable \texttt{t}. +selects the 32 most significants bits of the variable \texttt{t}. This version +sucesses the BigCrush of the TestU01 battery [P. L’ecuyer and + R. Simard. Testu01]. \section{Efficient prng based on chaotic iterations on GPU} @@ -80,7 +407,6 @@ On passe le BigCrush\\ On donne des temps de générations sur GPU/CPU\\ On donne des temps de générations de nombre sur GPU puis on rappatrie sur CPU / CPU ? bof bof, on verra -\section{Lyapunov} \section{Conclusion} \bibliographystyle{plain}