非线性动力系统的 NIST Suite 测试

Question

在我的以下代码中，我正在运行一个洛伦兹混沌方程，从中我将得到 xs、ys 和 zs 方面的随机数

import numpy as np
def lorenz(x, y, z, a=10,b=8/3,c=28 ):
            x_dot = a*(y -x) 
            y_dot = - y +c*x - x*z
            z_dot = -b*z + x*y
            return x_dot, y_dot, z_dot
       
dt = 0.01
num_steps =  10000
# Need one more for the initial values
xs = np.empty(num_steps + 1)
ys = np.empty(num_steps + 1)
zs = np.empty(num_steps + 1)


# Set initial values
xs[0], ys[0], zs[0]= (1,1,1)
# Step through "time", calculating the partial derivatives at the current point
# and using them to estimate the next point
for i in range(num_steps):
    x_dot, y_dot, z_dot= lorenz(xs[i], ys[i], zs[i])
    xs[i + 1] = xs[i] + (x_dot * dt)
    ys[i + 1] = ys[i] + (y_dot * dt)
    zs[i + 1] = zs[i] + (z_dot * dt)

我实际上是在尝试使用下面的代码通过 NIST 800 测试随机数生成测试的 xs、ys 和 zs 值

from __future__ import print_function

import math
from fractions import Fraction
from scipy.special import gamma, gammainc, gammaincc
# from gamma_functions import *
import numpy
import cmath
import random

#ones_table = [bin(i)[2:].count('1') for i in range(256)]
def count_ones_zeroes(bits):
    ones = 0
    zeroes = 0
    for bit in bits:
        if (bit == 1):
            ones += 1
        else:
            zeroes += 1
    return (zeroes,ones)

def runs_test(bits):
    n = len(bits)
    zeroes,ones = count_ones_zeroes(bits)

    prop = float(ones)/float(n)
    print("  prop ",prop)

    tau = 2.0/math.sqrt(n)
    print("  tau ",tau)

    if abs(prop-0.5) > tau:
        return (False,0.0,None)

    vobs = 1.0
    for i in range(n-1):
        if bits[i] != bits[i+1]:
            vobs += 1.0

    print("  vobs ",vobs)
      
    p = math.erfc(abs(vobs - (2.0*n*prop*(1.0-prop)))/(2.0*math.sqrt(2.0*n)*prop*(1-prop) ))
    success = (p >= 0.01)
    return (success,p,None)
print(runs_test(xs))
#%%
from __future__ import print_function

import math

def count_ones_zeroes(bits):
    ones = 0
    zeroes = 0
    for bit in bits:
        if (bit == 1):
            ones += 1
        else:
            zeroes += 1
    return (zeroes,ones)

def monobit_test(bits):
    n = len(bits)
    
    zeroes,ones = count_ones_zeroes(bits)
    s = abs(ones-zeroes)
    print("  Ones count   = %d" % ones)
    print("  Zeroes count = %d" % zeroes)
    
    p = math.erfc(float(s)/(math.sqrt(float(n)) * math.sqrt(2.0)))
    
    success = (p >= 0.01)
    return (success,p,None)
print(runs_test(xs))

我得到的输出是假的，即

输出： 道具 0.00019998000199980003 头 0.01999900007499375 （假，0.0，无）

我现在该怎么办？

Answer 1

洛伦兹系统是混乱的，不是随机的。您很好地实现了微分方程求解器，但似乎count_ones_zeroes 并没有像它的名字所暗示的那样做，至少在您提供的数据上没有。在xs 上，它返回(zeroes, ones) = (9999, 2)，这不是您想要的。代码检查xs 数组中的值，即x 值（例如8.2）对1，但x 是介于-20 和20 之间的浮点数，因此它通常为非1，并且将被计为0. 只有x==1会被计为1。

在 python 中，int/int 会产生浮点数，因此无需将其强制转换为浮点数，这与例如C 或 C++，因此您可以编写 prop = ones/n，而不是 prop = float(ones)/float(n)，类似的语句适用于 +、- 和 *