roaring 0.12.0

//! Ported from CRoaring and arXiv:1709.07821
//! Lemire et al, Roaring Bitmaps: Implementation of an Optimized Software Library
//!
//! Prior work: Schlegel et al., Fast Sorted-Set Intersection using SIMD Instructions
//!
//! Rust port notes:
//! Like CRoaring, the intersection and difference use the x86 PCMPISTRM instruction (falling back
//! to PCMPESTRM while a zero may be present) where SSE4.2 is available. Other SIMD levels use a
//! portable all-pairs comparison (see `matrix_cmp_u16`) instead.
//!
//! The public functions select a SIMD level at runtime with `fearless_simd::dispatch!`, and call
//! implementations generic over the SIMD level. Those are annotated with `#[simd]`, so that they
//! are compiled with the matching target features enabled.
//!
//! The small helpers they call are `#[inline(always)]` rather than `#[simd]`: they get the target
//! features by being inlined, and `#[simd]` does not force inlining, which can leave them out of
//! line in the hot loops.

use super::scalar;
use crate::bitmap::store::array_store::visitor::BinaryOperationVisitor;
use core::hint::select_unpredictable;
use fearless_simd::prelude::*;
use fearless_simd::{dispatch, mask16x8, u16x8, Level};
use fearless_simd_macros::simd;

/// The number of lanes in a `u16x8`
const LANES: usize = 8;

/// The SIMD level to use.
///
/// With `std`, this is detected at runtime on x86/x86-64 (and cached by `fearless_simd`).
/// Otherwise, it is determined by the target features enabled at compile time.
#[inline]
fn level() -> Level {
    #[cfg(feature = "std")]
    return Level::new();
    #[cfg(not(feature = "std"))]
    return Level::baseline();
}

// a one-pass union algorithm
pub fn or(lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    dispatch!(level(), simd => or_impl(simd, lhs, rhs, visitor))
}

#[simd]
fn or_impl<S: Simd>(simd: S, lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    // De-duplicates `slice` in place
    // Returns the end index of the deduplicated slice.
    // elements after the return value are not guaranteed to be unique or in order
    #[inline]
    fn dedup(slice: &mut [u16]) -> usize {
        let mut pos: usize = 1;
        for i in 1..slice.len() {
            if slice[i] != slice[i - 1] {
                slice[pos] = slice[i];
                pos += 1;
            }
        }
        pos
    }

    // returns `new`, with a mask of its values which were not already written,
    // assuming that the previously written vector was `old`
    #[inline(always)]
    fn handle_vector<S: Simd>(old: u16x8<S>, new: u16x8<S>) -> (u16x8<S>, mask16x8<S>) {
        // `[old[7], new[0], ..., new[6]]`
        let tmp: u16x8<S> = old.slide::<7>(new);
        (new, tmp.simd_ne(new))
    }

    let (lhs_chunks, lhs_remain) = lhs.as_chunks::<LANES>();
    let (rhs_chunks, rhs_remain) = rhs.as_chunks::<LANES>();
    let ([lhs_first, ..], [rhs_first, ..]) = (lhs_chunks, rhs_chunks) else {
        scalar::or(lhs, rhs, visitor);
        return;
    };

    let v_a = u16x8::load_array_ref(simd, lhs_first);
    let v_b = u16x8::load_array_ref(simd, rhs_first);
    let (mut v_min, mut v_max) = bitonic_merge(v_a, v_b);

    let mut i = 1;
    let mut j = 1;
    let (v, m) = handle_vector(u16x8::splat(simd, u16::MAX), v_min);
    visitor.visit_vector(v, m);
    let mut v_prev = v_min;
    while let (Some(lhs_chunk), Some(rhs_chunk)) = (lhs_chunks.get(i), rhs_chunks.get(j)) {
        let take_lhs = lhs_chunk[0] <= rhs_chunk[0];
        let (src, side) = select_unpredictable(take_lhs, (lhs_chunk, &mut i), (rhs_chunk, &mut j));
        let v = u16x8::load_array_ref(simd, src);
        *side += 1;
        (v_min, v_max) = bitonic_merge(v_max, v);
        let (v, m) = handle_vector(v_prev, v_min);
        visitor.visit_vector(v, m);
        v_prev = v_min;
    }

    debug_assert!(i == lhs_chunks.len() || j == rhs_chunks.len());

    // we finish the rest off using a scalar algorithm
    // could be improved?
    //
    // copy the small end on a tmp buffer
    let mut buffer: [u16; 16] = [0; 16];
    let (v, m) = handle_vector(v_prev, v_max);
    v.compress(m).store_array(buffer.first_chunk_mut().unwrap());
    let mut rem = count_set_lanes(m);

    let (tail_a, tail_b) = if i == lhs_chunks.len() {
        (lhs_remain, &rhs[LANES * j..])
    } else {
        (rhs_remain, &lhs[LANES * i..])
    };

    buffer[rem..rem + tail_a.len()].copy_from_slice(tail_a);
    rem += tail_a.len();

    if rem == 0 {
        visitor.visit_slice(tail_b)
    } else {
        buffer[..rem].sort_unstable();
        rem = dedup(&mut buffer[..rem]);
        scalar::or(&buffer[..rem], tail_b, visitor);
    }
}

