NNAlgorithm/opt.py at master · FloydZ/NNAlgorithm · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
#!/usr/bin/env python3
import math


def HH(i: float):
    if i == 1.0 or i == 0.0:
        return 0.0

    if i > 1.0 or i < 0.0:
        print("error: ", i)
        raise ValueError

    return -(i * math.log2(i) + (1 - i) * math.log2(1 - i))


def H1(value: float):
    if value == 1.0:
        return 0.5

    # approximate inverse binary entropy function
    steps = [0.1, 0.01, 0.001, 0.0001, 0.00001, 0.000001, 0.0000001, 0.00000001, 0.0000000001, 0.0000000000001, 0.000000000000001]
    r = 0.000000000000000000000000000000001
    if value < 0.:
        value = 0.

    for step in steps:
        i = r
        while (i + step < 1.0) and (HH(i) < value):
            i += step

        if step > i:
            break
        r = i - step

    return r


def calc_q_(n: int, w_: float, d_: float):
    """
    compute probability of
    """
    assert w_ < n
    assert d_ < n
    wm1 = 1 - w_
    return wm1 * HH((d_ - (w_ / 2)) / wm1)


def compute_optimal_params(n: int, lam, w, r=0):
    """
    taken from https://arxiv.org/pdf/2102.02597.pdf

    lambda = list size
    gamma = weight to diff
    delta = weight to match exactly on each bucket_window (relative to k)
    k = number of coords per window
    :return
        r: number of sublimbs
        N: number of leaves
        d: weight on each window
        k: window size

    """
    # convert to relative numbers
    if type(lam) == int or lam > 30:
        lam = math.log2(lam)/n

    if type(w) == int:
        w = w / n

    # only set r if explicity wanted
    if r == 0:
        r = n / (math.log2(n))**2

    assert(lam <= 1.)

    delta_star = H1(1. - lam)
    limit = 2.*delta_star * (1. - delta_star)
    if w > limit:
        d = 1/2. * (1. - math.sqrt(1. - 2.*w))
    else:
        d = delta_star

    q = calc_q_(n, w, d)
    assert q <= 1.

    N = int(n / q)

    k = n/r

    return r, N, d*k, k


def compute_time(n, lam, w, r, N, d):
    """
    compute the expected runtime in log2
    """

    if type(lam) == int or lam > 30:
        lam = math.log2(lam)/n

    w = w
    if type(w) == int:
        w = w / n

    delta_star = H1(1. - lam)
    w_star = 2*delta_star * (1. - delta_star)
    print(w, w_star, lam)
    if w <= w_star:
        theta = (1. - w) * (1.0 - HH((delta_star - w/2.)/(1. - w)))
    else:
        theta = 2*lam + HH(w) - 1.

    return theta*n


# TODO compute intermediate list size
n = 256
lam = 1 << 14
w = 16
r, N, d, k = compute_optimal_params(n, lam, w)
print("r", r, "N", N, "d", d)

time = compute_time(n, lam, w, r, N, d)
print("theta", time, 2**time)