Struct WaveletMatrix 

pub struct WaveletMatrix<T> {
    len: usize,
    bit_length: usize,
    zeros: Vec<usize>,
    bit_vectors: Vec<BitVector>,
    quad_vectors: Vec<WaveletMatrixQuadVector>,
    compress: VecCompress<T>,
    backend: SimdBackend,
}

Fields

len: usize

bit_length: usize

zeros: Vec<usize>

bit_vectors: Vec<BitVector>

quad_vectors: Vec<WaveletMatrixQuadVector>

compress: VecCompress<T>

backend: SimdBackend

Implementations

impl<T> WaveletMatrix<T>

where T: Ord + Clone,

pub fn new(v: Vec<T>) -> Self

Examples found in repository

crates/library_checker/src/data_structure/static_range_frequency.rs (line 20)

pub fn static_range_frequency_wavelet_matrix(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, a: [usize; n]);
    let wm = WaveletMatrix::new(a);
    for _ in 0..q {
        sc!(l, r, x: usize);
        let ans = wm.rank(x, l..r);
        pp!(ans);
    }
}

crates/library_checker/src/data_structure/range_kth_smallest.rs (line 8)

pub fn range_kth_smallest(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, a: [usize; n], queries: [(usize, usize, usize); iter q]);
    let wm = WaveletMatrix::new(a);
    let results = wm.quantile_batch(queries.map(|(l, r, k)| (l..r, k)));
    pp!(@lf @it results);
}

crates/library_checker/src/data_structure/static_range_sum_with_upper_bound.rs (line 26)

pub fn static_range_sum_with_upper_bound_wavelet_matrix(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, a: [i64; n]);
    let weights = a.clone();
    let wm = WaveletMatrix::new(a);
    let fold = wm.build_fold::<AdditiveOperation<i64>>(&weights);
    for _ in 0..q {
        sc!(l, r, x: i64);
        let (count, sum) = fold.fold_lessthan_with_count(x + 1, l..r);
        pp!(count, sum);
    }
}

crates/library_checker/src/data_structure/rectangle_sum.rs (line 20)

pub fn rectangle_sum(reader: impl Read, writer: impl Write) {
    prepare_io!(buffered; reader, writer);
    sc!(n, q, mut xyw: [(u32, u32, i64); n], queries: [(u32, u32, u32, u32); q]);
    xyw.radix_sort_by_key(|&(x, ..)| x);
    let xs: Vec<_> = xyw.iter().map(|&(x, ..)| x).collect();
    let search = StaticSearch::from_sorted(&xs);
    let endpoints: Vec<_> = queries.iter().flat_map(|&(l, _, r, _)| [l, r]).collect();
    let mut positions = vec![0; endpoints.len()];
    search.lower_bound_batch(&endpoints, &mut positions);
    let ys = xyw.iter().map(|&(_, y, _)| y).collect();
    let weights: Vec<_> = xyw.iter().map(|&(_, _, w)| w).collect();
    let wm = WaveletMatrix::new(ys);
    let fold = wm.build_fold::<AdditiveOperation<i64>>(&weights);
    let result = fold.fold_lessthan_batch(
        queries
            .into_iter()
            .zip(positions.as_chunks::<2>().0)
            .flat_map(|((_, d, _, u), &[l, r])| [(d, l..r), (u, l..r)]),
    );
    for &[lower, upper] in result.as_chunks::<2>().0 {
        pp!(upper - lower);
    }
}

crates/library_checker/src/data_structure/point_add_rectangle_sum.rs (line 39)

pub fn point_add_rectangle_sum(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, xyw: [(u32, u32, u64); iter n]);
    let mut points: Vec<_> = xyw.map(|(x, y, w)| (x, y, w as i64)).collect();
    sc!(queries: [Query; q]);
    points.extend(queries.iter().filter_map(|&query| match query {
        Query::Add { x, y, .. } => Some((x, y, 0)),
        Query::Sum { .. } => None,
    }));
    let mut order: Vec<_> = (0..points.len()).collect();
    order.radix_sort_by_key(|&i| points[i].0);
    let mut positions = vec![0; points.len()];
    let mut xs = Vec::with_capacity(points.len());
    let mut ys = Vec::with_capacity(points.len());
    let mut weights = Vec::with_capacity(points.len());
    for (i, &point) in order.iter().enumerate() {
        positions[point] = i;
        let (x, y, w) = points[point];
        xs.push(x);
        ys.push(y);
        weights.push(w);
    }
    let wm = WaveletMatrix::new(ys);
    let mut fold: WaveletMatrixPointAdd<_, AdditiveOperation<i64>> = wm.build_point_add(&weights);

    let mut point = n;
    for query in queries {
        match query {
            Query::Add { w, .. } => {
                fold.update(positions[point], w as i64);
                point += 1;
            }
            Query::Sum { l, d, r, u } => {
                let l = xs.partition_point(|&x| x < l);
                let r = xs.partition_point(|&x| x < r);
                pp!(fold.fold_range(d..u, l..r));
            }
        }
    }
}

crates/competitive/src/data_structure/wavelet_matrix.rs (line 597)

    pub fn new_with_init<F>(v: Vec<T>, mut f: F) -> Self
    where
        F: FnMut(usize, usize, T),
    {
        let this = Self::new(v.clone());
        if !this.quad_vectors.is_empty() {
            let bits = usize::BITS as usize
                - this.compress.size().saturating_sub(1).leading_zeros() as usize;
            for (mut k, value) in v.into_iter().enumerate() {
                if this.bit_length > bits {
                    f(this.bit_length - 1, k, value.clone());
                }
                for (level, vector) in this.quad_vectors.iter().enumerate() {
                    let d = (this.quad_vectors.len() - level - 1) * 2;
                    let (digit, rank) = vector.access_rank(k);
                    if d + 1 < bits {
                        let block = &vector.blocks[k / 64];
                        let high = block.rank[2] as usize
                            + block.rank[3] as usize
                            + (block.hi & !(u64::MAX << (k % 64))).count_ones() as usize;
                        let middle = if digit & 2 == 0 {
                            k - high
                        } else {
                            this.zeros[this.level(d + 1)] + high
                        };
                        f(d + 1, middle, value.clone());
                    }
                    k = vector.starts[digit] + rank;
                    f(d, k, value.clone());
                }
            }
            return this;
        }
        for (mut k, value) in v.into_iter().enumerate() {
            for d in (0..this.bit_length).rev() {
                let level = this.level(d);
                let (bit, rank1) = this.bit_vectors[level].access_rank1(k);
                k = if bit {
                    this.zeros[level] + rank1
                } else {
                    k - rank1
                };
                f(d, k, value.clone());
            }
        }
        this
    }

fn pack_words<I: Copy>(indices: &[I], bit: impl Fn(I) -> bool) -> Vec<u64>

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 469)

    pub fn new(v: Vec<T>) -> Self {
        if v.len() <= u32::MAX as usize {
            #[cfg(target_arch = "x86_64")]
            let backend = super::simd_backend();
            Self::from_values(
                v,
                |i| i as u32,
                |i| i as usize,
                |indices, d| {
                    #[cfg(target_arch = "x86_64")]
                    if backend != super::SimdBackend::Scalar && is_x86_feature_detected!("avx2") {
                        // SAFETY: AVX2 is available.
                        return unsafe { simd::pack_words(indices, d) };
                    }
                    Self::pack_words(indices, |i| (i >> d) & 1 != 0)
                },
                |indices, words, zeros, next| {
                    #[cfg(target_arch = "x86_64")]
                    if backend != super::SimdBackend::Scalar && is_x86_feature_detected!("avx2") {
                        // SAFETY: AVX2 is available, and the partition buffers have equal length.
                        unsafe { simd::partition_avx2(indices, words, zeros, next) };
                        return;
                    }
                    Self::partition(indices, words, zeros, next);
                },
            )
        } else {
            Self::from_values(
                v,
                |i| i,
                |i| i,
                |indices, d| Self::pack_words(indices, |i| (i >> d) & 1 != 0),
                Self::partition,
            )
        }
    }

fn partition<I: Copy>(indices: &[I], words: &[u64], one: usize, next: &mut [I])

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 478)

    pub fn new(v: Vec<T>) -> Self {
        if v.len() <= u32::MAX as usize {
            #[cfg(target_arch = "x86_64")]
            let backend = super::simd_backend();
            Self::from_values(
                v,
                |i| i as u32,
                |i| i as usize,
                |indices, d| {
                    #[cfg(target_arch = "x86_64")]
                    if backend != super::SimdBackend::Scalar && is_x86_feature_detected!("avx2") {
                        // SAFETY: AVX2 is available.
                        return unsafe { simd::pack_words(indices, d) };
                    }
                    Self::pack_words(indices, |i| (i >> d) & 1 != 0)
                },
                |indices, words, zeros, next| {
                    #[cfg(target_arch = "x86_64")]
                    if backend != super::SimdBackend::Scalar && is_x86_feature_detected!("avx2") {
                        // SAFETY: AVX2 is available, and the partition buffers have equal length.
                        unsafe { simd::partition_avx2(indices, words, zeros, next) };
                        return;
                    }
                    Self::partition(indices, words, zeros, next);
                },
            )
        } else {
            Self::from_values(
                v,
                |i| i,
                |i| i,
                |indices, d| Self::pack_words(indices, |i| (i >> d) & 1 != 0),
                Self::partition,
            )
        }
    }

fn from_values<I: Copy>( v: Vec<T>, code: impl Fn(usize) -> I, index: impl Fn(I) -> usize, pack: impl Fn(&[I], usize) -> Vec<u64>, partition: impl Fn(&[I], &[u64], usize, &mut [I]), ) -> Self

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (lines 459-480)

    pub fn new(v: Vec<T>) -> Self {
        if v.len() <= u32::MAX as usize {
            #[cfg(target_arch = "x86_64")]
            let backend = super::simd_backend();
            Self::from_values(
                v,
                |i| i as u32,
                |i| i as usize,
                |indices, d| {
                    #[cfg(target_arch = "x86_64")]
                    if backend != super::SimdBackend::Scalar && is_x86_feature_detected!("avx2") {
                        // SAFETY: AVX2 is available.
                        return unsafe { simd::pack_words(indices, d) };
                    }
                    Self::pack_words(indices, |i| (i >> d) & 1 != 0)
                },
                |indices, words, zeros, next| {
                    #[cfg(target_arch = "x86_64")]
                    if backend != super::SimdBackend::Scalar && is_x86_feature_detected!("avx2") {
                        // SAFETY: AVX2 is available, and the partition buffers have equal length.
                        unsafe { simd::partition_avx2(indices, words, zeros, next) };
                        return;
                    }
                    Self::partition(indices, words, zeros, next);
                },
            )
        } else {
            Self::from_values(
                v,
                |i| i,
                |i| i,
                |indices, d| Self::pack_words(indices, |i| (i >> d) & 1 != 0),
                Self::partition,
            )
        }
    }

pub fn new_with_init<F>(v: Vec<T>, f: F) -> Self

where F: FnMut(usize, usize, T),