pub fn and(lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    dispatch!(level(), simd => and_impl(simd, lhs, rhs, visitor))
}

#[simd]
fn and_impl<S: Simd>(simd: S, lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    let (lhs_chunks, _) = lhs.as_chunks::<LANES>();
    let (rhs_chunks, _) = rhs.as_chunks::<LANES>();

    let mut i = 0;
    let mut j = 0;
    if let ([lhs_first, ..], [rhs_first, ..]) = (lhs_chunks, rhs_chunks) {
        let mut lhs_chunk = lhs_first;
        let mut rhs_chunk = rhs_first;
        let mut v_a = u16x8::load_array_ref(simd, lhs_chunk);
        let mut v_b = u16x8::load_array_ref(simd, rhs_chunk);
        loop {
            let mask = matrix_cmp_mask(v_a, v_b, lhs_chunk[0] == 0 || rhs_chunk[0] == 0);
            visitor.visit_vector(v_a, mask);

            let a_max = *lhs_chunk.last().unwrap();
            let b_max = *rhs_chunk.last().unwrap();
            if a_max <= b_max {
                i += 1;
                match lhs_chunks.get(i) {
                    Some(next) => lhs_chunk = next,
                    None => break,
                }
                v_a = u16x8::load_array_ref(simd, lhs_chunk);
            }
            if b_max <= a_max {
                j += 1;
                match rhs_chunks.get(j) {
                    Some(next) => rhs_chunk = next,
                    None => break,
                }
                v_b = u16x8::load_array_ref(simd, rhs_chunk);
            }
        }
    }

    // intersect the tail using scalar intersection
    scalar::and(&lhs[LANES * i..], &rhs[LANES * j..], visitor);
}

// a one-pass xor algorithm
pub fn xor(lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    dispatch!(level(), simd => xor_impl(simd, lhs, rhs, visitor))
}

