1 Introduction

Chaotic pseudo-random number generators (cPRNGs), in particular the logistic map [1], have been considered as a means for secure PRNGs. Experiments [11], [13] revealed that the period lengths are not satisfactory in relation to state size. When studying the descriptions of these experiments, we found that they are missing some implementation details that also influence the period length, even for fixed parameter settings. Hence, replicating those studies experiences some challenges, which we report and from which we derive recommendations for describing experiments including numeric calculations. We also underline experimentally that the implementation details such as chosen number representation, rounding mode and order of evaluation in a formula influence the results.

Beyond the typical number representations integer, i. e. fixed-point, and (binary) floating-point, we next present number representation formats in-between those two that use radix higher than 2 and possibly mixed radix. We demonstrate experimentally that those formats enable to trade accuracy and storage. Thus, our investigations also exhibit relations to the field of approximate computing. Moreover, we investigate if exploiting the symmetry of the logistic map function, i. e. f ( x ) = f ( 1 − x ) for x ∈ [ 0 ; 1 ] can be mirrored in the number representation in order to achieve longer period lengths.

The remainder of this article is structured as follows. Section 2 gives information on pseudo-random number generation for security applications, especially when using chaotic functions, summarizes basic information on number representations for rational numbers, and discusses related work. Section 3 presents the challenges that we encountered in experiments with chaotic pseudo-random number generation based on logistic map and reports on own experiments to demonstrate those challenges. In Section 4 we present number formats between standard floating-point and fixed-point representations, as well as experimental results that indicate that those formats allow a tradeoff between storage size and achievable period lengths. Section 5 investigates how symmetry in the logistic map function can be exploited in number representations, while Section 6 concludes with recommendations on configuration details to be reported about numerical experiments, and an outlook on future work.

2 Background

3 Challenges

In this section, we will investigate challenges for replication of execution results, when evaluating the formula a · x · ( 1 − x ). Those challenges illustrate that even on such microscopic scale, replication is difficult, and that therefore, replication of experiments on a larger scale must work with the assumption that certain differences will be inevitable.

4 Comparing different number formats

Fixed-point and floating-point representations are just the opposite ends of a whole spectrum of possible number formats. In the following, we describe a generic scheme to specify different number formats. To do so, we partition the interval [0;1) into a sequence of d intervals I 0 to I d − 1 , where I 0 = [ 0 ; 2 e 0 ), I i = [ 2 e i − 1 ; 2 e i ) for i ≥ 1, and the e i are a strictly increasing sequence of integers, with e d − 1 = 0 to enforce that interval I d − 1 ends with 1. In each interval, we have (up to) 2 m representations that are equi-spaced within the interval by representing multiples of 2 e i − m , i. e. the represented number is

Thus, to store a number in this format, we need e = ⌈ log 2 d ⌉ bits to represent the interval index, and m bits to represent the mantissa (or m − 1 in case of binary floating-point representations, see below).

We start with fixed-point representations. They only comprise one interval I 0 , i. e. d = 1, and thus the number of exponent bits is zero. The mantissa bits represent a multiple of 2 − m , i. e. the represented number is r ( 0 , M ) = ⟨ 0 . M ⟩ = ⟨ M ⟩ · 2 − m .

In single precision floating-point representations, d = 127, and interval I 0 contains the denormalized numbers. The exponents form a sequence e i = i − 126, i. e. successive exponents only differ by 1. To ensure that r ( i , M ) ≥ 2 e i − 1 for i ≥ 1, i. e. that r ( i , M ) ∈ I i , it follows from Eq. (2) that ⟨ 0 . M ⟩ ≥ 0.5, i. e. that the uppermost bit of mantissa M must be set (so-called implicit 1). Thus, only m − 1 bits must be stored, or put otherwise, if we store m bits, we really work with a m + 1-bit mantissa.

Another example of a possible format is a separation in d = 4 intervals, i. e. using e = 2 exponent bits, with ( e i ) i = ( − 32 , − 16 , − 10 , 0 ). Thus, no particular regularity is necessary beyond the restrictions above. The mantissa is shorter than for fixed-point representations, but longer than for floating-point representations, so a possible advantage in cycle length might arise.

If all differences e i − e i − 1 are equal, then we can consider the representation to be designed with respect to higher radix 2 e i − e i − 1 , if they are not equal, they can be considered multi-radix representations.