fn level(&self, d: usize) -> usize

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 616)

    pub fn new_with_init<F>(v: Vec<T>, mut f: F) -> Self
    where
        F: FnMut(usize, usize, T),
    {
        let this = Self::new(v.clone());
        if !this.quad_vectors.is_empty() {
            let bits = usize::BITS as usize
                - this.compress.size().saturating_sub(1).leading_zeros() as usize;
            for (mut k, value) in v.into_iter().enumerate() {
                if this.bit_length > bits {
                    f(this.bit_length - 1, k, value.clone());
                }
                for (level, vector) in this.quad_vectors.iter().enumerate() {
                    let d = (this.quad_vectors.len() - level - 1) * 2;
                    let (digit, rank) = vector.access_rank(k);
                    if d + 1 < bits {
                        let block = &vector.blocks[k / 64];
                        let high = block.rank[2] as usize
                            + block.rank[3] as usize
                            + (block.hi & !(u64::MAX << (k % 64))).count_ones() as usize;
                        let middle = if digit & 2 == 0 {
                            k - high
                        } else {
                            this.zeros[this.level(d + 1)] + high
                        };
                        f(d + 1, middle, value.clone());
                    }
                    k = vector.starts[digit] + rank;
                    f(d, k, value.clone());
                }
            }
            return this;
        }
        for (mut k, value) in v.into_iter().enumerate() {
            for d in (0..this.bit_length).rev() {
                let level = this.level(d);
                let (bit, rank1) = this.bit_vectors[level].access_rank1(k);
                k = if bit {
                    this.zeros[level] + rank1
                } else {
                    k - rank1
                };
                f(d, k, value.clone());
            }
        }
        this
    }

    fn level(&self, d: usize) -> usize {
        self.bit_length - 1 - d
    }

    fn rank1(&self, level: usize, k: usize) -> usize {
        self.bit_vectors[level].rank1(k)
    }

    fn reorder<U>(&self, level: usize, current: Vec<U>) -> Vec<U> {
        assert_eq!(current.len(), self.len);
        let mut next = Vec::with_capacity(self.len);
        next.resize_with(self.len, MaybeUninit::uninit);
        let mut zero = 0;
        let mut one = self.zeros[level];
        let mut current = current.into_iter();
        for block in self.bit_vectors[level].blocks() {
            let count = current.len().min(64);
            if block.bits == 0 || block.bits == u64::MAX {
                let offset = if block.bits == 0 { &mut zero } else { &mut one };
                for (slot, value) in next[*offset..*offset + count]
                    .iter_mut()
                    .zip(current.by_ref().take(count))
                {
                    slot.write(value);
                }
                *offset += count;
            } else {
                for (i, value) in current.by_ref().take(count).enumerate() {
                    let bit = (block.bits >> i) & 1 != 0;
                    next[if bit { one } else { zero }].write(value);
                    zero += !bit as usize;
                    one += bit as usize;
                }
            }
        }
        // SAFETY: the partition counts fill every slot once, and `MaybeUninit<U>` has `U`'s layout.
        unsafe {
            let mut next = mem::ManuallyDrop::new(next);
            Vec::from_raw_parts(next.as_mut_ptr().cast(), next.len(), next.capacity())
        }
    }

    fn range_by_index(&self, idx: usize, mut range: Range<usize>) -> Range<usize> {
        if !self.quad_vectors.is_empty() {
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                if range.is_empty() {
                    break;
                }
                let digit = (idx >> ((self.quad_vectors.len() - level - 1) * 2)) & 3;
                range = vector.starts[digit] + vector.rank(digit, range.start)
                    ..vector.starts[digit] + vector.rank(digit, range.end);
            }
            return range;
        }
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            if range.is_empty() {
                break;
            }
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            if ((idx >> d) & 1) != 0 {
                range.start = self.zeros[level] + start1;
                range.end = self.zeros[level] + end1;
            } else {
                range.start -= start1;
                range.end -= end1;
            }
        }
        range
    }

    fn batch<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        #[cfg(target_arch = "x86_64")]
        if !cfg!(target_feature = "popcnt") && is_x86_feature_detected!("popcnt") {
            // SAFETY: POPCNT is checked above.
            unsafe {
                self.batch_popcnt::<OP>(states, count);
            }
            return;
        }
        self.batch_inner::<OP>(states, count);
    }

    #[cfg(target_arch = "x86_64")]
    #[target_feature(enable = "popcnt")]
    unsafe fn batch_popcnt<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        self.batch_inner::<OP>(states, count);
    }

    #[inline(always)]
    fn batch_inner<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        if OP == RANK_LESSTHAN && self.compress.size() <= 1 {
            for state in &mut states[..count] {
                state[3] = if state[2] == 0 {
                    0
                } else {
                    state[1] - state[0]
                };
            }
            return;
        }
        if !self.quad_vectors.is_empty() && matches!(OP, ACCESS | RANK | QUANTILE) {
            #[cfg(target_arch = "x86_64")]
            if OP == QUANTILE && count >= 8 && self.backend == super::SimdBackend::Avx512 {
                // SAFETY: checked batch ranges stay within the quad blocks after each
                // stable partition. Construction restricts quad counters to u32, and
                // the cached backend includes AVX512F and VPOPCNTDQ.
                unsafe {
                    simd::quad_avx512(&self.quad_vectors, &mut states[..count.next_multiple_of(8)]);
                }
                return;
            }

            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                #[cfg(target_arch = "x86_64")]
                let next = self
                    .quad_vectors
                    .get(level + 1)
                    .filter(|_| self.len >= 1048576 || (OP == ACCESS && self.len >= 262144));
                for state in &mut states[..count] {
                    if OP == ACCESS {
                        let (digit, rank) = vector.access_rank(state[0]);
                        state[0] = vector.starts[digit] + rank;
                        state[3] = state[3] * 4 + digit;
                    } else if OP == RANK {
                        let digit = (state[2] >> d) & 3;
                        state[0] = vector.starts[digit] + vector.rank(digit, state[0]);
                        state[1] = vector.starts[digit] + vector.rank(digit, state[1]);
                    } else {
                        let start = vector.ranks(state[0]);
                        let end = vector.ranks(state[1]);
                        let prefix = [
                            0,
                            end[0] - start[0],
                            end[0] + end[1] - start[0] - start[1],
                            state[1] - state[0] - (end[3] - start[3]),
                        ];
                        let digit = (state[2] >= prefix[1]) as usize
                            + (state[2] >= prefix[2]) as usize
                            + (state[2] >= prefix[3]) as usize;
                        state[2] -= prefix[digit];
                        state[0] = vector.starts[digit] + start[digit];
                        state[1] = vector.starts[digit] + end[digit];
                        state[3] = state[3] * 4 + digit;
                    }
                    #[cfg(target_arch = "x86_64")]
                    if let Some(next) = next {
                        // SAFETY: stable partitions keep both endpoints within the next vector.
                        unsafe {
                            std::arch::x86_64::_mm_prefetch::<{ std::arch::x86_64::_MM_HINT_T0 }>(
                                next.blocks.as_ptr().add(state[0] / 64).cast(),
                            );
                            if OP != ACCESS {
                                std::arch::x86_64::_mm_prefetch::<{ std::arch::x86_64::_MM_HINT_T0 }>(
                                    next.blocks.as_ptr().add(state[1] / 64).cast(),
                                );
                            }
                        }
                    }
                }
            }
            return;
        }
        let first = (matches!(OP, ACCESS | RANK | QUANTILE)
            && self.compress.size().is_power_of_two()) as usize;
        #[cfg(target_arch = "x86_64")]
        if OP == RANK_LESSTHAN && count >= 8 && self.len <= u32::MAX as usize {
            // SAFETY: the caller checked each initial position. Rank transitions remain in
            // 0..=len, including the sentinel block. Padding uses position zero. The
            // backend was detected at construction, and x86-64 blocks contain two u64s.
            unsafe {
                match self.backend {
                    super::SimdBackend::Avx512 => {
                        simd::rank_lessthan_avx512(
                            &self.bit_vectors[first..],
                            &self.zeros[first..],
                            &mut states[..count.next_multiple_of(8)],
                        );
                        return;
                    }
                    super::SimdBackend::Avx2 if self.len < 1048576 => {
                        simd::rank_lessthan_avx2(
                            &self.bit_vectors[first..],
                            &self.zeros[first..],
                            &mut states[..count.next_multiple_of(4)],
                        );
                        return;
                    }
                    _ => {}
                }
            }
        }
        for d in (0..self.bit_length - first).rev() {
            let level = self.level(d);
            for state in &mut states[..count] {
                let (bit, start1) = self.bit_vectors[level].access_rank1(state[0]);
                let start0 = state[0] - start1;
                let end1 = if OP == ACCESS {
                    0
                } else {
                    self.rank1(level, state[1])
                };
                let end0 = if OP == ACCESS { 0 } else { state[1] - end1 };
                let count0 = if OP == ACCESS { 0 } else { end0 - start0 };
                let bit = match OP {
                    ACCESS => bit,
                    RANK | RANK_LESSTHAN => (state[2] >> d) & 1 != 0,
                    _ => state[2] >= count0,
                };
                state[0] = if bit {
                    self.zeros[level] + start1
                } else {
                    start0
                };
                if OP != ACCESS {
                    state[1] = if bit { self.zeros[level] + end1 } else { end0 };
                }
                if OP == ACCESS || OP == QUANTILE {
                    state[3] |= (bit as usize) << d;
                }
                if OP == QUANTILE {
                    state[2] -= if bit { count0 } else { 0 };
                }
                if OP == RANK_LESSTHAN {
                    state[3] += if bit { count0 } else { 0 };
                }
            }
        }
    }

    /// get k-th value
    pub fn access(&self, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            let mut index = 0;
            for vector in &self.quad_vectors {
                let (digit, rank) = vector.access_rank(k);
                index = index * 4 + digit;
                k = vector.starts[digit] + rank;
            }
            return self.compress.values()[index].clone();
        }
        let mut idx = 0;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let (bit, rank1) = self.bit_vectors[level].access_rank1(k);
            idx |= (bit as usize) << d;
            k = if bit {
                self.zeros[level] + rank1
            } else {
                k - rank1
            };
        }
        self.compress.values()[idx].clone()
    }

    /// Returns the values at `indices` in input order.
    pub fn access_batch(&self, indices: impl IntoIterator<Item = usize>) -> Vec<T> {
        let indices: Vec<_> = indices.into_iter().collect();
        let mut result = Vec::with_capacity(indices.len());
        for indices in indices.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, &index) in states.iter_mut().zip(indices) {
                assert!(index < self.len);
                state[0] = index;
            }
            self.batch::<ACCESS>(&mut states, indices.len());
            result.extend(
                states[..indices.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// the number of val in range
    pub fn rank(&self, val: T, range: Range<usize>) -> usize {
        match self.compress.index_exact(&val) {
            Some(idx) => self.range_by_index(idx, range).len(),
            None => 0,
        }
    }

    /// Returns the number of exact matches for each `(value, range)` query.
    pub fn rank_batch(&self, queries: impl IntoIterator<Item = (T, Range<usize>)>) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                if let Some(index) = self.compress.index_exact(value) {
                    *state = [range.start, range.end, index, 0];
                }
            }
            self.batch::<RANK>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| state[1] - state[0]),
            );
        }
        result
    }

    /// index of k-th val
    pub fn select(&self, val: T, k: usize) -> Option<usize> {
        let idx = self.compress.index_exact(&val)?;
        let range = self.range_by_index(idx, 0..self.len);
        if range.len() <= k {
            return None;
        }
        let mut i = range.start + k;
        if !self.quad_vectors.is_empty() {
            for (level, vector) in self.quad_vectors.iter().enumerate().rev() {
                let digit = (idx >> ((self.quad_vectors.len() - level - 1) * 2)) & 3;
                i = vector.select(digit, i - vector.starts[digit]);
            }
            return Some(i);
        }
        for level in (self.compress.size().is_power_of_two() as usize..self.bit_length).rev() {
            if i >= self.zeros[level] {
                i = self.bit_vectors[level]
                    .select1(i - self.zeros[level])
                    .unwrap();
            } else {
                i = self.bit_vectors[level].select0(i).unwrap();
            }
        }
        Some(i)
    }

    /// get k-th smallest value in range
    pub fn quantile(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            return self.quad_quantile(range, k, 0, 0);
        }
        let mut idx = 0;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let z = end0 - start0;
            let bit = z <= k;
            k -= if bit { z } else { 0 };
            idx |= (bit as usize) << d;
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        self.compress.values()[idx].clone()
    }

    #[inline(always)]
    fn quad_quantile(
        &self,
        mut range: Range<usize>,
        mut k: usize,
        level: usize,
        mut index: usize,
    ) -> T {
        for vector in &self.quad_vectors[level..] {
            let start = vector.ranks(range.start);
            let end = vector.ranks(range.end);
            let mut digit = 0;
            while digit < 3 && k >= end[digit] - start[digit] {
                k -= end[digit] - start[digit];
                digit += 1;
            }
            index = index * 4 + digit;
            range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
        }
        self.compress.values()[index].clone()
    }

    pub fn quantile_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<usize>, usize)>,
    ) -> Vec<T> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (range, k)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len && *k < range.len());
                *state = [range.start, range.end, *k, 0];
            }
            if let [(range, k)] = queries {
                result.push(self.quantile(range.clone(), *k));
                continue;
            }
            self.batch::<QUANTILE>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// Returns the requested order statistics of one range. `ranks` must be sorted
    /// in nondecreasing order and every rank must be less than the range length.
    pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T> {
        assert!(range.start <= range.end && range.end <= self.len);
        assert!(ranks.windows(2).all(|pair| pair[0] <= pair[1]));
        assert!(ranks.last().is_none_or(|&k| k < range.len()));
        if ranks.is_empty() {
            return Vec::new();
        }
        if ranks[0] == ranks[ranks.len() - 1] {
            return vec![self.quantile(range, ranks[0]); ranks.len()];
        }
        let use_batch = ranks.len() < 1024;
        #[cfg(target_arch = "x86_64")]
        let use_batch = use_batch || self.backend == super::SimdBackend::Avx512;
        // Sparse requests amortize less of the shared traversal's branching work.
        if use_batch && ranks.len() < self.compress.size().div_ceil(16) {
            return self.quantile_batch(ranks.iter().map(|&k| (range.clone(), k)));
        }
        let mut result = Vec::with_capacity(ranks.len());
        if !self.quad_vectors.is_empty() {
            let mut stack = vec![(range, 0, 0, 0, ranks)];
            while let Some((range, base, index, level, ranks)) = stack.pop() {
                if level == self.quad_vectors.len() {
                    result.resize(
                        result.len() + ranks.len(),
                        self.compress.values()[index].clone(),
                    );
                    continue;
                }
                if ranks.len() == 1 {
                    result.push(self.quad_quantile(range, ranks[0] - base, level, index));
                    continue;
                }
                let vector = &self.quad_vectors[level];
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let mut split = base + range.len();
                let mut requested = ranks.len();
                for digit in (0..4).rev() {
                    split -= end[digit] - start[digit];
                    let boundary = ranks[..requested].partition_point(|&k| k < split);
                    if boundary < requested {
                        stack.push((
                            vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit],
                            split,
                            index * 4 + digit,
                            level + 1,
                            &ranks[boundary..requested],
                        ));
                    }
                    requested = boundary;
                }
            }
            return result;
        }
        let mut stack = vec![(range, 0, 0, self.bit_length, ranks)];
        while let Some((range, base, index, bits, ranks)) = stack.pop() {
            if bits == 0 {
                result.resize(
                    result.len() + ranks.len(),
                    self.compress.values()[index].clone(),
                );
                continue;
            }
            let d = bits - 1;
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let split = base + end0 - start0;
            let boundary = ranks.partition_point(|&k| k < split);
            if boundary < ranks.len() {
                stack.push((
                    self.zeros[level] + start1..self.zeros[level] + end1,
                    split,
                    index | (1 << d),
                    d,
                    &ranks[boundary..],
                ));
            }
            if boundary != 0 {
                stack.push((start0..end0, base, index, d, &ranks[..boundary]));
            }
        }
        result
    }

    /// get k-th smallest value out of range
    pub fn quantile_outer(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            let mut outer = 0..self.len;
            let mut index = 0;
            for vector in &self.quad_vectors {
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let outer_start = vector.ranks(outer.start);
                let outer_end = vector.ranks(outer.end);
                let mut digit = 0;
                while digit < 3 {
                    let count = outer_end[digit] - outer_start[digit] - (end[digit] - start[digit]);
                    if k < count {
                        break;
                    }
                    k -= count;
                    digit += 1;
                }
                index = index * 4 + digit;
                range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
                outer = vector.starts[digit] + outer_start[digit]
                    ..vector.starts[digit] + outer_end[digit];
            }
            return self.compress.values()[index].clone();
        }
        let mut idx = 0;
        let mut orange = 0..self.len;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let range_start1 = self.rank1(level, range.start);
            let range_end1 = self.rank1(level, range.end);
            let outer_start1 = self.rank1(level, orange.start);
            let outer_end1 = self.rank1(level, orange.end);
            let range_start0 = range.start - range_start1;
            let range_end0 = range.end - range_end1;
            let outer_start0 = orange.start - outer_start1;
            let outer_end0 = orange.end - outer_end1;
            let z = (outer_end0 - outer_start0) - (range_end0 - range_start0);
            if z <= k {
                k -= z;
                idx |= 1 << d;
                range.start = self.zeros[level] + range_start1;
                range.end = self.zeros[level] + range_end1;
                orange.start = self.zeros[level] + outer_start1;
                orange.end = self.zeros[level] + outer_end1;
            } else {
                range.start = range_start0;
                range.end = range_end0;
                orange.start = outer_start0;
                orange.end = outer_end0;
            }
        }
        self.compress.values()[idx].clone()
    }

    /// the number of value less than val in range
    pub fn rank_lessthan(&self, val: T, range: Range<usize>) -> usize {
        let index = self.compress.index_lower_bound(&val);
        if index == self.compress.size() {
            return range.len();
        }
        if self.len >= 262144 && !self.quad_vectors.is_empty() {
            self.quad_rank_lessthan_index(index, range, 0)
        } else {
            let bits = usize::BITS as usize
                - self.compress.size().saturating_sub(1).leading_zeros() as usize;
            self.rank_lessthan_index(index, range, bits)
        }
    }

    /// Returns the number of values below each query's threshold.
    pub fn rank_lessthan_batch(
        &self,
        queries: impl IntoIterator<Item = (T, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                *state = [
                    range.start,
                    range.end,
                    self.compress.index_lower_bound(value),
                    0,
                ];
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len());
            result.extend(states[..queries.len()].iter().map(|state| state[3]));
        }
        result
    }

    fn quad_rank_lessthan_index(
        &self,
        index: usize,
        mut range: Range<usize>,
        start_level: usize,
    ) -> usize {
        let mut result = 0;
        for (level, vector) in self.quad_vectors.iter().enumerate().skip(start_level) {
            let d = (self.quad_vectors.len() - level - 1) * 2;
            if d + 2 <= index.trailing_zeros() as usize {
                break;
            }
            let digit = (index >> d) & 3;
            let ranks = |position: usize| {
                let block = &vector.blocks[position / 64];
                let mask = !(u64::MAX << (position % 64));
                let low = if digit & 1 == 0 { !block.lo } else { block.lo };
                let high = if digit & 2 == 0 { !block.hi } else { block.hi };
                let exact = block.rank[digit] as usize + (low & high & mask).count_ones() as usize;
                let less = match digit {
                    0 => 0,
                    1 => {
                        block.rank[0] as usize
                            + (!(block.lo | block.hi) & mask).count_ones() as usize
                    }
                    2 => {
                        block.rank[0] as usize
                            + block.rank[1] as usize
                            + (!block.hi & mask).count_ones() as usize
                    }
                    _ => position - exact,
                };
                (exact, less)
            };
            let (start, start_less) = ranks(range.start);
            let (end, end_less) = ranks(range.end);
            result += end_less - start_less;
            range = vector.starts[digit] + start..vector.starts[digit] + end;
        }
        result
    }

    fn rank_lessthan_index(&self, idx: usize, mut range: Range<usize>, bits: usize) -> usize {
        let mut res = 0;
        for d in (idx.trailing_zeros() as usize..bits).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let bit = (idx >> d) & 1 != 0;
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            res += if bit { end0 - start0 } else { 0 };
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        res
    }

    /// the number of valrange in range
    pub fn rank_range(&self, valrange: Range<T>, mut range: Range<usize>) -> usize {
        let lower = self.compress.index_lower_bound(&valrange.start);
        let upper = self.compress.index_lower_bound(&valrange.end);
        if lower >= upper {
            return 0;
        }
        if !self.quad_vectors.is_empty() {
            if upper == self.compress.size() {
                return range.len() - self.quad_rank_lessthan_index(lower, range, 0);
            }
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let low = (lower >> d) & 3;
                let high = (upper >> d) & 3;
                if low == high {
                    range = vector.starts[low] + vector.rank(low, range.start)
                        ..vector.starts[low] + vector.rank(low, range.end);
                    continue;
                }
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let count: usize = (low..high).map(|digit| end[digit] - start[digit]).sum();
                return count
                    - self.quad_rank_lessthan_index(
                        lower,
                        vector.starts[low] + start[low]..vector.starts[low] + end[low],
                        level + 1,
                    )
                    + self.quad_rank_lessthan_index(
                        upper,
                        vector.starts[high] + start[high]..vector.starts[high] + end[high],
                        level + 1,
                    );
            }
            return 0;
        }
        for d in (0..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let zero_range = range.start - start1..range.end - end1;
            let one_range = self.zeros[level] + start1..self.zeros[level] + end1;
            if ((lower ^ upper) >> d) & 1 != 0 {
                return zero_range.len() - self.rank_lessthan_index(lower, zero_range, d)
                    + self.rank_lessthan_index(upper, one_range, d);
            }
            range = if (lower >> d) & 1 != 0 {
                one_range
            } else {
                zero_range
            };
        }
        0
    }

    /// Returns each query's count in a half-open value range.
    pub fn rank_range_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<T>, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(8) {
            let mut states = [[0; 4]; 16];
            for (i, (values, range)) in queries.iter().enumerate() {
                assert!(range.start <= range.end && range.end <= self.len);
                let lower = self.compress.index_lower_bound(&values.start);
                let upper = self.compress.index_lower_bound(&values.end);
                if lower < upper {
                    states[i * 2] = [range.start, range.end, lower, 0];
                    states[i * 2 + 1] = [range.start, range.end, upper, 0];
                }
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len() * 2);
            result.extend(
                states[..queries.len() * 2]
                    .as_chunks::<2>()
                    .0
                    .iter()
                    .map(|pair| pair[1][3] - pair[0][3]),
            );
        }
        result
    }

    pub fn query_less_than<F>(&self, val: T, mut range: Range<usize>, mut f: F)
    where
        F: FnMut(usize, Range<usize>),
    {
        let idx = self.compress.index_lower_bound(&val);
        if !self.quad_vectors.is_empty() {
            if idx == self.compress.size() && idx.is_power_of_two() {
                f(self.bit_length - 1, range);
                return;
            }
            for (level, vector) in self.quad_vectors.iter().enumerate().take(
                self.quad_vectors
                    .len()
                    .saturating_sub(idx.trailing_zeros() as usize / 2),
            ) {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let digit = (idx >> d) & 3;
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                if digit & 2 != 0 {
                    f(d + 1, start[0] + start[1]..end[0] + end[1]);
                }
                if digit & 1 != 0 {
                    let zero = digit & 2;
                    f(
                        d,
                        vector.starts[zero] + start[zero]..vector.starts[zero] + end[zero],
                    );
                }
                range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
            }
            return;
        }
        for d in (idx.trailing_zeros() as usize..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            if ((idx >> d) & 1) != 0 {
                f(d, start0..end0);
                range.start = self.zeros[level] + start1;
                range.end = self.zeros[level] + end1;
            } else {
                range.start = start0;
                range.end = end0;
            }
        }
    }

    pub fn build_fold<M>(&self, weights: &[M::T]) -> WaveletMatrixFold<'_, T, M>
    where
        M: AbelianGroup,
    {
        assert_eq!(weights.len(), self.len);
        let mut offsets = Vec::with_capacity(self.bit_length);
        let mut prefix = Vec::with_capacity(self.zeros.iter().map(|&zero| zero + 1).sum());
        let mut current: Vec<M::T> = weights.to_vec();
        for level in 0..self.bit_length {
            current = self.reorder(level, current);
            offsets.push(prefix.len());
            let mut acc = M::unit();
            prefix.push(acc.clone());
            for w in &current[..self.zeros[level]] {
                acc = M::operate(&acc, w);
                prefix.push(acc.clone());
            }
        }
        WaveletMatrixFold {
            wavelet_matrix: self,
            prefix,
            offsets,
        }
    }

    pub fn build_point_add<M>(&self, weights: &[M::T]) -> WaveletMatrixPointAdd<'_, T, M>
    where
        M: AbelianGroup,
    {
        assert_eq!(weights.len(), self.len);
        let mut current = weights.to_vec();
        let mut bits = Vec::with_capacity(self.bit_length);
        for level in 0..self.bit_length {
            current = self.reorder(level, current);
            bits.push(BinaryIndexedTree::from_slice(&current[..self.zeros[level]]));
        }
        WaveletMatrixPointAdd {
            wavelet_matrix: self,
            bits,
        }
    }
}