#[simd]
fn xor_impl<S: Simd>(simd: S, lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    /// De-duplicates `slice` in place, removing _both_ duplicates
    /// Returns the end index of the xor-ed slice.
    /// elements after the return value are not guaranteed to be unique or in order
    #[inline]
    fn xor_slice(slice: &mut [u16]) -> usize {
        let mut pos: usize = 1;
        for i in 1..slice.len() {
            if slice[i] != slice[i - 1] {
                slice[pos] = slice[i];
                pos += 1;
            } else {
                pos -= 1; // it is identical to previous, delete it
            }
        }
        pos
    }

    // returns the vector to write, and a mask of the values to write from it,
    // omitting repeated values assuming that previously written vector was "old"
    #[inline(always)]
    fn handle_vector<S: Simd>(old: u16x8<S>, new: u16x8<S>) -> (u16x8<S>, mask16x8<S>) {
        // `[old[6], old[7], new[0], ..., new[5]]`
        let tmp1: u16x8<S> = old.slide::<6>(new);
        // `[old[7], new[0], ..., new[6]]`
        let tmp2: u16x8<S> = old.slide::<7>(new);
        let ne_l: mask16x8<S> = tmp2.simd_ne(tmp1);
        let ne_r: mask16x8<S> = tmp2.simd_ne(new);
        (tmp2, ne_l & ne_r)
    }

    let (lhs_chunks, lhs_remain) = lhs.as_chunks::<LANES>();
    let (rhs_chunks, rhs_remain) = rhs.as_chunks::<LANES>();
    let ([lhs_first, ..], [rhs_first, ..]) = (lhs_chunks, rhs_chunks) else {
        scalar::xor(lhs, rhs, visitor);
        return;
    };

    let v_a = u16x8::load_array_ref(simd, lhs_first);
    let v_b = u16x8::load_array_ref(simd, rhs_first);
    let (mut v_min, mut v_max) = bitonic_merge(v_a, v_b);

    let mut i = 1;
    let mut j = 1;
    let (v, m) = handle_vector(u16x8::splat(simd, u16::MAX), v_min);
    visitor.visit_vector(v, m);
    let mut v_prev = v_min;
    while let (Some(lhs_chunk), Some(rhs_chunk)) = (lhs_chunks.get(i), rhs_chunks.get(j)) {
        let take_lhs = lhs_chunk[0] <= rhs_chunk[0];
        let (src, side) = select_unpredictable(take_lhs, (lhs_chunk, &mut i), (rhs_chunk, &mut j));
        let v = u16x8::load_array_ref(simd, src);
        *side += 1;
        (v_min, v_max) = bitonic_merge(v_max, v);
        let (v, m) = handle_vector(v_prev, v_min);
        visitor.visit_vector(v, m);
        v_prev = v_min;
    }

    debug_assert!(i == lhs_chunks.len() || j == rhs_chunks.len());

    // we finish the rest off using a scalar algorithm
    // could be improved?
    // conditionally stores the last value of laststore as well as all but the
    // last value of vecMax,
    let mut buffer: [u16; 17] = [0; 17];
    // remaining size
    let (v, m) = handle_vector(v_prev, v_max);
    v.compress(m).store_array(buffer.first_chunk_mut().unwrap());
    let mut rem = count_set_lanes(m);

    // Store `v_max` rather than reading its lanes with indexing, which makes LLVM split
    // `v_max` into pieces in the merge loop above
    let mut arr_max = [0u16; LANES];
    v_max.store_array(&mut arr_max);
    let vec7 = arr_max[7];
    let vec6 = arr_max[6];
    if vec6 != vec7 {
        buffer[rem] = vec7;
        rem += 1;
    }

    let (tail_a, tail_b) = if i == lhs_chunks.len() {
        (lhs_remain, &rhs[LANES * j..])
    } else {
        (rhs_remain, &lhs[LANES * i..])
    };

    buffer[rem..rem + tail_a.len()].copy_from_slice(tail_a);
    rem += tail_a.len();

    if rem == 0 {
        visitor.visit_slice(tail_b)
    } else {
        buffer[..rem].sort_unstable();
        rem = xor_slice(&mut buffer[..rem]);
        scalar::xor(&buffer[..rem], tail_b, visitor);
    }
}

pub fn sub(lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    dispatch!(level(), simd => sub_impl(simd, lhs, rhs, visitor))
}

