range_minimum_query 0.2.0

Range Minimum Query (RMQ) is used on arrays to find the position of an element with the minimum value between two specified indices.
Documentation 
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
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254

pub const ONES_STEP_9: u64 = 1 << 0 | 1 << 9 | 1 << 18 | 1 << 27 | 1 << 36 | 1 << 45 | 1 << 54;

pub const INV_COUNT_STEP_9: u64 = 1 << 54 | 2 << 45 | 3 << 36 | 4 << 27 | 5 << 18 | 6 << 9 | 7;

pub const MSBS_STEP_9: u64 = 0x100 * ONES_STEP_9;

pub const fn uleq_step_9(x: u64, y: u64) -> u64 {
    (((((y | MSBS_STEP_9) - (x & !MSBS_STEP_9)) | (x ^ y)) ^ (x & !y)) & MSBS_STEP_9) >> 8
}

const PS_OVERFLOW: [u64; 65] = [
    0x8080_8080_8080_8080,
    0x7f7f_7f7f_7f7f_7f7f,
    0x7e7e_7e7e_7e7e_7e7e,
    0x7d7d_7d7d_7d7d_7d7d,
    0x7c7c_7c7c_7c7c_7c7c,
    0x7b7b_7b7b_7b7b_7b7b,
    0x7a7a_7a7a_7a7a_7a7a,
    0x7979_7979_7979_7979,
    0x7878_7878_7878_7878,
    0x7777_7777_7777_7777,
    0x7676_7676_7676_7676,
    0x7575_7575_7575_7575,
    0x7474_7474_7474_7474,
    0x7373_7373_7373_7373,
    0x7272_7272_7272_7272,
    0x7171_7171_7171_7171,
    0x7070_7070_7070_7070,
    0x6f6f_6f6f_6f6f_6f6f,
    0x6e6e_6e6e_6e6e_6e6e,
    0x6d6d_6d6d_6d6d_6d6d,
    0x6c6c_6c6c_6c6c_6c6c,
    0x6b6b_6b6b_6b6b_6b6b,
    0x6a6a_6a6a_6a6a_6a6a,
    0x6969_6969_6969_6969,
    0x6868_6868_6868_6868,
    0x6767_6767_6767_6767,
    0x6666_6666_6666_6666,
    0x6565_6565_6565_6565,
    0x6464_6464_6464_6464,
    0x6363_6363_6363_6363,
    0x6262_6262_6262_6262,
    0x6161_6161_6161_6161,
    0x6060_6060_6060_6060,
    0x5f5f_5f5f_5f5f_5f5f,
    0x5e5e_5e5e_5e5e_5e5e,
    0x5d5d_5d5d_5d5d_5d5d,
    0x5c5c_5c5c_5c5c_5c5c,
    0x5b5b_5b5b_5b5b_5b5b,
    0x5a5a_5a5a_5a5a_5a5a,
    0x5959_5959_5959_5959,
    0x5858_5858_5858_5858,
    0x5757_5757_5757_5757,
    0x5656_5656_5656_5656,
    0x5555_5555_5555_5555,
    0x5454_5454_5454_5454,
    0x5353_5353_5353_5353,
    0x5252_5252_5252_5252,
    0x5151_5151_5151_5151,
    0x5050_5050_5050_5050,
    0x4f4f_4f4f_4f4f_4f4f,
    0x4e4e_4e4e_4e4e_4e4e,
    0x4d4d_4d4d_4d4d_4d4d,
    0x4c4c_4c4c_4c4c_4c4c,
    0x4b4b_4b4b_4b4b_4b4b,
    0x4a4a_4a4a_4a4a_4a4a,
    0x4949_4949_4949_4949,
    0x4848_4848_4848_4848,
    0x4747_4747_4747_4747,
    0x4646_4646_4646_4646,
    0x4545_4545_4545_4545,
    0x4444_4444_4444_4444,
    0x4343_4343_4343_4343,
    0x4242_4242_4242_4242,
    0x4141_4141_4141_4141,
    0x4040_4040_4040_4040,
];