pub struct WaveletMatrixPointAdd<'a, T, M>
where
    T: Ord + Clone,
    M: AbelianGroup,
{
    wavelet_matrix: &'a WaveletMatrix<T>,
    bits: Vec<BinaryIndexedTree<M>>,
}

impl<'a, T, M> WaveletMatrixPointAdd<'a, T, M>
where
    T: Ord + Clone,
    M: AbelianGroup,
{
    pub fn update(&mut self, mut index: usize, value: M::T) {
        debug_assert!(index < self.wavelet_matrix.len);
        for d in (0..self.wavelet_matrix.bit_length).rev() {
            let level = self.wavelet_matrix.level(d);
            let (bit, rank1) = self.wavelet_matrix.bit_vectors[level].access_rank1(index);
            if bit {
                index = self.wavelet_matrix.zeros[level] + rank1;
            } else {
                index -= rank1;
                self.bits[level].update(index, value.clone());
            }
        }
    }

    pub fn fold_lessthan(&self, value: T, range: Range<usize>) -> M::T {
        let mut result = M::unit();
        self.wavelet_matrix
            .query_less_than(value, range, |d, range| {
                M::operate_assign(
                    &mut result,
                    &self.bits[self.wavelet_matrix.level(d)].fold_abelian(range.start, range.end),
                );
            });
        result
    }

    pub fn fold_range(&self, values: Range<T>, range: Range<usize>) -> M::T {
        let lower = self
            .wavelet_matrix
            .compress
            .index_lower_bound(&values.start);
        let upper = self.wavelet_matrix.compress.index_lower_bound(&values.end);
        if lower >= upper {
            return M::unit();
        }
        let mut range = range;
        for d in (0..self.wavelet_matrix.bit_length).rev() {
            let level = self.wavelet_matrix.level(d);
            let start1 = self.wavelet_matrix.bit_vectors[level].rank1(range.start);
            let end1 = self.wavelet_matrix.bit_vectors[level].rank1(range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            if ((lower >> d) & 1) == ((upper >> d) & 1) {
                if ((lower >> d) & 1) == 0 {
                    range = start0..end0;
                } else {
                    range = self.wavelet_matrix.zeros[level] + start1
                        ..self.wavelet_matrix.zeros[level] + end1;
                }
                continue;
            }
            let zero_range = start0..end0;
            let one_range =
                self.wavelet_matrix.zeros[level] + start1..self.wavelet_matrix.zeros[level] + end1;
            let lower_sum = self.fold_lessthan_index(lower, zero_range.clone(), d);
            let upper_sum = self.fold_lessthan_index(upper, one_range, d);
            let zero_sum = self.bits[level].fold_abelian(zero_range.start, zero_range.end);
            let mut result = M::rinv_operate(&zero_sum, &lower_sum);
            M::operate_assign(&mut result, &upper_sum);
            return result;
        }
        M::unit()
    }

    fn fold_lessthan_index(&self, idx: usize, mut range: Range<usize>, bits: usize) -> M::T {
        let mut result = M::unit();
        for d in (idx.trailing_zeros() as usize..bits).rev() {
            let level = self.wavelet_matrix.level(d);
            let start1 = self.wavelet_matrix.bit_vectors[level].rank1(range.start);
            let end1 = self.wavelet_matrix.bit_vectors[level].rank1(range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            if ((idx >> d) & 1) != 0 {
                M::operate_assign(&mut result, &self.bits[level].fold_abelian(start0, end0));
                range.start = self.wavelet_matrix.zeros[level] + start1;
                range.end = self.wavelet_matrix.zeros[level] + end1;
            } else {
                range.start = start0;
                range.end = end0;
            }
        }
        result
    }
}

#[derive(Debug, Clone)]
pub struct WaveletMatrixFold<'a, T, M>
where
    T: Ord + Clone,
    M: AbelianGroup,
{
    wavelet_matrix: &'a WaveletMatrix<T>,
    prefix: Vec<M::T>,
    offsets: Vec<usize>,
}

impl<'a, T, M> WaveletMatrixFold<'a, T, M>
where
    T: Ord + Clone,
    M: AbelianGroup,
{
    pub fn fold_lessthan(&self, val: T, range: Range<usize>) -> M::T {
        self.fold_lessthan_with_count(val, range).1
    }

    pub fn fold_lessthan_with_count(&self, val: T, range: Range<usize>) -> (usize, M::T) {
        debug_assert!(range.end <= self.wavelet_matrix.len);
        let [result] = self.fold_lessthan_indices_with_count(
            [self.wavelet_matrix.compress.index_lower_bound(&val)],
            [range],
            self.wavelet_matrix.bit_length,
        );
        result
    }

    pub fn fold_range(&self, valrange: Range<T>, range: Range<usize>) -> M::T {
        self.fold_range_with_count(valrange, range).1
    }

    pub fn fold_range_with_count(
        &self,
        valrange: Range<T>,
        mut range: Range<usize>,
    ) -> (usize, M::T) {
        debug_assert!(range.end <= self.wavelet_matrix.len);
        let lower = self
            .wavelet_matrix
            .compress
            .index_lower_bound(&valrange.start);
        let upper = self
            .wavelet_matrix
            .compress
            .index_lower_bound(&valrange.end);
        if lower >= upper {
            return (0, M::unit());
        }
        for d in (0..self.wavelet_matrix.bit_length).rev() {
            let level = self.wavelet_matrix.level(d);
            let start1 = self.wavelet_matrix.bit_vectors[level].rank1(range.start);
            let end1 = self.wavelet_matrix.bit_vectors[level].rank1(range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            if ((lower >> d) & 1) == ((upper >> d) & 1) {
                if ((lower >> d) & 1) == 0 {
                    range = start0..end0;
                } else {
                    range = self.wavelet_matrix.zeros[level] + start1
                        ..self.wavelet_matrix.zeros[level] + end1;
                }
                continue;
            }
            let zero_range = start0..end0;
            let one_range =
                self.wavelet_matrix.zeros[level] + start1..self.wavelet_matrix.zeros[level] + end1;
            let [(lower_count, lower_sum), (upper_count, upper_sum)] = self
                .fold_lessthan_indices_with_count(
                    [lower, upper],
                    [zero_range.clone(), one_range],
                    d,
                );
            let zero_sum = self.range_sum(level, zero_range.clone());
            return (
                zero_range.len() - lower_count + upper_count,
                M::operate(&M::rinv_operate(&zero_sum, &lower_sum), &upper_sum),
            );
        }
        (0, M::unit())
    }

    #[inline]
    fn range_sum(&self, level: usize, range: Range<usize>) -> M::T {
        let offset = self.offsets[level];
        M::rinv_operate(
            &self.prefix[offset + range.end],
            &self.prefix[offset + range.start],
        )
    }

    fn fold_lessthan_indices_with_count<const N: usize>(
        &self,
        indices: [usize; N],
        mut ranges: [Range<usize>; N],
        bits: usize,
    ) -> [(usize, M::T); N] {
        let mut results = std::array::from_fn(|_| (0, M::unit()));
        let last = indices
            .iter()
            .map(|index| index.trailing_zeros() as usize)
            .min()
            .unwrap_or(bits);
        for d in (last..bits).rev() {
            let level = self.wavelet_matrix.level(d);
            for ((&index, range), (count, sum)) in indices.iter().zip(&mut ranges).zip(&mut results)
            {
                let start1 = self.wavelet_matrix.bit_vectors[level].rank1(range.start);
                let end1 = self.wavelet_matrix.bit_vectors[level].rank1(range.end);
                let start0 = range.start - start1;
                let end0 = range.end - end1;
                if ((index >> d) & 1) != 0 {
                    *count += end0 - start0;
                    M::operate_assign(sum, &self.range_sum(level, start0..end0));
                    range.start = self.wavelet_matrix.zeros[level] + start1;
                    range.end = self.wavelet_matrix.zeros[level] + end1;
                } else {
                    range.start = start0;
                    range.end = end0;
                }
            }
        }
        results
    }

fn rank1(&self, level: usize, k: usize) -> usize

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 700)

    fn range_by_index(&self, idx: usize, mut range: Range<usize>) -> Range<usize> {
        if !self.quad_vectors.is_empty() {
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                if range.is_empty() {
                    break;
                }
                let digit = (idx >> ((self.quad_vectors.len() - level - 1) * 2)) & 3;
                range = vector.starts[digit] + vector.rank(digit, range.start)
                    ..vector.starts[digit] + vector.rank(digit, range.end);
            }
            return range;
        }
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            if range.is_empty() {
                break;
            }
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            if ((idx >> d) & 1) != 0 {
                range.start = self.zeros[level] + start1;
                range.end = self.zeros[level] + end1;
            } else {
                range.start -= start1;
                range.end -= end1;
            }
        }
        range
    }

    fn batch<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        #[cfg(target_arch = "x86_64")]
        if !cfg!(target_feature = "popcnt") && is_x86_feature_detected!("popcnt") {
            // SAFETY: POPCNT is checked above.
            unsafe {
                self.batch_popcnt::<OP>(states, count);
            }
            return;
        }
        self.batch_inner::<OP>(states, count);
    }

    #[cfg(target_arch = "x86_64")]
    #[target_feature(enable = "popcnt")]
    unsafe fn batch_popcnt<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        self.batch_inner::<OP>(states, count);
    }

    #[inline(always)]
    fn batch_inner<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        if OP == RANK_LESSTHAN && self.compress.size() <= 1 {
            for state in &mut states[..count] {
                state[3] = if state[2] == 0 {
                    0
                } else {
                    state[1] - state[0]
                };
            }
            return;
        }
        if !self.quad_vectors.is_empty() && matches!(OP, ACCESS | RANK | QUANTILE) {
            #[cfg(target_arch = "x86_64")]
            if OP == QUANTILE && count >= 8 && self.backend == super::SimdBackend::Avx512 {
                // SAFETY: checked batch ranges stay within the quad blocks after each
                // stable partition. Construction restricts quad counters to u32, and
                // the cached backend includes AVX512F and VPOPCNTDQ.
                unsafe {
                    simd::quad_avx512(&self.quad_vectors, &mut states[..count.next_multiple_of(8)]);
                }
                return;
            }

            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                #[cfg(target_arch = "x86_64")]
                let next = self
                    .quad_vectors
                    .get(level + 1)
                    .filter(|_| self.len >= 1048576 || (OP == ACCESS && self.len >= 262144));
                for state in &mut states[..count] {
                    if OP == ACCESS {
                        let (digit, rank) = vector.access_rank(state[0]);
                        state[0] = vector.starts[digit] + rank;
                        state[3] = state[3] * 4 + digit;
                    } else if OP == RANK {
                        let digit = (state[2] >> d) & 3;
                        state[0] = vector.starts[digit] + vector.rank(digit, state[0]);
                        state[1] = vector.starts[digit] + vector.rank(digit, state[1]);
                    } else {
                        let start = vector.ranks(state[0]);
                        let end = vector.ranks(state[1]);
                        let prefix = [
                            0,
                            end[0] - start[0],
                            end[0] + end[1] - start[0] - start[1],
                            state[1] - state[0] - (end[3] - start[3]),
                        ];
                        let digit = (state[2] >= prefix[1]) as usize
                            + (state[2] >= prefix[2]) as usize
                            + (state[2] >= prefix[3]) as usize;
                        state[2] -= prefix[digit];
                        state[0] = vector.starts[digit] + start[digit];
                        state[1] = vector.starts[digit] + end[digit];
                        state[3] = state[3] * 4 + digit;
                    }
                    #[cfg(target_arch = "x86_64")]
                    if let Some(next) = next {
                        // SAFETY: stable partitions keep both endpoints within the next vector.
                        unsafe {
                            std::arch::x86_64::_mm_prefetch::<{ std::arch::x86_64::_MM_HINT_T0 }>(
                                next.blocks.as_ptr().add(state[0] / 64).cast(),
                            );
                            if OP != ACCESS {
                                std::arch::x86_64::_mm_prefetch::<{ std::arch::x86_64::_MM_HINT_T0 }>(
                                    next.blocks.as_ptr().add(state[1] / 64).cast(),
                                );
                            }
                        }
                    }
                }
            }
            return;
        }
        let first = (matches!(OP, ACCESS | RANK | QUANTILE)
            && self.compress.size().is_power_of_two()) as usize;
        #[cfg(target_arch = "x86_64")]
        if OP == RANK_LESSTHAN && count >= 8 && self.len <= u32::MAX as usize {
            // SAFETY: the caller checked each initial position. Rank transitions remain in
            // 0..=len, including the sentinel block. Padding uses position zero. The
            // backend was detected at construction, and x86-64 blocks contain two u64s.
            unsafe {
                match self.backend {
                    super::SimdBackend::Avx512 => {
                        simd::rank_lessthan_avx512(
                            &self.bit_vectors[first..],
                            &self.zeros[first..],
                            &mut states[..count.next_multiple_of(8)],
                        );
                        return;
                    }
                    super::SimdBackend::Avx2 if self.len < 1048576 => {
                        simd::rank_lessthan_avx2(
                            &self.bit_vectors[first..],
                            &self.zeros[first..],
                            &mut states[..count.next_multiple_of(4)],
                        );
                        return;
                    }
                    _ => {}
                }
            }
        }
        for d in (0..self.bit_length - first).rev() {
            let level = self.level(d);
            for state in &mut states[..count] {
                let (bit, start1) = self.bit_vectors[level].access_rank1(state[0]);
                let start0 = state[0] - start1;
                let end1 = if OP == ACCESS {
                    0
                } else {
                    self.rank1(level, state[1])
                };
                let end0 = if OP == ACCESS { 0 } else { state[1] - end1 };
                let count0 = if OP == ACCESS { 0 } else { end0 - start0 };
                let bit = match OP {
                    ACCESS => bit,
                    RANK | RANK_LESSTHAN => (state[2] >> d) & 1 != 0,
                    _ => state[2] >= count0,
                };
                state[0] = if bit {
                    self.zeros[level] + start1
                } else {
                    start0
                };
                if OP != ACCESS {
                    state[1] = if bit { self.zeros[level] + end1 } else { end0 };
                }
                if OP == ACCESS || OP == QUANTILE {
                    state[3] |= (bit as usize) << d;
                }
                if OP == QUANTILE {
                    state[2] -= if bit { count0 } else { 0 };
                }
                if OP == RANK_LESSTHAN {
                    state[3] += if bit { count0 } else { 0 };
                }
            }
        }
    }

    /// get k-th value
    pub fn access(&self, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            let mut index = 0;
            for vector in &self.quad_vectors {
                let (digit, rank) = vector.access_rank(k);
                index = index * 4 + digit;
                k = vector.starts[digit] + rank;
            }
            return self.compress.values()[index].clone();
        }
        let mut idx = 0;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let (bit, rank1) = self.bit_vectors[level].access_rank1(k);
            idx |= (bit as usize) << d;
            k = if bit {
                self.zeros[level] + rank1
            } else {
                k - rank1
            };
        }
        self.compress.values()[idx].clone()
    }

    /// Returns the values at `indices` in input order.
    pub fn access_batch(&self, indices: impl IntoIterator<Item = usize>) -> Vec<T> {
        let indices: Vec<_> = indices.into_iter().collect();
        let mut result = Vec::with_capacity(indices.len());
        for indices in indices.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, &index) in states.iter_mut().zip(indices) {
                assert!(index < self.len);
                state[0] = index;
            }
            self.batch::<ACCESS>(&mut states, indices.len());
            result.extend(
                states[..indices.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// the number of val in range
    pub fn rank(&self, val: T, range: Range<usize>) -> usize {
        match self.compress.index_exact(&val) {
            Some(idx) => self.range_by_index(idx, range).len(),
            None => 0,
        }
    }

    /// Returns the number of exact matches for each `(value, range)` query.
    pub fn rank_batch(&self, queries: impl IntoIterator<Item = (T, Range<usize>)>) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                if let Some(index) = self.compress.index_exact(value) {
                    *state = [range.start, range.end, index, 0];
                }
            }
            self.batch::<RANK>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| state[1] - state[0]),
            );
        }
        result
    }

    /// index of k-th val
    pub fn select(&self, val: T, k: usize) -> Option<usize> {
        let idx = self.compress.index_exact(&val)?;
        let range = self.range_by_index(idx, 0..self.len);
        if range.len() <= k {
            return None;
        }
        let mut i = range.start + k;
        if !self.quad_vectors.is_empty() {
            for (level, vector) in self.quad_vectors.iter().enumerate().rev() {
                let digit = (idx >> ((self.quad_vectors.len() - level - 1) * 2)) & 3;
                i = vector.select(digit, i - vector.starts[digit]);
            }
            return Some(i);
        }
        for level in (self.compress.size().is_power_of_two() as usize..self.bit_length).rev() {
            if i >= self.zeros[level] {
                i = self.bit_vectors[level]
                    .select1(i - self.zeros[level])
                    .unwrap();
            } else {
                i = self.bit_vectors[level].select0(i).unwrap();
            }
        }
        Some(i)
    }

    /// get k-th smallest value in range
    pub fn quantile(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            return self.quad_quantile(range, k, 0, 0);
        }
        let mut idx = 0;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let z = end0 - start0;
            let bit = z <= k;
            k -= if bit { z } else { 0 };
            idx |= (bit as usize) << d;
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        self.compress.values()[idx].clone()
    }

    #[inline(always)]
    fn quad_quantile(
        &self,
        mut range: Range<usize>,
        mut k: usize,
        level: usize,
        mut index: usize,
    ) -> T {
        for vector in &self.quad_vectors[level..] {
            let start = vector.ranks(range.start);
            let end = vector.ranks(range.end);
            let mut digit = 0;
            while digit < 3 && k >= end[digit] - start[digit] {
                k -= end[digit] - start[digit];
                digit += 1;
            }
            index = index * 4 + digit;
            range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
        }
        self.compress.values()[index].clone()
    }

    pub fn quantile_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<usize>, usize)>,
    ) -> Vec<T> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (range, k)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len && *k < range.len());
                *state = [range.start, range.end, *k, 0];
            }
            if let [(range, k)] = queries {
                result.push(self.quantile(range.clone(), *k));
                continue;
            }
            self.batch::<QUANTILE>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// Returns the requested order statistics of one range. `ranks` must be sorted
    /// in nondecreasing order and every rank must be less than the range length.
    pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T> {
        assert!(range.start <= range.end && range.end <= self.len);
        assert!(ranks.windows(2).all(|pair| pair[0] <= pair[1]));
        assert!(ranks.last().is_none_or(|&k| k < range.len()));
        if ranks.is_empty() {
            return Vec::new();
        }
        if ranks[0] == ranks[ranks.len() - 1] {
            return vec![self.quantile(range, ranks[0]); ranks.len()];
        }
        let use_batch = ranks.len() < 1024;
        #[cfg(target_arch = "x86_64")]
        let use_batch = use_batch || self.backend == super::SimdBackend::Avx512;
        // Sparse requests amortize less of the shared traversal's branching work.
        if use_batch && ranks.len() < self.compress.size().div_ceil(16) {
            return self.quantile_batch(ranks.iter().map(|&k| (range.clone(), k)));
        }
        let mut result = Vec::with_capacity(ranks.len());
        if !self.quad_vectors.is_empty() {
            let mut stack = vec![(range, 0, 0, 0, ranks)];
            while let Some((range, base, index, level, ranks)) = stack.pop() {
                if level == self.quad_vectors.len() {
                    result.resize(
                        result.len() + ranks.len(),
                        self.compress.values()[index].clone(),
                    );
                    continue;
                }
                if ranks.len() == 1 {
                    result.push(self.quad_quantile(range, ranks[0] - base, level, index));
                    continue;
                }
                let vector = &self.quad_vectors[level];
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let mut split = base + range.len();
                let mut requested = ranks.len();
                for digit in (0..4).rev() {
                    split -= end[digit] - start[digit];
                    let boundary = ranks[..requested].partition_point(|&k| k < split);
                    if boundary < requested {
                        stack.push((
                            vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit],
                            split,
                            index * 4 + digit,
                            level + 1,
                            &ranks[boundary..requested],
                        ));
                    }
                    requested = boundary;
                }
            }
            return result;
        }
        let mut stack = vec![(range, 0, 0, self.bit_length, ranks)];
        while let Some((range, base, index, bits, ranks)) = stack.pop() {
            if bits == 0 {
                result.resize(
                    result.len() + ranks.len(),
                    self.compress.values()[index].clone(),
                );
                continue;
            }
            let d = bits - 1;
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let split = base + end0 - start0;
            let boundary = ranks.partition_point(|&k| k < split);
            if boundary < ranks.len() {
                stack.push((
                    self.zeros[level] + start1..self.zeros[level] + end1,
                    split,
                    index | (1 << d),
                    d,
                    &ranks[boundary..],
                ));
            }
            if boundary != 0 {
                stack.push((start0..end0, base, index, d, &ranks[..boundary]));
            }
        }
        result
    }

    /// get k-th smallest value out of range
    pub fn quantile_outer(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            let mut outer = 0..self.len;
            let mut index = 0;
            for vector in &self.quad_vectors {
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let outer_start = vector.ranks(outer.start);
                let outer_end = vector.ranks(outer.end);
                let mut digit = 0;
                while digit < 3 {
                    let count = outer_end[digit] - outer_start[digit] - (end[digit] - start[digit]);
                    if k < count {
                        break;
                    }
                    k -= count;
                    digit += 1;
                }
                index = index * 4 + digit;
                range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
                outer = vector.starts[digit] + outer_start[digit]
                    ..vector.starts[digit] + outer_end[digit];
            }
            return self.compress.values()[index].clone();
        }
        let mut idx = 0;
        let mut orange = 0..self.len;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let range_start1 = self.rank1(level, range.start);
            let range_end1 = self.rank1(level, range.end);
            let outer_start1 = self.rank1(level, orange.start);
            let outer_end1 = self.rank1(level, orange.end);
            let range_start0 = range.start - range_start1;
            let range_end0 = range.end - range_end1;
            let outer_start0 = orange.start - outer_start1;
            let outer_end0 = orange.end - outer_end1;
            let z = (outer_end0 - outer_start0) - (range_end0 - range_start0);
            if z <= k {
                k -= z;
                idx |= 1 << d;
                range.start = self.zeros[level] + range_start1;
                range.end = self.zeros[level] + range_end1;
                orange.start = self.zeros[level] + outer_start1;
                orange.end = self.zeros[level] + outer_end1;
            } else {
                range.start = range_start0;
                range.end = range_end0;
                orange.start = outer_start0;
                orange.end = outer_end0;
            }
        }
        self.compress.values()[idx].clone()
    }

    /// the number of value less than val in range
    pub fn rank_lessthan(&self, val: T, range: Range<usize>) -> usize {
        let index = self.compress.index_lower_bound(&val);
        if index == self.compress.size() {
            return range.len();
        }
        if self.len >= 262144 && !self.quad_vectors.is_empty() {
            self.quad_rank_lessthan_index(index, range, 0)
        } else {
            let bits = usize::BITS as usize
                - self.compress.size().saturating_sub(1).leading_zeros() as usize;
            self.rank_lessthan_index(index, range, bits)
        }
    }

    /// Returns the number of values below each query's threshold.
    pub fn rank_lessthan_batch(
        &self,
        queries: impl IntoIterator<Item = (T, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                *state = [
                    range.start,
                    range.end,
                    self.compress.index_lower_bound(value),
                    0,
                ];
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len());
            result.extend(states[..queries.len()].iter().map(|state| state[3]));
        }
        result
    }

    fn quad_rank_lessthan_index(
        &self,
        index: usize,
        mut range: Range<usize>,
        start_level: usize,
    ) -> usize {
        let mut result = 0;
        for (level, vector) in self.quad_vectors.iter().enumerate().skip(start_level) {
            let d = (self.quad_vectors.len() - level - 1) * 2;
            if d + 2 <= index.trailing_zeros() as usize {
                break;
            }
            let digit = (index >> d) & 3;
            let ranks = |position: usize| {
                let block = &vector.blocks[position / 64];
                let mask = !(u64::MAX << (position % 64));
                let low = if digit & 1 == 0 { !block.lo } else { block.lo };
                let high = if digit & 2 == 0 { !block.hi } else { block.hi };
                let exact = block.rank[digit] as usize + (low & high & mask).count_ones() as usize;
                let less = match digit {
                    0 => 0,
                    1 => {
                        block.rank[0] as usize
                            + (!(block.lo | block.hi) & mask).count_ones() as usize
                    }
                    2 => {
                        block.rank[0] as usize
                            + block.rank[1] as usize
                            + (!block.hi & mask).count_ones() as usize
                    }
                    _ => position - exact,
                };
                (exact, less)
            };
            let (start, start_less) = ranks(range.start);
            let (end, end_less) = ranks(range.end);
            result += end_less - start_less;
            range = vector.starts[digit] + start..vector.starts[digit] + end;
        }
        result
    }

    fn rank_lessthan_index(&self, idx: usize, mut range: Range<usize>, bits: usize) -> usize {
        let mut res = 0;
        for d in (idx.trailing_zeros() as usize..bits).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let bit = (idx >> d) & 1 != 0;
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            res += if bit { end0 - start0 } else { 0 };
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        res
    }

    /// the number of valrange in range
    pub fn rank_range(&self, valrange: Range<T>, mut range: Range<usize>) -> usize {
        let lower = self.compress.index_lower_bound(&valrange.start);
        let upper = self.compress.index_lower_bound(&valrange.end);
        if lower >= upper {
            return 0;
        }
        if !self.quad_vectors.is_empty() {
            if upper == self.compress.size() {
                return range.len() - self.quad_rank_lessthan_index(lower, range, 0);
            }
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let low = (lower >> d) & 3;
                let high = (upper >> d) & 3;
                if low == high {
                    range = vector.starts[low] + vector.rank(low, range.start)
                        ..vector.starts[low] + vector.rank(low, range.end);
                    continue;
                }
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let count: usize = (low..high).map(|digit| end[digit] - start[digit]).sum();
                return count
                    - self.quad_rank_lessthan_index(
                        lower,
                        vector.starts[low] + start[low]..vector.starts[low] + end[low],
                        level + 1,
                    )
                    + self.quad_rank_lessthan_index(
                        upper,
                        vector.starts[high] + start[high]..vector.starts[high] + end[high],
                        level + 1,
                    );
            }
            return 0;
        }
        for d in (0..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let zero_range = range.start - start1..range.end - end1;
            let one_range = self.zeros[level] + start1..self.zeros[level] + end1;
            if ((lower ^ upper) >> d) & 1 != 0 {
                return zero_range.len() - self.rank_lessthan_index(lower, zero_range, d)
                    + self.rank_lessthan_index(upper, one_range, d);
            }
            range = if (lower >> d) & 1 != 0 {
                one_range
            } else {
                zero_range
            };
        }
        0
    }

    /// Returns each query's count in a half-open value range.
    pub fn rank_range_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<T>, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(8) {
            let mut states = [[0; 4]; 16];
            for (i, (values, range)) in queries.iter().enumerate() {
                assert!(range.start <= range.end && range.end <= self.len);
                let lower = self.compress.index_lower_bound(&values.start);
                let upper = self.compress.index_lower_bound(&values.end);
                if lower < upper {
                    states[i * 2] = [range.start, range.end, lower, 0];
                    states[i * 2 + 1] = [range.start, range.end, upper, 0];
                }
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len() * 2);
            result.extend(
                states[..queries.len() * 2]
                    .as_chunks::<2>()
                    .0
                    .iter()
                    .map(|pair| pair[1][3] - pair[0][3]),
            );
        }
        result
    }

    pub fn query_less_than<F>(&self, val: T, mut range: Range<usize>, mut f: F)
    where
        F: FnMut(usize, Range<usize>),
    {
        let idx = self.compress.index_lower_bound(&val);
        if !self.quad_vectors.is_empty() {
            if idx == self.compress.size() && idx.is_power_of_two() {
                f(self.bit_length - 1, range);
                return;
            }
            for (level, vector) in self.quad_vectors.iter().enumerate().take(
                self.quad_vectors
                    .len()
                    .saturating_sub(idx.trailing_zeros() as usize / 2),
            ) {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let digit = (idx >> d) & 3;
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                if digit & 2 != 0 {
                    f(d + 1, start[0] + start[1]..end[0] + end[1]);
                }
                if digit & 1 != 0 {
                    let zero = digit & 2;
                    f(
                        d,
                        vector.starts[zero] + start[zero]..vector.starts[zero] + end[zero],
                    );
                }
                range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
            }
            return;
        }
        for d in (idx.trailing_zeros() as usize..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            if ((idx >> d) & 1) != 0 {
                f(d, start0..end0);
                range.start = self.zeros[level] + start1;
                range.end = self.zeros[level] + end1;
            } else {
                range.start = start0;
                range.end = end0;
            }
        }
    }