Hence, for given total storage size of t = e + m bits, there are different numbers d = 2 e of possible exponents, leading to different mantissa lengths m, and for each d the exponents can have a number of possible distances leading to different number formats, which will in turn lead to different cycle lengths in the logistic map application. In this configuration space, we can search for the representation which is best with respect to a particular application.

Alternatively, one can fix the required performance, i. e. cycle length in our application, and search for the representation with shortest storage size that achieves that performance, taking up the view of approximate computing [19].

Moreover, from all combinations of performance and storage size, the user can draw a plot and identify the Pareto front, i. e. the curve of Pareto-optimal combinations, which is the (typically small) number of combinations that the user can choose from depending on her needs [3].

In general, if e i − e i − 1 ≥ m, then all representations of Eq. (2) (if not all mantissa bits are 0) are in I i . Otherwise, the uppermost e i − e i − 1 mantissa bits may not all be 0 simultaneously. As there is no similar trick as the implicit 1 if e i − e i − 1 ≥ 2, a fraction of 2 e i − e i − 1 of the representations does not belong to I i but to I j with j < i, i. e. there is a certain redundancy, i. e. reduction in the total number of distinct representations.

To explore the spectrum of possible number formats, we have implemented a custom library where d, M and the sequence of e i can be specified for a representation type.^[1]

We investigate total lengths t from 32 down to 24, with mantissa lengths ranging from t (resulting in e = 0) to 24 (resulting in e = t − 24). For each combination of e and m we explore sequences of e i with regular distances e i − e i − 1 from the set {1,2,4,6,8,10}. So far, only round-to-zero is supported and multiplication with a = 3 + ε is performed last, by first adding x ( 1 − x ) to itself, by computing ε x ( 1 − x ) and adding x ( 1 − x ), and finally adding the two intermediate results.

For each number representation, we choose again 63 different values of a as before, and compute the minimum, maximum and average cycle lengths (plus standard deviation).

Like in previous experiments, we see that standard deviation is quite large (ranging 49 % to 69 % of the average) and does not get better if we take larger number of samples. For almost all (31 of 33) combinations of e and m, we find that e i − e i − 1 = 1 leads to largest average cycle lengths. In the two exceptions, e i − e i − 1 = 2 performs better: in one case ( e = 1, m = 29) the difference is only 0.32 %, yet in the other case ( e = 1, m 3 0) the difference is 30.9 %, i. e. notably, mainly driven by a really large maximum cycle length achieved. Thus that case can be considered an outlier.

In the following, we further investigate the e , m combinations with distance leading to maximum average cycle length. Figure 1 depicts that the average cycle length is mainly determined by the number of mantissa bits, for each m it is largely constant over the different exponent lengths e, and the curves are ordered by growing m. Please note that the cycle lengths for e = 0 are taken from Table 2 where different rounding mode and multiplication order are applied, so that the step from e = 0 to e = 1 is not completely fair, which might explain the unusual shape for m = 31.

Figure 2 depicts that for a fixed exponent length e, the average cycle length grows with total length t = e + m, i. e. with mantissa length m, and the curves are ordered by shrinking exponent length e. The only exception is a crossing of the curves for e = 1 and e = 2 at the right end of the curves. Figure 3 completes this view by plotting curves for each total length t = e + m. We see that each curve decreases with growing e, i. e. shrinking m, and the curves are ordered by t.

Figure 1

Dependence of average cycle length on exponent length e for different mantissa lengths m.

Figure 2

Dependence of average cycle length on mantissa length m for different exponent lengths e.

Figure 3

Dependence of average cycle length on exponent length e for different total lengths t = e + m.

In summary, we find that for our application, distances e i − e i − 1 = 1 are best, as in binary floating-point representations, and that fixed-point representation leads indeed to longest average cycle lengths, even when we consider the large space of representations between (binary) floating-point and fixed-point representations.

Still, our view unifies several proposals for number representations from literature, cf. Section 2.3, and our implementation allows to find the optimal representation, e. g. for a reconfigurable hardware device with fixed application, with low effort.

5 Extending number formats for symmetry