const LT_SEL: [u8; 2048] = [
    0, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    6, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    7, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    6, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    5, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0, 4, 0, 1, 0, 2, 0, 1, 0, 3, 0, 1, 0, 2, 0, 1, 0,
    0, 0, 0, 1, 0, 2, 2, 1, 0, 3, 3, 1, 3, 2, 2, 1, 0, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    0, 5, 5, 1, 5, 2, 2, 1, 5, 3, 3, 1, 3, 2, 2, 1, 5, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    0, 6, 6, 1, 6, 2, 2, 1, 6, 3, 3, 1, 3, 2, 2, 1, 6, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    6, 5, 5, 1, 5, 2, 2, 1, 5, 3, 3, 1, 3, 2, 2, 1, 5, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    0, 7, 7, 1, 7, 2, 2, 1, 7, 3, 3, 1, 3, 2, 2, 1, 7, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    7, 5, 5, 1, 5, 2, 2, 1, 5, 3, 3, 1, 3, 2, 2, 1, 5, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    7, 6, 6, 1, 6, 2, 2, 1, 6, 3, 3, 1, 3, 2, 2, 1, 6, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    6, 5, 5, 1, 5, 2, 2, 1, 5, 3, 3, 1, 3, 2, 2, 1, 5, 4, 4, 1, 4, 2, 2, 1, 4, 3, 3, 1, 3, 2, 2, 1,
    0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 3, 0, 3, 3, 2, 0, 0, 0, 4, 0, 4, 4, 2, 0, 4, 4, 3, 4, 3, 3, 2,
    0, 0, 0, 5, 0, 5, 5, 2, 0, 5, 5, 3, 5, 3, 3, 2, 0, 5, 5, 4, 5, 4, 4, 2, 5, 4, 4, 3, 4, 3, 3, 2,
    0, 0, 0, 6, 0, 6, 6, 2, 0, 6, 6, 3, 6, 3, 3, 2, 0, 6, 6, 4, 6, 4, 4, 2, 6, 4, 4, 3, 4, 3, 3, 2,
    0, 6, 6, 5, 6, 5, 5, 2, 6, 5, 5, 3, 5, 3, 3, 2, 6, 5, 5, 4, 5, 4, 4, 2, 5, 4, 4, 3, 4, 3, 3, 2,
    0, 0, 0, 7, 0, 7, 7, 2, 0, 7, 7, 3, 7, 3, 3, 2, 0, 7, 7, 4, 7, 4, 4, 2, 7, 4, 4, 3, 4, 3, 3, 2,
    0, 7, 7, 5, 7, 5, 5, 2, 7, 5, 5, 3, 5, 3, 3, 2, 7, 5, 5, 4, 5, 4, 4, 2, 5, 4, 4, 3, 4, 3, 3, 2,
    0, 7, 7, 6, 7, 6, 6, 2, 7, 6, 6, 3, 6, 3, 3, 2, 7, 6, 6, 4, 6, 4, 4, 2, 6, 4, 4, 3, 4, 3, 3, 2,
    7, 6, 6, 5, 6, 5, 5, 2, 6, 5, 5, 3, 5, 3, 3, 2, 6, 5, 5, 4, 5, 4, 4, 2, 5, 4, 4, 3, 4, 3, 3, 2,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3, 0, 0, 0, 0, 0, 0, 0, 4, 0, 0, 0, 4, 0, 4, 4, 3,
    0, 0, 0, 0, 0, 0, 0, 5, 0, 0, 0, 5, 0, 5, 5, 3, 0, 0, 0, 5, 0, 5, 5, 4, 0, 5, 5, 4, 5, 4, 4, 3,
    0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 6, 0, 6, 6, 3, 0, 0, 0, 6, 0, 6, 6, 4, 0, 6, 6, 4, 6, 4, 4, 3,
    0, 0, 0, 6, 0, 6, 6, 5, 0, 6, 6, 5, 6, 5, 5, 3, 0, 6, 6, 5, 6, 5, 5, 4, 6, 5, 5, 4, 5, 4, 4, 3,
    0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 3, 0, 0, 0, 7, 0, 7, 7, 4, 0, 7, 7, 4, 7, 4, 4, 3,
    0, 0, 0, 7, 0, 7, 7, 5, 0, 7, 7, 5, 7, 5, 5, 3, 0, 7, 7, 5, 7, 5, 5, 4, 7, 5, 5, 4, 5, 4, 4, 3,
    0, 0, 0, 7, 0, 7, 7, 6, 0, 7, 7, 6, 7, 6, 6, 3, 0, 7, 7, 6, 7, 6, 6, 4, 7, 6, 6, 4, 6, 4, 4, 3,
    0, 7, 7, 6, 7, 6, 6, 5, 7, 6, 6, 5, 6, 5, 5, 3, 7, 6, 6, 5, 6, 5, 5, 4, 6, 5, 5, 4, 5, 4, 4, 3,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 4,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 5, 0, 0, 0, 0, 0, 0, 0, 5, 0, 0, 0, 5, 0, 5, 5, 4,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 6, 0, 6, 6, 4,
    0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 6, 0, 6, 6, 5, 0, 0, 0, 6, 0, 6, 6, 5, 0, 6, 6, 5, 6, 5, 5, 4,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 4,
    0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 5, 0, 0, 0, 7, 0, 7, 7, 5, 0, 7, 7, 5, 7, 5, 5, 4,
    0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 6, 0, 0, 0, 7, 0, 7, 7, 6, 0, 7, 7, 6, 7, 6, 6, 4,
    0, 0, 0, 7, 0, 7, 7, 6, 0, 7, 7, 6, 7, 6, 6, 5, 0, 7, 7, 6, 7, 6, 6, 5, 7, 6, 6, 5, 6, 5, 5, 4,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 5,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 0, 0, 0, 0, 6, 0, 0, 0, 6, 0, 6, 6, 5,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 5,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 6,
    0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 6, 0, 0, 0, 7, 0, 7, 7, 6, 0, 7, 7, 6, 7, 6, 6, 5,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 6,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 0, 0, 0, 0, 7, 0, 0, 0, 7, 0, 7, 7, 6,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 7,
];