fn reorder<U>(&self, level: usize, current: Vec<U>) -> Vec<U>

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 1440)

    pub fn build_fold<M>(&self, weights: &[M::T]) -> WaveletMatrixFold<'_, T, M>
    where
        M: AbelianGroup,
    {
        assert_eq!(weights.len(), self.len);
        let mut offsets = Vec::with_capacity(self.bit_length);
        let mut prefix = Vec::with_capacity(self.zeros.iter().map(|&zero| zero + 1).sum());
        let mut current: Vec<M::T> = weights.to_vec();
        for level in 0..self.bit_length {
            current = self.reorder(level, current);
            offsets.push(prefix.len());
            let mut acc = M::unit();
            prefix.push(acc.clone());
            for w in &current[..self.zeros[level]] {
                acc = M::operate(&acc, w);
                prefix.push(acc.clone());
            }
        }
        WaveletMatrixFold {
            wavelet_matrix: self,
            prefix,
            offsets,
        }
    }

    pub fn build_point_add<M>(&self, weights: &[M::T]) -> WaveletMatrixPointAdd<'_, T, M>
    where
        M: AbelianGroup,
    {
        assert_eq!(weights.len(), self.len);
        let mut current = weights.to_vec();
        let mut bits = Vec::with_capacity(self.bit_length);
        for level in 0..self.bit_length {
            current = self.reorder(level, current);
            bits.push(BinaryIndexedTree::from_slice(&current[..self.zeros[level]]));
        }
        WaveletMatrixPointAdd {
            wavelet_matrix: self,
            bits,
        }
    }