We note that the logistic map as defined in Eq. (1) is symmetric with middle 0.5, i. e. for z ≤ 0.5 we have f ( 0.5 − z ) = f ( 0.5 + z ) because for x = 0.5 − z we get 1 − x = 0.5 + z (and vice versa) and from Eq. (1) we get f ( x ) = f ( 1 − x ). Hence we can extend previous number formats to respect this symmetry. If z is any of the previous representations, only restricted to [0;0.5] instead of [0;1], then ( s , z ) with x = 0.5 + ( − 1 ) s · z is a representation for [0;1] where representation densities are symmetric around 0.5, and 1 − x can be represented by ( 1 − s , z ), i. e. by inverting the “sign” bit s. The product x · ( 1 − x ) can be expressed as 0.25 − z 2 independent of s, and

If z ˜ = a / 4 + a · z 2 ≥ 0.5, then z ′ is non-negative, and f ( x ) ≥ 0.5 can be represented by ( 0 , z ′ ). Otherwise, − z ′ is non-negative and f ( x ) < 0.5 can be represented by ( 1 , − z ′ ). Thus, one first computes z ˜, and depending on whether z ˜ ≥ 0.5 or not, one computes z ˜ − 0.5 or 0.5 − z ˜ to determine z ′ or − z ′ , respectively. If ( s , z ) can be represented together as in floating-point format, then one can simply compute z ′ ∈ [ − 0.5 ; + 0.5 ]. Please note that in this case almost twice as many values of IEEE754 floats can be used: on the one hand one value less from the exponent range is used (but still close to half of the values), but on the other hand the sign bit is used.

To see if “symmetric” floating-point representations are advantageous, we did experiments similar to those in Section 3.4: for each a tested, we find the largest cycle reachable from 32 starting points, and we compute minimum, maximum, and average cycle length over 32 values of a starting from 4 − 2 − 12 and equi-spaced with distance 2 − 17 . We only use rounding to nearest, but compare different orders of computation. For ordinary single precision (32 bit) IEEE754 floating-point numbers we distinguish multiplication with a last, or x last, or 1 − x last. For symmetric single precision floating-point numbers we distinguish multiplication with a last, i. e. computation of z 2 first, or multiplication with z last, i. e. computation of a · z first, and we distinguish subtraction of 0.5 last, or subtraction of a · z 2 last. Moreover, we compute minimum, maximum and average over all possibilities.

The results can be seen in Table 3. There is no clear winner. With respect to average cycle length found, ordinary floating-point representation with multiplication order “a last” is best, but the best variant of symmetric floating-point representations is much better in minimum (and maximum) cycle length found. Moreover, the standard deviations are very large compared to the averages, so that the averages are not very meaningful, and minimum and maximum lengths might provide a better insight. Please note that increasing the number of trials over which averages are computed does not help here, due to the chaotic nature of the logistic map for different values of a close to 4. For 1024 equi-spaced values of a starting from 4 − 2 − 12 , the standard deviation does not shrink substantially, e. g. from 1610 to 1523 for symmetric float representations and multiplication with a last. Yet, in some variants standard deviation even increases, e. g. from 1630 to 1688 for conventional floats and multiplication with 1 − x last. The biggest difference is that the average cycle length for ordinary floats where multiplication with a is done last reduces to 2459.7, so that the range of averages is much closer together.

6 Conclusions

We have presented challenges that occur with respect to replicability in security or cryptographic applications that involve numeric parts originally devised for real numbers. As example, we have chosen the logistic map, which has been one proposal for a chaotic pseudo-random number generator. We have demonstrated by experiments the relevance of considering the differences that possibly result from not arbitrating the ambiguities underlying the above challenges.

As similar problems might occur in other security applications as well, we derive the following recommendations. First, the rounding mode used in the original experiments should be found out by the experimenters and should be reported. Second, the order of computations used in the original experiments should be reported as well, and possibly should be modified by the first experimenters to check if notable differences in result can be the consequence. Finally, different number representations should be considered by experimenters, and results on their influence and reasons for the final choice of representation should be reported.

Moreover, we have investigated higher and mixed radix representations in between fixed-point and binary floating-point representations, plus symmetric representations, and demonstrated how their tradeoff between storage size and performance (in the sense of period length in the application) can be used by users to tailor their application to their needs, thus connecting this investigation to the field of approximate computing.

Future work will comprise a deeper analysis of our experimental results. Also, the scope of numerical challenges will be extended from chaotic random numbers to other branches of IT security. Moreover, the introduction of redundancy by higher radix representations, e. g. by skipping normalization [2], opens the path for information hiding, e. g. storage covert channels [18], in such representations, which we plan to investigate further.