#[cfg(target_feature = "bmi2")]
pub unsafe fn select1_in_u64_bmi2(a: u64, k: u64) -> u64 {
    let i = 1 << k;
    let y = core::arch::x86_64::_pdep_u64(i, a);
    core::arch::x86_64::_tzcnt_u64(y)
}

// taken from sdsl-lite
pub const fn select1_in_u64_slow(x: u64, mut i: u64) -> u64 {
    i += 1; // both are 1 based
    let mut s = x; // s = sum
    s = s - ((s >> 1) & 0x5555_5555_5555_5555);
    s = (s & 0x3333_3333_3333_3333) + ((s >> 2) & 0x3333_3333_3333_3333);
    s = (s + (s >> 4)) & 0x0F0F_0F0F_0F0F_0F0F;
    s = s.wrapping_mul(0x0101_0101_0101_0101);
    let b = s + PS_OVERFLOW[i as usize]; //&0x8080808080808080;// add something to the partial sums to cause overflow
    i = (i - 1) << 8;
    if (b & 0x0000_0000_8000_0000) != 0 {
        // byte <=3
        if (b & 0x0000_0000_0000_8000) != 0 {
            //byte <= 1
            if (b & 0x0000_0000_0000_0080) != 0 {
                LT_SEL[((x & 0xFF) + i) as usize] as u64
            } else {
                8 + LT_SEL[((((x >> 8) & 0xFF) + i - ((s & 0xFF) << 8)) & 0x7FF) as usize] as u64
                //byte 1;
            }
        } else {
            //byte >1
            if (b & 0x0000_0000_0080_0000) != 0 {
                //byte <=2
                16 + LT_SEL[((((x >> 16) & 0xFF) + i - (s & 0xFF00)) & 0x7FF) as usize] as u64
            //byte 2;
            } else {
                24 + LT_SEL[((((x >> 24) & 0xFF) + i - ((s >> 8) & 0xFF00)) & 0x7FF) as usize]
                    as u64
                //byte 3;
            }
        }
    } else {
        //  byte > 3
        if (b & 0x0000_8000_0000_0000) != 0 {
            // byte <=5
            if (b & 0x0000_0080_0000_0000) != 0 {
                //byte <=4
                32 + LT_SEL[((((x >> 32) & 0xFF) + i - ((s >> 16) & 0xFF00)) & 0x7FF) as usize]
                    as u64
            //byte 4;
            } else {
                40 + LT_SEL[((((x >> 40) & 0xFF) + i - ((s >> 24) & 0xFF00)) & 0x7FF) as usize]
                    as u64
                //byte 5;
            }
        } else {
            // byte >5
            if (b & 0x0080_0000_0000_0000) != 0 {
                //byte<=6
                48 + LT_SEL[((((x >> 48) & 0xFF) + i - ((s >> 32) & 0xFF00)) & 0x7FF) as usize]
                    as u64
            //byte 6;
            } else {
                56 + LT_SEL[((((x >> 56) & 0xFF) + i - ((s >> 40) & 0xFF00)) & 0x7FF) as usize]
                    as u64
                //byte 7;
            }
        }
    }
}

pub fn select1_in_u64(a: u64, k: u64) -> u64 {
    #[cfg(target_feature = "bmi2")]
    {
        unsafe { select1_in_u64_bmi2(a, k) }
    }
    #[cfg(not(target_feature = "bmi2"))]
    {
        select1_in_u64_slow(a, k)
    }
}

#[cfg(test)]
mod tests {
    use bitvec::prelude::*;
    use proptest::prelude::*;

    proptest! {
        #[test]
        fn select1_in_u64(data: u64) {
            let bv = BitVec::<_, Lsb0>::from_element(data);
            for (found, one_pos) in bv.iter_ones().enumerate() {
                let comp_pos = super::select1_in_u64(data,found as u64);
                assert_eq!(one_pos,comp_pos as usize);
            }
        }
    }

    proptest! {
        #[test]
        #[cfg(target_feature = "bmi2")]
        fn select1_in_u64_slow_eq_fast(data: u64) {
            let bv = BitVec::<_, Lsb0>::from_element(data);
            for (found, _one_pos) in bv.iter_ones().enumerate() {
                let comp_pos_fast = unsafe { super::select1_in_u64_bmi2(data,found as u64) };
                let comp_pos_slow = super::select1_in_u64_slow(data,found as u64);
                assert_eq!(comp_pos_fast,comp_pos_slow);
            }
        }
    }
}