fn range_by_index(&self, idx: usize, range: Range<usize>) -> Range<usize>

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 921)

    pub fn rank(&self, val: T, range: Range<usize>) -> usize {
        match self.compress.index_exact(&val) {
            Some(idx) => self.range_by_index(idx, range).len(),
            None => 0,
        }
    }

    /// Returns the number of exact matches for each `(value, range)` query.
    pub fn rank_batch(&self, queries: impl IntoIterator<Item = (T, Range<usize>)>) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                if let Some(index) = self.compress.index_exact(value) {
                    *state = [range.start, range.end, index, 0];
                }
            }
            self.batch::<RANK>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| state[1] - state[0]),
            );
        }
        result
    }

    /// index of k-th val
    pub fn select(&self, val: T, k: usize) -> Option<usize> {
        let idx = self.compress.index_exact(&val)?;
        let range = self.range_by_index(idx, 0..self.len);
        if range.len() <= k {
            return None;
        }
        let mut i = range.start + k;
        if !self.quad_vectors.is_empty() {
            for (level, vector) in self.quad_vectors.iter().enumerate().rev() {
                let digit = (idx >> ((self.quad_vectors.len() - level - 1) * 2)) & 3;
                i = vector.select(digit, i - vector.starts[digit]);
            }
            return Some(i);
        }
        for level in (self.compress.size().is_power_of_two() as usize..self.bit_length).rev() {
            if i >= self.zeros[level] {
                i = self.bit_vectors[level]
                    .select1(i - self.zeros[level])
                    .unwrap();
            } else {
                i = self.bit_vectors[level].select0(i).unwrap();
            }
        }
        Some(i)
    }