#[simd]
fn sub_impl<S: Simd>(simd: S, lhs: &[u16], rhs: &[u16], visitor: &mut impl BinaryOperationVisitor) {
    // we handle the degenerate cases
    if lhs.is_empty() {
        return;
    } else if rhs.is_empty() {
        visitor.visit_slice(lhs);
        return;
    }

    let (lhs_chunks, _) = lhs.as_chunks::<LANES>();
    let (rhs_chunks, _) = rhs.as_chunks::<LANES>();

    let mut i = 0;
    let mut j = 0;
    let mut rhs_tail = rhs;
    if let ([lhs_first, ..], [rhs_first, ..]) = (lhs_chunks, rhs_chunks) {
        let mut lhs_chunk = lhs_first;
        let mut rhs_chunk = rhs_first;
        let mut v_a = u16x8::load_array_ref(simd, lhs_chunk);
        let mut v_b = u16x8::load_array_ref(simd, rhs_chunk);
        // we have a running mask which indicates which values from a have been
        // spotted in b, these don't get written out.
        let mut runningmask_a_found_in_b = mask16x8::splat(simd, false);
        loop {
            // a_found_in_b will contain a mask indicate for each entry in A
            // whether it is seen in B
            let a_found_in_b = matrix_cmp_mask(v_a, v_b, lhs_chunk[0] == 0 || rhs_chunk[0] == 0);
            runningmask_a_found_in_b |= a_found_in_b;
            // we always compare the last values of A and B
            let a_max = *lhs_chunk.last().unwrap();
            let b_max = *rhs_chunk.last().unwrap();
            if a_max <= b_max {
                // Ok. In this code path, we are ready to write our v_a
                // because there is no need to read more from B, they will
                // all be large values.
                let bitmask_belongs_to_difference = !runningmask_a_found_in_b;
                visitor.visit_vector(v_a, bitmask_belongs_to_difference);
                i += 1;
                match lhs_chunks.get(i) {
                    Some(next) => lhs_chunk = next,
                    None => break,
                }
                runningmask_a_found_in_b = mask16x8::splat(simd, false);
                v_a = u16x8::load_array_ref(simd, lhs_chunk);
            }
            if b_max <= a_max {
                // in this code path, the current v_b has become useless
                j += 1;
                match rhs_chunks.get(j) {
                    Some(next) => rhs_chunk = next,
                    None => break,
                }
                v_b = u16x8::load_array_ref(simd, rhs_chunk);
            }
        }

        debug_assert!(i == lhs_chunks.len() || j == rhs_chunks.len());
        rhs_tail = &rhs[LANES * j..];

        // End of main vectorized loop
        // At this point either all of lhs's chunks are processed, which is the end of the
        // vectorized processing, or all of rhs's chunks are and we are not done processing
        // `lhs_chunk`... so we need to finish it off.
        if i < lhs_chunks.len() {
            if !rhs_tail.is_empty() {
                let mut buffer: [u16; 8] = [0; 8]; // buffer to do a masked load
                buffer[..rhs_tail.len()].copy_from_slice(rhs_tail);
                // Ensure the buffer is filled with a value we should remove: we do not want to
                // end up trying to remove zero values which aren't actually in rhs
                buffer[rhs_tail.len()..].fill(rhs_tail[0]);
                v_b = u16x8::load_array_ref(simd, &buffer);
                let a_found_in_b = matrix_cmp_mask(v_a, v_b, lhs_chunk[0] == 0 || rhs_tail[0] == 0);
                runningmask_a_found_in_b |= a_found_in_b;
                // Read from `lhs_chunk` (which `v_a` was loaded from): reading a lane of `v_a`
                // makes LLVM split `v_a` into pieces in the loop above
                let max_va = *lhs_chunk.last().unwrap();
                let used_rhs = rhs_tail.partition_point(|&b| b <= max_va);
                rhs_tail = &rhs_tail[used_rhs..];
            }
            let bitmask_belongs_to_difference = !runningmask_a_found_in_b;
            visitor.visit_vector(v_a, bitmask_belongs_to_difference);
            i += 1;
        }
    }

    // do the tail using scalar code
    scalar::sub(&lhs[LANES * i..], rhs_tail, visitor);
}

/// Compare all lanes in `a` to all lanes in `b`
///
/// Returns result mask will be set if any lane at `a[i]` is in any lane of `b`
///
/// ### Example
/// ```ignore
/// let a = u16x8::from_slice(simd, &[1, 2, 3, 4, 32, 33, 34, 35]);
/// let b = u16x8::from_slice(simd, &[2, 4, 6, 8, 10, 12, 14, 16]);
/// let result = matrix_cmp_u16(a, b);
/// assert_eq!(result.to_bitmask(), 0b0000_1010);
/// ```
#[inline(always)]
fn matrix_cmp_u16<S: Simd>(a: u16x8<S>, b: u16x8<S>) -> mask16x8<S> {
    a.simd_eq(b)
        | a.simd_eq(b.rotate_elements_left::<1>())
        | a.simd_eq(b.rotate_elements_left::<2>())
        | a.simd_eq(b.rotate_elements_left::<3>())
        | a.simd_eq(b.rotate_elements_left::<4>())
        | a.simd_eq(b.rotate_elements_left::<5>())
        | a.simd_eq(b.rotate_elements_left::<6>())
        | a.simd_eq(b.rotate_elements_left::<7>())
}

