|
1 | 1 | module openzeppelin_math::vector; |
2 | 2 |
|
| 3 | +use openzeppelin_math::macros; |
| 4 | +use openzeppelin_math::rounding::RoundingMode; |
| 5 | + |
| 6 | +// === Errors === |
| 7 | + |
| 8 | +#[error(code = 0)] |
| 9 | +const EMedianOfEmptyVector: vector<u8> = "Median of empty vector is undefined"; |
| 10 | + |
3 | 11 | // === Public Functions === |
4 | 12 |
|
5 | 13 | /// Sort an unsigned integer vector in-place using the quicksort algorithm. |
@@ -194,3 +202,199 @@ public macro fun quick_sort_by<$T>($vec: &mut vector<$T>, $le: |&$T, &$T| -> boo |
194 | 202 | }; |
195 | 203 | }; |
196 | 204 | } |
| 205 | + |
| 206 | +/// Compute the median of an unsigned integer vector. |
| 207 | +/// |
| 208 | +/// For odd-length vectors, the median is the central order statistic and |
| 209 | +/// `$rounding_mode` has no effect on the result. |
| 210 | +/// |
| 211 | +/// For even-length vectors, the median is the arithmetic mean of the two central |
| 212 | +/// values, rounded according to `$rounding_mode`. |
| 213 | +/// |
| 214 | +/// The input vector is not mutated. |
| 215 | +/// |
| 216 | +/// Median selection uses quickselect over a `u256` working vector instead of |
| 217 | +/// sorting the full input. Benchmarks show similar cost for tiny vectors and |
| 218 | +/// substantially lower cost than sorting as input size grows. |
| 219 | +/// |
| 220 | +/// #### Generics |
| 221 | +/// - `$Int`: Any unsigned integer type (`u8`, `u16`, `u32`, `u64`, `u128`, or `u256`). |
| 222 | +/// |
| 223 | +/// #### Parameters |
| 224 | +/// - `$vec`: Reference to the vector whose median is desired. |
| 225 | +/// - `$rounding_mode`: Rounding strategy, applied only when the length is even. |
| 226 | +/// |
| 227 | +/// #### Returns |
| 228 | +/// - The median of `$vec`. |
| 229 | +/// |
| 230 | +/// #### Aborts |
| 231 | +/// - `EMedianOfEmptyVector` if `$vec` is empty. |
| 232 | +/// |
| 233 | +/// #### Example |
| 234 | +/// ```move |
| 235 | +/// use openzeppelin_math::rounding |
| 236 | +/// use openzeppelin_math::vector; |
| 237 | +/// |
| 238 | +/// let v = vector[3u64, 1, 4, 1, 5, 9, 2, 6]; |
| 239 | +/// let m = vector::median!(&v, rounding::down()); |
| 240 | +/// // m == 3 |
| 241 | +/// ``` |
| 242 | +public macro fun median<$Int>($vec: &vector<$Int>, $rounding_mode: RoundingMode): $Int { |
| 243 | + let vec = $vec; |
| 244 | + let vec_u256 = vec.map_ref!(|x| *x as u256); |
| 245 | + median_u256(vec_u256, $rounding_mode) as $Int |
| 246 | +} |
| 247 | + |
| 248 | +// === Concrete Public Functions === |
| 249 | + |
| 250 | +/// Compute the median of a `u256` vector with configurable rounding. |
| 251 | +/// |
| 252 | +/// This function consumes and partially reorders `vec`. |
| 253 | +/// |
| 254 | +/// Median selection uses quickselect instead of sorting the full input. |
| 255 | +/// Benchmarks show similar cost for tiny vectors and substantially lower cost |
| 256 | +/// than sorting as input size grows. |
| 257 | +/// |
| 258 | +/// #### Parameters |
| 259 | +/// - `vec`: Vector whose median is desired (consumed). |
| 260 | +/// - `rounding_mode`: Rounding strategy, applied only when the length is even. |
| 261 | +/// |
| 262 | +/// #### Returns |
| 263 | +/// - The median of `vec`. |
| 264 | +/// |
| 265 | +/// #### Aborts |
| 266 | +/// - `EMedianOfEmptyVector` if `vec` is empty. |
| 267 | +/// |
| 268 | +/// #### Example |
| 269 | +/// ```move |
| 270 | +/// use openzeppelin_math::rounding |
| 271 | +/// use openzeppelin_math::vector; |
| 272 | +/// |
| 273 | +/// let m = vector::median_u256(vector[10u256, 2, 8, 4], rounding::down()); |
| 274 | +/// // m == 6 |
| 275 | +/// ``` |
| 276 | +public fun median_u256(mut vec: vector<u256>, rounding_mode: RoundingMode): u256 { |
| 277 | + let len = vec.length(); |
| 278 | + assert!(len != 0, EMedianOfEmptyVector); |
| 279 | + |
| 280 | + let mid = len / 2; |
| 281 | + if (len % 2 == 1) { |
| 282 | + select_k_u256(&mut vec, mid) |
| 283 | + } else { |
| 284 | + let upper = select_k_u256(&mut vec, mid); |
| 285 | + // Selecting index `mid` leaves all preceding elements less than or equal |
| 286 | + // to `upper`, so the lower central statistic is the prefix maximum. |
| 287 | + let lower = max_prefix_u256(&vec, mid); |
| 288 | + macros::average!(lower, upper, rounding_mode) |
| 289 | + } |
| 290 | +} |
| 291 | + |
| 292 | +// === Internal Functions === |
| 293 | + |
| 294 | +// Return the kth order statistic, partially partitioning `vec` in place. |
| 295 | +// |
| 296 | +// Quickselect repeatedly partitions the active range around a pivot and keeps |
| 297 | +// only the side that can contain index `k`. This avoids sorting the whole vector |
| 298 | +// while still placing the selected value at the same index it would occupy in |
| 299 | +// sorted order. |
| 300 | +fun select_k_u256(vec: &mut vector<u256>, k: u64): u256 { |
| 301 | + let mut start = 0; |
| 302 | + let mut end = vec.length(); |
| 303 | + |
| 304 | + loop { |
| 305 | + if (start + 1 >= end) return vec[k]; |
| 306 | + |
| 307 | + // Use insertion sort as the quickselect base case for small active ranges. |
| 308 | + if (end - start <= 10) { |
| 309 | + insertion_sort_u256_range(vec, start, end); |
| 310 | + return vec[k] |
| 311 | + }; |
| 312 | + |
| 313 | + let pivot_index = end - 1; |
| 314 | + let mid = start + (end - start) / 2; |
| 315 | + |
| 316 | + // Median-of-three pivot selection using start, middle, and end - 1. |
| 317 | + // The median of those three values is moved into `pivot_index`. |
| 318 | + if (vec[mid] <= vec[start]) { |
| 319 | + vec.swap(start, mid); |
| 320 | + }; |
| 321 | + if (vec[pivot_index] <= vec[start]) { |
| 322 | + vec.swap(start, pivot_index) |
| 323 | + }; |
| 324 | + if (vec[mid] <= vec[pivot_index]) { |
| 325 | + vec.swap(mid, pivot_index); |
| 326 | + }; |
| 327 | + |
| 328 | + let mut lt = start; |
| 329 | + let mut i = start; |
| 330 | + let mut gt = pivot_index; |
| 331 | + |
| 332 | + // Three-way partition around the pivot at `pivot_index`. |
| 333 | + // |
| 334 | + // During the scan: |
| 335 | + // - [start, lt) contains values less than the pivot. |
| 336 | + // - [lt, i) contains values equal to the pivot. |
| 337 | + // - [i, gt) is unprocessed. |
| 338 | + // - [gt, pivot_index) contains values greater than the pivot. |
| 339 | + while (i < gt) { |
| 340 | + if (vec[i] <= vec[pivot_index]) { |
| 341 | + if (vec[pivot_index] <= vec[i]) { |
| 342 | + i = i + 1; |
| 343 | + } else { |
| 344 | + vec.swap(lt, i); |
| 345 | + lt = lt + 1; |
| 346 | + i = i + 1; |
| 347 | + } |
| 348 | + } else { |
| 349 | + gt = gt - 1; |
| 350 | + vec.swap(i, gt); |
| 351 | + }; |
| 352 | + }; |
| 353 | + |
| 354 | + // Move the pivot next to the equal region. After this, [lt, eq_end) |
| 355 | + // contains exactly the pivot-equivalent values. |
| 356 | + vec.swap(gt, pivot_index); |
| 357 | + let eq_end = gt + 1; |
| 358 | + |
| 359 | + // Discard the partitions that cannot contain the kth order statistic. |
| 360 | + if (k < lt) { |
| 361 | + end = lt; |
| 362 | + } else if (k >= eq_end) { |
| 363 | + start = eq_end; |
| 364 | + } else { |
| 365 | + return vec[k] |
| 366 | + }; |
| 367 | + } |
| 368 | +} |
| 369 | + |
| 370 | +// Sort a half-open sub-range [start, end). Used as the quickselect base case |
| 371 | +// for small active ranges. |
| 372 | +fun insertion_sort_u256_range(vec: &mut vector<u256>, start: u64, end: u64) { |
| 373 | + let mut i = start + 1; |
| 374 | + while (i < end) { |
| 375 | + let mut j = i; |
| 376 | + while (j != start && vec[j] < vec[j - 1]) { |
| 377 | + vec.swap(j, j - 1); |
| 378 | + j = j - 1; |
| 379 | + }; |
| 380 | + i = i + 1; |
| 381 | + }; |
| 382 | +} |
| 383 | + |
| 384 | +// Return the maximum value in vec[0..end). |
| 385 | +// |
| 386 | +// For even-length medians, selecting the upper middle index partitions the |
| 387 | +// vector enough to guarantee every earlier element is less than or equal to it, |
| 388 | +// but that prefix is not sorted. The lower middle value is therefore the maximum |
| 389 | +// of the prefix. |
| 390 | +fun max_prefix_u256(vec: &vector<u256>, end: u64): u256 { |
| 391 | + let mut max = vec[0]; |
| 392 | + let mut i = 1; |
| 393 | + while (i < end) { |
| 394 | + if (max < vec[i]) { |
| 395 | + max = vec[i]; |
| 396 | + }; |
| 397 | + i = i + 1; |
| 398 | + }; |
| 399 | + max |
| 400 | +} |
0 commit comments