fn batch<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize)

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 908)

    pub fn access_batch(&self, indices: impl IntoIterator<Item = usize>) -> Vec<T> {
        let indices: Vec<_> = indices.into_iter().collect();
        let mut result = Vec::with_capacity(indices.len());
        for indices in indices.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, &index) in states.iter_mut().zip(indices) {
                assert!(index < self.len);
                state[0] = index;
            }
            self.batch::<ACCESS>(&mut states, indices.len());
            result.extend(
                states[..indices.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// the number of val in range
    pub fn rank(&self, val: T, range: Range<usize>) -> usize {
        match self.compress.index_exact(&val) {
            Some(idx) => self.range_by_index(idx, range).len(),
            None => 0,
        }
    }

    /// Returns the number of exact matches for each `(value, range)` query.
    pub fn rank_batch(&self, queries: impl IntoIterator<Item = (T, Range<usize>)>) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                if let Some(index) = self.compress.index_exact(value) {
                    *state = [range.start, range.end, index, 0];
                }
            }
            self.batch::<RANK>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| state[1] - state[0]),
            );
        }
        result
    }

    /// index of k-th val
    pub fn select(&self, val: T, k: usize) -> Option<usize> {
        let idx = self.compress.index_exact(&val)?;
        let range = self.range_by_index(idx, 0..self.len);
        if range.len() <= k {
            return None;
        }
        let mut i = range.start + k;
        if !self.quad_vectors.is_empty() {
            for (level, vector) in self.quad_vectors.iter().enumerate().rev() {
                let digit = (idx >> ((self.quad_vectors.len() - level - 1) * 2)) & 3;
                i = vector.select(digit, i - vector.starts[digit]);
            }
            return Some(i);
        }
        for level in (self.compress.size().is_power_of_two() as usize..self.bit_length).rev() {
            if i >= self.zeros[level] {
                i = self.bit_vectors[level]
                    .select1(i - self.zeros[level])
                    .unwrap();
            } else {
                i = self.bit_vectors[level].select0(i).unwrap();
            }
        }
        Some(i)
    }

    /// get k-th smallest value in range
    pub fn quantile(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            return self.quad_quantile(range, k, 0, 0);
        }
        let mut idx = 0;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let z = end0 - start0;
            let bit = z <= k;
            k -= if bit { z } else { 0 };
            idx |= (bit as usize) << d;
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        self.compress.values()[idx].clone()
    }

    #[inline(always)]
    fn quad_quantile(
        &self,
        mut range: Range<usize>,
        mut k: usize,
        level: usize,
        mut index: usize,
    ) -> T {
        for vector in &self.quad_vectors[level..] {
            let start = vector.ranks(range.start);
            let end = vector.ranks(range.end);
            let mut digit = 0;
            while digit < 3 && k >= end[digit] - start[digit] {
                k -= end[digit] - start[digit];
                digit += 1;
            }
            index = index * 4 + digit;
            range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
        }
        self.compress.values()[index].clone()
    }

    pub fn quantile_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<usize>, usize)>,
    ) -> Vec<T> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (range, k)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len && *k < range.len());
                *state = [range.start, range.end, *k, 0];
            }
            if let [(range, k)] = queries {
                result.push(self.quantile(range.clone(), *k));
                continue;
            }
            self.batch::<QUANTILE>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// Returns the requested order statistics of one range. `ranks` must be sorted
    /// in nondecreasing order and every rank must be less than the range length.
    pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T> {
        assert!(range.start <= range.end && range.end <= self.len);
        assert!(ranks.windows(2).all(|pair| pair[0] <= pair[1]));
        assert!(ranks.last().is_none_or(|&k| k < range.len()));
        if ranks.is_empty() {
            return Vec::new();
        }
        if ranks[0] == ranks[ranks.len() - 1] {
            return vec![self.quantile(range, ranks[0]); ranks.len()];
        }
        let use_batch = ranks.len() < 1024;
        #[cfg(target_arch = "x86_64")]
        let use_batch = use_batch || self.backend == super::SimdBackend::Avx512;
        // Sparse requests amortize less of the shared traversal's branching work.
        if use_batch && ranks.len() < self.compress.size().div_ceil(16) {
            return self.quantile_batch(ranks.iter().map(|&k| (range.clone(), k)));
        }
        let mut result = Vec::with_capacity(ranks.len());
        if !self.quad_vectors.is_empty() {
            let mut stack = vec![(range, 0, 0, 0, ranks)];
            while let Some((range, base, index, level, ranks)) = stack.pop() {
                if level == self.quad_vectors.len() {
                    result.resize(
                        result.len() + ranks.len(),
                        self.compress.values()[index].clone(),
                    );
                    continue;
                }
                if ranks.len() == 1 {
                    result.push(self.quad_quantile(range, ranks[0] - base, level, index));
                    continue;
                }
                let vector = &self.quad_vectors[level];
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let mut split = base + range.len();
                let mut requested = ranks.len();
                for digit in (0..4).rev() {
                    split -= end[digit] - start[digit];
                    let boundary = ranks[..requested].partition_point(|&k| k < split);
                    if boundary < requested {
                        stack.push((
                            vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit],
                            split,
                            index * 4 + digit,
                            level + 1,
                            &ranks[boundary..requested],
                        ));
                    }
                    requested = boundary;
                }
            }
            return result;
        }
        let mut stack = vec![(range, 0, 0, self.bit_length, ranks)];
        while let Some((range, base, index, bits, ranks)) = stack.pop() {
            if bits == 0 {
                result.resize(
                    result.len() + ranks.len(),
                    self.compress.values()[index].clone(),
                );
                continue;
            }
            let d = bits - 1;
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let split = base + end0 - start0;
            let boundary = ranks.partition_point(|&k| k < split);
            if boundary < ranks.len() {
                stack.push((
                    self.zeros[level] + start1..self.zeros[level] + end1,
                    split,
                    index | (1 << d),
                    d,
                    &ranks[boundary..],
                ));
            }
            if boundary != 0 {
                stack.push((start0..end0, base, index, d, &ranks[..boundary]));
            }
        }
        result
    }

    /// get k-th smallest value out of range
    pub fn quantile_outer(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            let mut outer = 0..self.len;
            let mut index = 0;
            for vector in &self.quad_vectors {
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let outer_start = vector.ranks(outer.start);
                let outer_end = vector.ranks(outer.end);
                let mut digit = 0;
                while digit < 3 {
                    let count = outer_end[digit] - outer_start[digit] - (end[digit] - start[digit]);
                    if k < count {
                        break;
                    }
                    k -= count;
                    digit += 1;
                }
                index = index * 4 + digit;
                range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
                outer = vector.starts[digit] + outer_start[digit]
                    ..vector.starts[digit] + outer_end[digit];
            }
            return self.compress.values()[index].clone();
        }
        let mut idx = 0;
        let mut orange = 0..self.len;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let range_start1 = self.rank1(level, range.start);
            let range_end1 = self.rank1(level, range.end);
            let outer_start1 = self.rank1(level, orange.start);
            let outer_end1 = self.rank1(level, orange.end);
            let range_start0 = range.start - range_start1;
            let range_end0 = range.end - range_end1;
            let outer_start0 = orange.start - outer_start1;
            let outer_end0 = orange.end - outer_end1;
            let z = (outer_end0 - outer_start0) - (range_end0 - range_start0);
            if z <= k {
                k -= z;
                idx |= 1 << d;
                range.start = self.zeros[level] + range_start1;
                range.end = self.zeros[level] + range_end1;
                orange.start = self.zeros[level] + outer_start1;
                orange.end = self.zeros[level] + outer_end1;
            } else {
                range.start = range_start0;
                range.end = range_end0;
                orange.start = outer_start0;
                orange.end = outer_end0;
            }
        }
        self.compress.values()[idx].clone()
    }

    /// the number of value less than val in range
    pub fn rank_lessthan(&self, val: T, range: Range<usize>) -> usize {
        let index = self.compress.index_lower_bound(&val);
        if index == self.compress.size() {
            return range.len();
        }
        if self.len >= 262144 && !self.quad_vectors.is_empty() {
            self.quad_rank_lessthan_index(index, range, 0)
        } else {
            let bits = usize::BITS as usize
                - self.compress.size().saturating_sub(1).leading_zeros() as usize;
            self.rank_lessthan_index(index, range, bits)
        }
    }

    /// Returns the number of values below each query's threshold.
    pub fn rank_lessthan_batch(
        &self,
        queries: impl IntoIterator<Item = (T, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                *state = [
                    range.start,
                    range.end,
                    self.compress.index_lower_bound(value),
                    0,
                ];
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len());
            result.extend(states[..queries.len()].iter().map(|state| state[3]));
        }
        result
    }

    fn quad_rank_lessthan_index(
        &self,
        index: usize,
        mut range: Range<usize>,
        start_level: usize,
    ) -> usize {
        let mut result = 0;
        for (level, vector) in self.quad_vectors.iter().enumerate().skip(start_level) {
            let d = (self.quad_vectors.len() - level - 1) * 2;
            if d + 2 <= index.trailing_zeros() as usize {
                break;
            }
            let digit = (index >> d) & 3;
            let ranks = |position: usize| {
                let block = &vector.blocks[position / 64];
                let mask = !(u64::MAX << (position % 64));
                let low = if digit & 1 == 0 { !block.lo } else { block.lo };
                let high = if digit & 2 == 0 { !block.hi } else { block.hi };
                let exact = block.rank[digit] as usize + (low & high & mask).count_ones() as usize;
                let less = match digit {
                    0 => 0,
                    1 => {
                        block.rank[0] as usize
                            + (!(block.lo | block.hi) & mask).count_ones() as usize
                    }
                    2 => {
                        block.rank[0] as usize
                            + block.rank[1] as usize
                            + (!block.hi & mask).count_ones() as usize
                    }
                    _ => position - exact,
                };
                (exact, less)
            };
            let (start, start_less) = ranks(range.start);
            let (end, end_less) = ranks(range.end);
            result += end_less - start_less;
            range = vector.starts[digit] + start..vector.starts[digit] + end;
        }
        result
    }

    fn rank_lessthan_index(&self, idx: usize, mut range: Range<usize>, bits: usize) -> usize {
        let mut res = 0;
        for d in (idx.trailing_zeros() as usize..bits).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let bit = (idx >> d) & 1 != 0;
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            res += if bit { end0 - start0 } else { 0 };
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        res
    }

    /// the number of valrange in range
    pub fn rank_range(&self, valrange: Range<T>, mut range: Range<usize>) -> usize {
        let lower = self.compress.index_lower_bound(&valrange.start);
        let upper = self.compress.index_lower_bound(&valrange.end);
        if lower >= upper {
            return 0;
        }
        if !self.quad_vectors.is_empty() {
            if upper == self.compress.size() {
                return range.len() - self.quad_rank_lessthan_index(lower, range, 0);
            }
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let low = (lower >> d) & 3;
                let high = (upper >> d) & 3;
                if low == high {
                    range = vector.starts[low] + vector.rank(low, range.start)
                        ..vector.starts[low] + vector.rank(low, range.end);
                    continue;
                }
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let count: usize = (low..high).map(|digit| end[digit] - start[digit]).sum();
                return count
                    - self.quad_rank_lessthan_index(
                        lower,
                        vector.starts[low] + start[low]..vector.starts[low] + end[low],
                        level + 1,
                    )
                    + self.quad_rank_lessthan_index(
                        upper,
                        vector.starts[high] + start[high]..vector.starts[high] + end[high],
                        level + 1,
                    );
            }
            return 0;
        }
        for d in (0..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let zero_range = range.start - start1..range.end - end1;
            let one_range = self.zeros[level] + start1..self.zeros[level] + end1;
            if ((lower ^ upper) >> d) & 1 != 0 {
                return zero_range.len() - self.rank_lessthan_index(lower, zero_range, d)
                    + self.rank_lessthan_index(upper, one_range, d);
            }
            range = if (lower >> d) & 1 != 0 {
                one_range
            } else {
                zero_range
            };
        }
        0
    }

    /// Returns each query's count in a half-open value range.
    pub fn rank_range_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<T>, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(8) {
            let mut states = [[0; 4]; 16];
            for (i, (values, range)) in queries.iter().enumerate() {
                assert!(range.start <= range.end && range.end <= self.len);
                let lower = self.compress.index_lower_bound(&values.start);
                let upper = self.compress.index_lower_bound(&values.end);
                if lower < upper {
                    states[i * 2] = [range.start, range.end, lower, 0];
                    states[i * 2 + 1] = [range.start, range.end, upper, 0];
                }
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len() * 2);
            result.extend(
                states[..queries.len() * 2]
                    .as_chunks::<2>()
                    .0
                    .iter()
                    .map(|pair| pair[1][3] - pair[0][3]),
            );
        }
        result
    }

unsafe fn batch_popcnt<const OP: u8>( &self, states: &mut [[usize; 4]; 16], count: usize, )

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 718)

    fn batch<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        #[cfg(target_arch = "x86_64")]
        if !cfg!(target_feature = "popcnt") && is_x86_feature_detected!("popcnt") {
            // SAFETY: POPCNT is checked above.
            unsafe {
                self.batch_popcnt::<OP>(states, count);
            }
            return;
        }
        self.batch_inner::<OP>(states, count);
    }

fn batch_inner<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize)

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 722)

    fn batch<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        #[cfg(target_arch = "x86_64")]
        if !cfg!(target_feature = "popcnt") && is_x86_feature_detected!("popcnt") {
            // SAFETY: POPCNT is checked above.
            unsafe {
                self.batch_popcnt::<OP>(states, count);
            }
            return;
        }
        self.batch_inner::<OP>(states, count);
    }

    #[cfg(target_arch = "x86_64")]
    #[target_feature(enable = "popcnt")]
    unsafe fn batch_popcnt<const OP: u8>(&self, states: &mut [[usize; 4]; 16], count: usize) {
        self.batch_inner::<OP>(states, count);
    }