/// Like [`matrix_cmp_u16`], with a hint for when PCMPISTRM can be used
///
/// `may_contain_zero` must be true if either vector may contain a zero lane: because the vectors
/// come from sorted arrays, that can only be the first lane of the first vector of either array.
///
/// Uses the SSE4.2 PCMPISTRM instruction where it is available, falling back to the slower
/// PCMPESTRM when there may be a zero (which PCMPISTRM treats as the end of the string).
#[inline(always)]
fn matrix_cmp_mask<S: Simd>(a: u16x8<S>, b: u16x8<S>, may_contain_zero: bool) -> mask16x8<S> {
    if !may_contain_zero {
        debug_assert!(a[0] != 0 && b[0] != 0)
    }
    #[cfg(any(target_arch = "x86", target_arch = "x86_64"))]
    if let Some(sse4_2) = a.token().level().as_sse4_2() {
        return x86::matrix_cmp_mask(sse4_2, a.into(), b.into(), may_contain_zero)
            .simd_into(a.simd);
    }
    matrix_cmp_u16(a, b)
}

#[cfg(any(target_arch = "x86", target_arch = "x86_64"))]
mod x86 {
    #[cfg(target_arch = "x86")]
    use core::arch::x86::*;
    #[cfg(target_arch = "x86_64")]
    use core::arch::x86_64::*;

    fearless_simd::kernel!(
        /// Returns a mask with lane `i` set if lane `i` of `a` is equal to any lane of `b`
        #[inline]
        pub(super) fn matrix_cmp_mask(
            sse4_2: Sse4_2,
            a: __m128i,
            b: __m128i,
            may_contain_zero: bool,
        ) -> __m128i {
            const MODE: i32 = _SIDD_UWORD_OPS | _SIDD_CMP_EQUAL_ANY | _SIDD_UNIT_MASK;
            if may_contain_zero {
                _mm_cmpestrm::<MODE>(b, 8, a, 8)
            } else {
                _mm_cmpistrm::<MODE>(b, a)
            }
        }
    );
}

/// Returns the number of lanes set in `mask`
#[inline(always)]
pub fn count_set_lanes<S: Simd>(mask: mask16x8<S>) -> usize {
    // Masking off the high bits lets LLVM prove the count is at most `LANES`, but for some
    // `to_bitmask` implementations (e.g. NEON) LLVM can't see that without help.
    (mask.to_bitmask() & 0xff).count_ones() as usize
}

/// Assuming that a and b are sorted, returns an array of the sorted output.
#[inline(always)]
fn bitonic_merge<S: Simd, V: SimdInt<S>>(mut a: V, mut b: V) -> (V, V) {
    const { assert!(V::LEN.is_power_of_two()) }

    // Conceptually, a ++ reverse(b) is a bitonic sequence.
    b = b.reverse();

    // Each bitonic stage compare-exchanges positions that differ in one index
    // bit, from the most significant down. Lane-wise min/max between `a` and
    // `b` handles the top bit (which vector a value is in). Interleaving is a
    // perfect shuffle: it rotates the index bits left by one, moving the next
    // bit to compare into the top position. After log2(2 * LEN) stages, the
    // bits have rotated all the way around, so the lanes are back in
    // ascending order.
    for _ in 0..V::LEN.trailing_zeros() + 1 {
        let lo = a.min(b);
        let hi = a.max(b);
        (a, b) = lo.interleave(hi);
    }

    (a, b)
}

#[cfg(test)]
mod test {
    use super::*;
    use crate::bitmap::store::array_store::visitor::VecWriter;
    use alloc::vec::Vec;
    use proptest::prelude::*;

    fn check_bitonic_merge<S: Simd, V: SimdInt<S>>(
        simd: S,
        lhs: &[V::Element],
        rhs: &[V::Element],
    ) {
        let mut expected = [lhs, rhs].concat();
        expected.sort_unstable();

        let (lo, hi) = bitonic_merge(V::from_slice(simd, lhs), V::from_slice(simd, rhs));
        assert_eq!(lo.as_slice(), &expected[..V::LEN], "lhs: {lhs:?}, rhs: {rhs:?}");
        assert_eq!(hi.as_slice(), &expected[V::LEN..], "lhs: {lhs:?}, rhs: {rhs:?}");
    }

    #[test]
    fn bitonic_merge_handles_every_distinct_lane_interleaving() {
        dispatch!(level(), simd => {
            for mask in 0..=u16::MAX {
                if mask.count_ones() as usize != LANES {
                    continue;
                }
                let mut lhs = Vec::with_capacity(LANES);
                let mut rhs = Vec::with_capacity(LANES);
                for value in 0..16 {
                    if mask & (1 << value) != 0 {
                        lhs.push(value);
                    } else {
                        rhs.push(value);
                    }
                }
                check_bitonic_merge::<_, u16x8<_>>(simd, &lhs, &rhs);
            }
        });
    }

    /// Checks the vectorized op produces the same result as the scalar op
    fn check_op(
        vector: impl Fn(&[u16], &[u16], &mut VecWriter),
        scalar: impl Fn(&[u16], &[u16], &mut VecWriter),
        lhs: &[u16],
        rhs: &[u16],
    ) {
        let mut expected = VecWriter::new(lhs.len() + rhs.len());
        scalar(lhs, rhs, &mut expected);
        let mut actual = VecWriter::new(lhs.len() + rhs.len());
        vector(lhs, rhs, &mut actual);
        assert_eq!(actual.into_inner(), expected.into_inner());
    }

    /// Checks all vectorized ops produce the same results as the scalar ops
    fn check_all(lhs: &[u16], rhs: &[u16]) {
        check_op(or, scalar::or, lhs, rhs);
        check_op(and, scalar::and, lhs, rhs);
        check_op(xor, scalar::xor, lhs, rhs);
        check_op(sub, scalar::sub, lhs, rhs);
    }

    proptest! {
        #[test]
        fn bitonic_merge_matches_sorted_concatenation(
            mut lhs in prop::array::uniform8(prop_oneof![any::<u16>(), 0u16..16]),
            mut rhs in prop::array::uniform8(prop_oneof![any::<u16>(), 0u16..16]),
        ) {
            lhs.sort_unstable();
            rhs.sort_unstable();
            dispatch!(level(), simd => {
                check_bitonic_merge::<_, u16x8<_>>(simd, &lhs, &rhs);
            });
        }

        #[test]
        fn vector_ops_match_scalar_dense(
            lhs in prop::collection::btree_set(0u16..256, 0..200),
            rhs in prop::collection::btree_set(0u16..256, 0..200),
        ) {
            check_all(&Vec::from_iter(lhs), &Vec::from_iter(rhs));
        }

        #[test]
        fn vector_ops_match_scalar_sparse(
            lhs in prop::collection::btree_set(any::<u16>(), 0..200),
            rhs in prop::collection::btree_set(any::<u16>(), 0..200),
        ) {
            check_all(&Vec::from_iter(lhs), &Vec::from_iter(rhs));
        }
    }

    #[test]
    fn vector_ops_match_scalar_at_lane_boundaries() {
        // Lengths around multiples of `LANES` exercise the scalar fallback, the vector loops,
        // and the partial vector tails; the offsets give disjoint, interleaved, and equal inputs.
        let lens = [0, 1, 7, 8, 9, 15, 16, 17, 24, 25];
        for len_l in lens {
            for len_r in lens {
                for offset in [0, 1, 4, 8, 64] {
                    let lhs: Vec<u16> = (0..len_l).map(|x| 2 * x).collect();
                    let rhs: Vec<u16> = (0..len_r).map(|x| 2 * x + offset).collect();
                    check_all(&lhs, &rhs);
                }
            }
        }
    }

    #[test]
    fn matrix_cmp_finds_lanes_of_a_in_b() {
        dispatch!(level(), simd => {
            let a = u16x8::from_slice(simd, &[1, 2, 3, 4, 32, 33, 34, 35]);
            let b = u16x8::from_slice(simd, &[2, 4, 6, 8, 10, 12, 14, 35]);
            assert_eq!(matrix_cmp_u16(a, b).to_bitmask(), 0b1000_1010);
            assert_eq!(matrix_cmp_mask(a, b, false).to_bitmask(), 0b1000_1010);
            assert_eq!(matrix_cmp_mask(a, b, true).to_bitmask(), 0b1000_1010);
            let a = u16x8::from_slice(simd, &[0, 2, 3, 4, 32, 33, 34, 35]);
            let b = u16x8::from_slice(simd, &[0, 4, 6, 8, 10, 12, 14, 35]);
            assert_eq!(matrix_cmp_mask(a, b, true).to_bitmask(), 0b1000_1001);
        });
    }
}