pub fn access(&self, k: usize) -> T

get k-th value

pub fn access_batch(&self, indices: impl IntoIterator<Item = usize>) -> Vec<T>

Returns the values at indices in input order.

pub fn rank(&self, val: T, range: Range<usize>) -> usize

the number of val in range

Examples found in repository

crates/library_checker/src/data_structure/static_range_frequency.rs (line 23)

pub fn static_range_frequency_wavelet_matrix(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, a: [usize; n]);
    let wm = WaveletMatrix::new(a);
    for _ in 0..q {
        sc!(l, r, x: usize);
        let ans = wm.rank(x, l..r);
        pp!(ans);
    }
}

pub fn rank_batch( &self, queries: impl IntoIterator<Item = (T, Range<usize>)>, ) -> Vec<usize>

Returns the number of exact matches for each (value, range) query.

pub fn select(&self, val: T, k: usize) -> Option<usize>

index of k-th val

pub fn quantile(&self, range: Range<usize>, k: usize) -> T

get k-th smallest value in range

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 1036)

    pub fn quantile_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<usize>, usize)>,
    ) -> Vec<T> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (range, k)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len && *k < range.len());
                *state = [range.start, range.end, *k, 0];
            }
            if let [(range, k)] = queries {
                result.push(self.quantile(range.clone(), *k));
                continue;
            }
            self.batch::<QUANTILE>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// Returns the requested order statistics of one range. `ranks` must be sorted
    /// in nondecreasing order and every rank must be less than the range length.
    pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T> {
        assert!(range.start <= range.end && range.end <= self.len);
        assert!(ranks.windows(2).all(|pair| pair[0] <= pair[1]));
        assert!(ranks.last().is_none_or(|&k| k < range.len()));
        if ranks.is_empty() {
            return Vec::new();
        }
        if ranks[0] == ranks[ranks.len() - 1] {
            return vec![self.quantile(range, ranks[0]); ranks.len()];
        }
        let use_batch = ranks.len() < 1024;
        #[cfg(target_arch = "x86_64")]
        let use_batch = use_batch || self.backend == super::SimdBackend::Avx512;
        // Sparse requests amortize less of the shared traversal's branching work.
        if use_batch && ranks.len() < self.compress.size().div_ceil(16) {
            return self.quantile_batch(ranks.iter().map(|&k| (range.clone(), k)));
        }
        let mut result = Vec::with_capacity(ranks.len());
        if !self.quad_vectors.is_empty() {
            let mut stack = vec![(range, 0, 0, 0, ranks)];
            while let Some((range, base, index, level, ranks)) = stack.pop() {
                if level == self.quad_vectors.len() {
                    result.resize(
                        result.len() + ranks.len(),
                        self.compress.values()[index].clone(),
                    );
                    continue;
                }
                if ranks.len() == 1 {
                    result.push(self.quad_quantile(range, ranks[0] - base, level, index));
                    continue;
                }
                let vector = &self.quad_vectors[level];
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let mut split = base + range.len();
                let mut requested = ranks.len();
                for digit in (0..4).rev() {
                    split -= end[digit] - start[digit];
                    let boundary = ranks[..requested].partition_point(|&k| k < split);
                    if boundary < requested {
                        stack.push((
                            vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit],
                            split,
                            index * 4 + digit,
                            level + 1,
                            &ranks[boundary..requested],
                        ));
                    }
                    requested = boundary;
                }
            }
            return result;
        }
        let mut stack = vec![(range, 0, 0, self.bit_length, ranks)];
        while let Some((range, base, index, bits, ranks)) = stack.pop() {
            if bits == 0 {
                result.resize(
                    result.len() + ranks.len(),
                    self.compress.values()[index].clone(),
                );
                continue;
            }
            let d = bits - 1;
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let split = base + end0 - start0;
            let boundary = ranks.partition_point(|&k| k < split);
            if boundary < ranks.len() {
                stack.push((
                    self.zeros[level] + start1..self.zeros[level] + end1,
                    split,
                    index | (1 << d),
                    d,
                    &ranks[boundary..],
                ));
            }
            if boundary != 0 {
                stack.push((start0..end0, base, index, d, &ranks[..boundary]));
            }
        }
        result
    }

fn quad_quantile( &self, range: Range<usize>, k: usize, level: usize, index: usize, ) -> T

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 978)

    pub fn quantile(&self, mut range: Range<usize>, mut k: usize) -> T {
        if !self.quad_vectors.is_empty() {
            return self.quad_quantile(range, k, 0, 0);
        }
        let mut idx = 0;
        for d in (0..self.bit_length - self.compress.size().is_power_of_two() as usize).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let z = end0 - start0;
            let bit = z <= k;
            k -= if bit { z } else { 0 };
            idx |= (bit as usize) << d;
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        self.compress.values()[idx].clone()
    }

    #[inline(always)]
    fn quad_quantile(
        &self,
        mut range: Range<usize>,
        mut k: usize,
        level: usize,
        mut index: usize,
    ) -> T {
        for vector in &self.quad_vectors[level..] {
            let start = vector.ranks(range.start);
            let end = vector.ranks(range.end);
            let mut digit = 0;
            while digit < 3 && k >= end[digit] - start[digit] {
                k -= end[digit] - start[digit];
                digit += 1;
            }
            index = index * 4 + digit;
            range = vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit];
        }
        self.compress.values()[index].clone()
    }

    pub fn quantile_batch(
        &self,
        queries: impl IntoIterator<Item = (Range<usize>, usize)>,
    ) -> Vec<T> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (range, k)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len && *k < range.len());
                *state = [range.start, range.end, *k, 0];
            }
            if let [(range, k)] = queries {
                result.push(self.quantile(range.clone(), *k));
                continue;
            }
            self.batch::<QUANTILE>(&mut states, queries.len());
            result.extend(
                states[..queries.len()]
                    .iter()
                    .map(|state| self.compress.values()[state[3]].clone()),
            );
        }
        result
    }

    /// Returns the requested order statistics of one range. `ranks` must be sorted
    /// in nondecreasing order and every rank must be less than the range length.
    pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T> {
        assert!(range.start <= range.end && range.end <= self.len);
        assert!(ranks.windows(2).all(|pair| pair[0] <= pair[1]));
        assert!(ranks.last().is_none_or(|&k| k < range.len()));
        if ranks.is_empty() {
            return Vec::new();
        }
        if ranks[0] == ranks[ranks.len() - 1] {
            return vec![self.quantile(range, ranks[0]); ranks.len()];
        }
        let use_batch = ranks.len() < 1024;
        #[cfg(target_arch = "x86_64")]
        let use_batch = use_batch || self.backend == super::SimdBackend::Avx512;
        // Sparse requests amortize less of the shared traversal's branching work.
        if use_batch && ranks.len() < self.compress.size().div_ceil(16) {
            return self.quantile_batch(ranks.iter().map(|&k| (range.clone(), k)));
        }
        let mut result = Vec::with_capacity(ranks.len());
        if !self.quad_vectors.is_empty() {
            let mut stack = vec![(range, 0, 0, 0, ranks)];
            while let Some((range, base, index, level, ranks)) = stack.pop() {
                if level == self.quad_vectors.len() {
                    result.resize(
                        result.len() + ranks.len(),
                        self.compress.values()[index].clone(),
                    );
                    continue;
                }
                if ranks.len() == 1 {
                    result.push(self.quad_quantile(range, ranks[0] - base, level, index));
                    continue;
                }
                let vector = &self.quad_vectors[level];
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let mut split = base + range.len();
                let mut requested = ranks.len();
                for digit in (0..4).rev() {
                    split -= end[digit] - start[digit];
                    let boundary = ranks[..requested].partition_point(|&k| k < split);
                    if boundary < requested {
                        stack.push((
                            vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit],
                            split,
                            index * 4 + digit,
                            level + 1,
                            &ranks[boundary..requested],
                        ));
                    }
                    requested = boundary;
                }
            }
            return result;
        }
        let mut stack = vec![(range, 0, 0, self.bit_length, ranks)];
        while let Some((range, base, index, bits, ranks)) = stack.pop() {
            if bits == 0 {
                result.resize(
                    result.len() + ranks.len(),
                    self.compress.values()[index].clone(),
                );
                continue;
            }
            let d = bits - 1;
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let split = base + end0 - start0;
            let boundary = ranks.partition_point(|&k| k < split);
            if boundary < ranks.len() {
                stack.push((
                    self.zeros[level] + start1..self.zeros[level] + end1,
                    split,
                    index | (1 << d),
                    d,
                    &ranks[boundary..],
                ));
            }
            if boundary != 0 {
                stack.push((start0..end0, base, index, d, &ranks[..boundary]));
            }
        }
        result
    }

pub fn quantile_batch( &self, queries: impl IntoIterator<Item = (Range<usize>, usize)>, ) -> Vec<T>

Examples found in repository

crates/library_checker/src/data_structure/range_kth_smallest.rs (line 9)

pub fn range_kth_smallest(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, a: [usize; n], queries: [(usize, usize, usize); iter q]);
    let wm = WaveletMatrix::new(a);
    let results = wm.quantile_batch(queries.map(|(l, r, k)| (l..r, k)));
    pp!(@lf @it results);
}

crates/competitive/src/data_structure/wavelet_matrix.rs (line 1066)

    pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T> {
        assert!(range.start <= range.end && range.end <= self.len);
        assert!(ranks.windows(2).all(|pair| pair[0] <= pair[1]));
        assert!(ranks.last().is_none_or(|&k| k < range.len()));
        if ranks.is_empty() {
            return Vec::new();
        }
        if ranks[0] == ranks[ranks.len() - 1] {
            return vec![self.quantile(range, ranks[0]); ranks.len()];
        }
        let use_batch = ranks.len() < 1024;
        #[cfg(target_arch = "x86_64")]
        let use_batch = use_batch || self.backend == super::SimdBackend::Avx512;
        // Sparse requests amortize less of the shared traversal's branching work.
        if use_batch && ranks.len() < self.compress.size().div_ceil(16) {
            return self.quantile_batch(ranks.iter().map(|&k| (range.clone(), k)));
        }
        let mut result = Vec::with_capacity(ranks.len());
        if !self.quad_vectors.is_empty() {
            let mut stack = vec![(range, 0, 0, 0, ranks)];
            while let Some((range, base, index, level, ranks)) = stack.pop() {
                if level == self.quad_vectors.len() {
                    result.resize(
                        result.len() + ranks.len(),
                        self.compress.values()[index].clone(),
                    );
                    continue;
                }
                if ranks.len() == 1 {
                    result.push(self.quad_quantile(range, ranks[0] - base, level, index));
                    continue;
                }
                let vector = &self.quad_vectors[level];
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let mut split = base + range.len();
                let mut requested = ranks.len();
                for digit in (0..4).rev() {
                    split -= end[digit] - start[digit];
                    let boundary = ranks[..requested].partition_point(|&k| k < split);
                    if boundary < requested {
                        stack.push((
                            vector.starts[digit] + start[digit]..vector.starts[digit] + end[digit],
                            split,
                            index * 4 + digit,
                            level + 1,
                            &ranks[boundary..requested],
                        ));
                    }
                    requested = boundary;
                }
            }
            return result;
        }
        let mut stack = vec![(range, 0, 0, self.bit_length, ranks)];
        while let Some((range, base, index, bits, ranks)) = stack.pop() {
            if bits == 0 {
                result.resize(
                    result.len() + ranks.len(),
                    self.compress.values()[index].clone(),
                );
                continue;
            }
            let d = bits - 1;
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            let split = base + end0 - start0;
            let boundary = ranks.partition_point(|&k| k < split);
            if boundary < ranks.len() {
                stack.push((
                    self.zeros[level] + start1..self.zeros[level] + end1,
                    split,
                    index | (1 << d),
                    d,
                    &ranks[boundary..],
                ));
            }
            if boundary != 0 {
                stack.push((start0..end0, base, index, d, &ranks[..boundary]));
            }
        }
        result
    }

pub fn quantiles_sorted(&self, range: Range<usize>, ranks: &[usize]) -> Vec<T>

Returns the requested order statistics of one range. ranks must be sorted in nondecreasing order and every rank must be less than the range length.

pub fn quantile_outer(&self, range: Range<usize>, k: usize) -> T

get k-th smallest value out of range

pub fn rank_lessthan(&self, val: T, range: Range<usize>) -> usize

the number of value less than val in range

pub fn rank_lessthan_batch( &self, queries: impl IntoIterator<Item = (T, Range<usize>)>, ) -> Vec<usize>

Returns the number of values below each query’s threshold.

fn quad_rank_lessthan_index( &self, index: usize, range: Range<usize>, start_level: usize, ) -> usize

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 1201)

    pub fn rank_lessthan(&self, val: T, range: Range<usize>) -> usize {
        let index = self.compress.index_lower_bound(&val);
        if index == self.compress.size() {
            return range.len();
        }
        if self.len >= 262144 && !self.quad_vectors.is_empty() {
            self.quad_rank_lessthan_index(index, range, 0)
        } else {
            let bits = usize::BITS as usize
                - self.compress.size().saturating_sub(1).leading_zeros() as usize;
            self.rank_lessthan_index(index, range, bits)
        }
    }

    /// Returns the number of values below each query's threshold.
    pub fn rank_lessthan_batch(
        &self,
        queries: impl IntoIterator<Item = (T, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                *state = [
                    range.start,
                    range.end,
                    self.compress.index_lower_bound(value),
                    0,
                ];
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len());
            result.extend(states[..queries.len()].iter().map(|state| state[3]));
        }
        result
    }

    fn quad_rank_lessthan_index(
        &self,
        index: usize,
        mut range: Range<usize>,
        start_level: usize,
    ) -> usize {
        let mut result = 0;
        for (level, vector) in self.quad_vectors.iter().enumerate().skip(start_level) {
            let d = (self.quad_vectors.len() - level - 1) * 2;
            if d + 2 <= index.trailing_zeros() as usize {
                break;
            }
            let digit = (index >> d) & 3;
            let ranks = |position: usize| {
                let block = &vector.blocks[position / 64];
                let mask = !(u64::MAX << (position % 64));
                let low = if digit & 1 == 0 { !block.lo } else { block.lo };
                let high = if digit & 2 == 0 { !block.hi } else { block.hi };
                let exact = block.rank[digit] as usize + (low & high & mask).count_ones() as usize;
                let less = match digit {
                    0 => 0,
                    1 => {
                        block.rank[0] as usize
                            + (!(block.lo | block.hi) & mask).count_ones() as usize
                    }
                    2 => {
                        block.rank[0] as usize
                            + block.rank[1] as usize
                            + (!block.hi & mask).count_ones() as usize
                    }
                    _ => position - exact,
                };
                (exact, less)
            };
            let (start, start_less) = ranks(range.start);
            let (end, end_less) = ranks(range.end);
            result += end_less - start_less;
            range = vector.starts[digit] + start..vector.starts[digit] + end;
        }
        result
    }

    fn rank_lessthan_index(&self, idx: usize, mut range: Range<usize>, bits: usize) -> usize {
        let mut res = 0;
        for d in (idx.trailing_zeros() as usize..bits).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let bit = (idx >> d) & 1 != 0;
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            res += if bit { end0 - start0 } else { 0 };
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        res
    }

    /// the number of valrange in range
    pub fn rank_range(&self, valrange: Range<T>, mut range: Range<usize>) -> usize {
        let lower = self.compress.index_lower_bound(&valrange.start);
        let upper = self.compress.index_lower_bound(&valrange.end);
        if lower >= upper {
            return 0;
        }
        if !self.quad_vectors.is_empty() {
            if upper == self.compress.size() {
                return range.len() - self.quad_rank_lessthan_index(lower, range, 0);
            }
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let low = (lower >> d) & 3;
                let high = (upper >> d) & 3;
                if low == high {
                    range = vector.starts[low] + vector.rank(low, range.start)
                        ..vector.starts[low] + vector.rank(low, range.end);
                    continue;
                }
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let count: usize = (low..high).map(|digit| end[digit] - start[digit]).sum();
                return count
                    - self.quad_rank_lessthan_index(
                        lower,
                        vector.starts[low] + start[low]..vector.starts[low] + end[low],
                        level + 1,
                    )
                    + self.quad_rank_lessthan_index(
                        upper,
                        vector.starts[high] + start[high]..vector.starts[high] + end[high],
                        level + 1,
                    );
            }
            return 0;
        }
        for d in (0..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let zero_range = range.start - start1..range.end - end1;
            let one_range = self.zeros[level] + start1..self.zeros[level] + end1;
            if ((lower ^ upper) >> d) & 1 != 0 {
                return zero_range.len() - self.rank_lessthan_index(lower, zero_range, d)
                    + self.rank_lessthan_index(upper, one_range, d);
            }
            range = if (lower >> d) & 1 != 0 {
                one_range
            } else {
                zero_range
            };
        }
        0
    }

fn rank_lessthan_index( &self, idx: usize, range: Range<usize>, bits: usize, ) -> usize

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (line 1205)

    pub fn rank_lessthan(&self, val: T, range: Range<usize>) -> usize {
        let index = self.compress.index_lower_bound(&val);
        if index == self.compress.size() {
            return range.len();
        }
        if self.len >= 262144 && !self.quad_vectors.is_empty() {
            self.quad_rank_lessthan_index(index, range, 0)
        } else {
            let bits = usize::BITS as usize
                - self.compress.size().saturating_sub(1).leading_zeros() as usize;
            self.rank_lessthan_index(index, range, bits)
        }
    }

    /// Returns the number of values below each query's threshold.
    pub fn rank_lessthan_batch(
        &self,
        queries: impl IntoIterator<Item = (T, Range<usize>)>,
    ) -> Vec<usize> {
        let queries: Vec<_> = queries.into_iter().collect();
        let mut result = Vec::with_capacity(queries.len());
        for queries in queries.chunks(16) {
            let mut states = [[0; 4]; 16];
            for (state, (value, range)) in states.iter_mut().zip(queries) {
                assert!(range.start <= range.end && range.end <= self.len);
                *state = [
                    range.start,
                    range.end,
                    self.compress.index_lower_bound(value),
                    0,
                ];
            }
            self.batch::<RANK_LESSTHAN>(&mut states, queries.len());
            result.extend(states[..queries.len()].iter().map(|state| state[3]));
        }
        result
    }

    fn quad_rank_lessthan_index(
        &self,
        index: usize,
        mut range: Range<usize>,
        start_level: usize,
    ) -> usize {
        let mut result = 0;
        for (level, vector) in self.quad_vectors.iter().enumerate().skip(start_level) {
            let d = (self.quad_vectors.len() - level - 1) * 2;
            if d + 2 <= index.trailing_zeros() as usize {
                break;
            }
            let digit = (index >> d) & 3;
            let ranks = |position: usize| {
                let block = &vector.blocks[position / 64];
                let mask = !(u64::MAX << (position % 64));
                let low = if digit & 1 == 0 { !block.lo } else { block.lo };
                let high = if digit & 2 == 0 { !block.hi } else { block.hi };
                let exact = block.rank[digit] as usize + (low & high & mask).count_ones() as usize;
                let less = match digit {
                    0 => 0,
                    1 => {
                        block.rank[0] as usize
                            + (!(block.lo | block.hi) & mask).count_ones() as usize
                    }
                    2 => {
                        block.rank[0] as usize
                            + block.rank[1] as usize
                            + (!block.hi & mask).count_ones() as usize
                    }
                    _ => position - exact,
                };
                (exact, less)
            };
            let (start, start_less) = ranks(range.start);
            let (end, end_less) = ranks(range.end);
            result += end_less - start_less;
            range = vector.starts[digit] + start..vector.starts[digit] + end;
        }
        result
    }

    fn rank_lessthan_index(&self, idx: usize, mut range: Range<usize>, bits: usize) -> usize {
        let mut res = 0;
        for d in (idx.trailing_zeros() as usize..bits).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let bit = (idx >> d) & 1 != 0;
            let start0 = range.start - start1;
            let end0 = range.end - end1;
            res += if bit { end0 - start0 } else { 0 };
            range.start = if bit {
                self.zeros[level] + start1
            } else {
                start0
            };
            range.end = if bit { self.zeros[level] + end1 } else { end0 };
        }
        res
    }

    /// the number of valrange in range
    pub fn rank_range(&self, valrange: Range<T>, mut range: Range<usize>) -> usize {
        let lower = self.compress.index_lower_bound(&valrange.start);
        let upper = self.compress.index_lower_bound(&valrange.end);
        if lower >= upper {
            return 0;
        }
        if !self.quad_vectors.is_empty() {
            if upper == self.compress.size() {
                return range.len() - self.quad_rank_lessthan_index(lower, range, 0);
            }
            for (level, vector) in self.quad_vectors.iter().enumerate() {
                let d = (self.quad_vectors.len() - level - 1) * 2;
                let low = (lower >> d) & 3;
                let high = (upper >> d) & 3;
                if low == high {
                    range = vector.starts[low] + vector.rank(low, range.start)
                        ..vector.starts[low] + vector.rank(low, range.end);
                    continue;
                }
                let start = vector.ranks(range.start);
                let end = vector.ranks(range.end);
                let count: usize = (low..high).map(|digit| end[digit] - start[digit]).sum();
                return count
                    - self.quad_rank_lessthan_index(
                        lower,
                        vector.starts[low] + start[low]..vector.starts[low] + end[low],
                        level + 1,
                    )
                    + self.quad_rank_lessthan_index(
                        upper,
                        vector.starts[high] + start[high]..vector.starts[high] + end[high],
                        level + 1,
                    );
            }
            return 0;
        }
        for d in (0..self.bit_length).rev() {
            let level = self.level(d);
            let start1 = self.rank1(level, range.start);
            let end1 = self.rank1(level, range.end);
            let zero_range = range.start - start1..range.end - end1;
            let one_range = self.zeros[level] + start1..self.zeros[level] + end1;
            if ((lower ^ upper) >> d) & 1 != 0 {
                return zero_range.len() - self.rank_lessthan_index(lower, zero_range, d)
                    + self.rank_lessthan_index(upper, one_range, d);
            }
            range = if (lower >> d) & 1 != 0 {
                one_range
            } else {
                zero_range
            };
        }
        0
    }

pub fn rank_range(&self, valrange: Range<T>, range: Range<usize>) -> usize

the number of valrange in range

pub fn rank_range_batch( &self, queries: impl IntoIterator<Item = (Range<T>, Range<usize>)>, ) -> Vec<usize>

Returns each query’s count in a half-open value range.

pub fn query_less_than<F>(&self, val: T, range: Range<usize>, f: F)

where F: FnMut(usize, Range<usize>),

Examples found in repository

crates/competitive/src/data_structure/wavelet_matrix.rs (lines 1505-1510)

    pub fn fold_lessthan(&self, value: T, range: Range<usize>) -> M::T {
        let mut result = M::unit();
        self.wavelet_matrix
            .query_less_than(value, range, |d, range| {
                M::operate_assign(
                    &mut result,
                    &self.bits[self.wavelet_matrix.level(d)].fold_abelian(range.start, range.end),
                );
            });
        result
    }

pub fn build_fold<M>(&self, weights: &[M::T]) -> WaveletMatrixFold<'_, T, M>

where M: AbelianGroup,

Examples found in repository

crates/library_checker/src/data_structure/static_range_sum_with_upper_bound.rs (line 27)

pub fn static_range_sum_with_upper_bound_wavelet_matrix(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, a: [i64; n]);
    let weights = a.clone();
    let wm = WaveletMatrix::new(a);
    let fold = wm.build_fold::<AdditiveOperation<i64>>(&weights);
    for _ in 0..q {
        sc!(l, r, x: i64);
        let (count, sum) = fold.fold_lessthan_with_count(x + 1, l..r);
        pp!(count, sum);
    }
}

crates/library_checker/src/data_structure/rectangle_sum.rs (line 21)

pub fn rectangle_sum(reader: impl Read, writer: impl Write) {
    prepare_io!(buffered; reader, writer);
    sc!(n, q, mut xyw: [(u32, u32, i64); n], queries: [(u32, u32, u32, u32); q]);
    xyw.radix_sort_by_key(|&(x, ..)| x);
    let xs: Vec<_> = xyw.iter().map(|&(x, ..)| x).collect();
    let search = StaticSearch::from_sorted(&xs);
    let endpoints: Vec<_> = queries.iter().flat_map(|&(l, _, r, _)| [l, r]).collect();
    let mut positions = vec![0; endpoints.len()];
    search.lower_bound_batch(&endpoints, &mut positions);
    let ys = xyw.iter().map(|&(_, y, _)| y).collect();
    let weights: Vec<_> = xyw.iter().map(|&(_, _, w)| w).collect();
    let wm = WaveletMatrix::new(ys);
    let fold = wm.build_fold::<AdditiveOperation<i64>>(&weights);
    let result = fold.fold_lessthan_batch(
        queries
            .into_iter()
            .zip(positions.as_chunks::<2>().0)
            .flat_map(|((_, d, _, u), &[l, r])| [(d, l..r), (u, l..r)]),
    );
    for &[lower, upper] in result.as_chunks::<2>().0 {
        pp!(upper - lower);
    }
}

pub fn build_point_add<M>( &self, weights: &[M::T], ) -> WaveletMatrixPointAdd<'_, T, M>

where M: AbelianGroup,

Examples found in repository

crates/library_checker/src/data_structure/point_add_rectangle_sum.rs (line 40)

pub fn point_add_rectangle_sum(reader: impl Read, writer: impl Write) {
    prepare_io!(reader, writer);
    sc!(n, q, xyw: [(u32, u32, u64); iter n]);
    let mut points: Vec<_> = xyw.map(|(x, y, w)| (x, y, w as i64)).collect();
    sc!(queries: [Query; q]);
    points.extend(queries.iter().filter_map(|&query| match query {
        Query::Add { x, y, .. } => Some((x, y, 0)),
        Query::Sum { .. } => None,
    }));
    let mut order: Vec<_> = (0..points.len()).collect();
    order.radix_sort_by_key(|&i| points[i].0);
    let mut positions = vec![0; points.len()];
    let mut xs = Vec::with_capacity(points.len());
    let mut ys = Vec::with_capacity(points.len());
    let mut weights = Vec::with_capacity(points.len());
    for (i, &point) in order.iter().enumerate() {
        positions[point] = i;
        let (x, y, w) = points[point];
        xs.push(x);
        ys.push(y);
        weights.push(w);
    }
    let wm = WaveletMatrix::new(ys);
    let mut fold: WaveletMatrixPointAdd<_, AdditiveOperation<i64>> = wm.build_point_add(&weights);

    let mut point = n;
    for query in queries {
        match query {
            Query::Add { w, .. } => {
                fold.update(positions[point], w as i64);
                point += 1;
            }
            Query::Sum { l, d, r, u } => {
                let l = xs.partition_point(|&x| x < l);
                let r = xs.partition_point(|&x| x < r);
                pp!(fold.fold_range(d..u, l..r));
            }
        }
    }
}

Trait Implementations

impl<T: Clone> Clone for WaveletMatrix<T>

fn clone(&self) -> Self

Returns a duplicate of the value.

fn clone_from(&mut self, source: &Self)

Performs copy-assignment from source.

impl<T: Debug> Debug for WaveletMatrix<T>

fn fmt(&self, f: &mut Formatter<'_>) -> Result

Formats the value using the given formatter.