Selection.java
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.numbers.arrays;
/**
* Select indices in array data.
*
* <p>Arranges elements such that indices {@code k} correspond to their correctly
* sorted value in the equivalent fully sorted array. For all indices {@code k}
* and any index {@code i}:
*
* <pre>{@code
* data[i < k] <= data[k] <= data[k < i]
* }</pre>
*
* <p>Examples:
*
* <pre>
* data [0, 1, 2, 1, 2, 5, 2, 3, 3, 6, 7, 7, 7, 7]
*
*
* k=4 : [0, 1, 2, 1], [2], [5, 2, 3, 3, 6, 7, 7, 7, 7]
* k=4,8 : [0, 1, 2, 1], [2], [3, 3, 2], [5], [6, 7, 7, 7, 7]
* </pre>
*
* <p>This implementation can select on multiple indices and will handle duplicate and
* unordered indices. The method detects ordered indices (with or without duplicates) and
* uses this during processing. Passing ordered indices is recommended if the order is already
* known; for example using uniform spacing through the array data, or to select the top and
* bottom {@code n} values from the data.
*
* <p>A quickselect adaptive method is used for single indices. This uses analysis of the
* partition sizes after each division to update the algorithm mode. If the partition
* containing the target does not sufficiently reduce in size then the algorithm is
* progressively changed to use partitions with guaranteed margins. This ensures a set fraction
* of data is eliminated each step and worse-case linear run time performance. This method can
* handle a range of indices {@code [ka, kb]} with a small separation by targeting the start of
* the range {@code ka} and then selecting the remaining elements {@code (ka, kb]} that are at
* the edge of the partition bounded by {@code ka}.
*
* <p>Multiple keys are partitioned collectively using an introsort method which only recurses
* into partitions containing indices. Excess recursion will trigger use of a heapselect
* on the remaining range of indices ensuring non-quadratic worse case performance. Any
* partition containing a single index, adjacent pair of indices, or range of indices with a
* small separation will use the quickselect adaptive method for single keys. Note that the
* maximum number of times that {@code n} indices can be split is {@code n - 1} before all
* indices are handled as singles.
*
* <p>Floating-point order
*
* <p>The {@code <} relation does not impose a total order on all floating-point values.
* This class respects the ordering imposed by {@link Double#compare(double, double)}.
* {@code -0.0} is treated as less than value {@code 0.0}; {@code Double.NaN} is
* considered greater than any other value; and all {@code Double.NaN} values are
* considered equal.
*
* <p>References
*
* <p>Quickselect is introduced in Hoare [1]. This selects an element {@code k} from {@code n}
* using repeat division of the data around a partition element, recursing into the
* partition that contains {@code k}.
*
* <p>Introsort/select is introduced in Musser [2]. This detects excess recursion in
* quicksort/select and reverts to a heapsort or linear select to achieve an improved worst
* case bound.
*
* <p>Use of sampling to identify a pivot that places {@code k} in the smaller partition is
* performed in the SELECT algorithm of Floyd and Rivest [3, 4].
*
* <p>A worst-case linear time algorithm PICK is described in Blum <i>et al</i> [5]. This uses
* the median of medians as a partition element for selection which ensures a minimum fraction of
* the elements are eliminated per iteration. This was extended to use an asymmetric pivot choice
* with efficient placement of the medians sample location in the QuickselectAdpative algorithm of
* Alexandrescu [6].
*
* <ol>
* <li>Hoare (1961)
* Algorithm 65: Find
* <a href="https://doi.org/10.1145%2F366622.366647">Comm. ACM. 4 (7): 321–322</a>
* <li>Musser (1999)
* Introspective Sorting and Selection Algorithms
* <a href="https://doi.org/10.1002/(SICI)1097-024X(199708)27:8%3C983::AID-SPE117%3E3.0.CO;2-%23">
* Software: Practice and Experience 27, 983-993.</a>
* <li>Floyd and Rivest (1975)
* Algorithm 489: The Algorithm SELECT—for Finding the ith Smallest of n elements.
* Comm. ACM. 18 (3): 173.
* <li>Kiwiel (2005)
* On Floyd and Rivest's SELECT algorithm.
* Theoretical Computer Science 347, 214-238.
* <li>Blum, Floyd, Pratt, Rivest, and Tarjan (1973)
* Time bounds for selection.
* <a href="https://doi.org/10.1016%2FS0022-0000%2873%2980033-9">
* Journal of Computer and System Sciences. 7 (4): 448–461</a>.
* <li>Alexandrescu (2016)
* Fast Deterministic Selection
* <a href="https://arxiv.org/abs/1606.00484">arXiv:1606.00484</a>.
* <li><a href="https://en.wikipedia.org/wiki/Quickselect">Quickselect (Wikipedia)</a>
* <li><a href="https://en.wikipedia.org/wiki/Introsort">Introsort (Wikipedia)</a>
* <li><a href="https://en.wikipedia.org/wiki/Introselect">Introselect (Wikipedia)</a>
* <li><a href="https://en.wikipedia.org/wiki/Floyd%E2%80%93Rivest_algorithm">Floyd-Rivest algorithm (Wikipedia)</a>
* <li><a href="https://en.wikipedia.org/wiki/Median_of_medians">Median of medians (Wikipedia)</a>
* </ol>
*
* @since 1.2
*/
public final class Selection {
/** No instances. */
private Selection() {}
/**
* Partition the array such that index {@code k} corresponds to its correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param k Index.
* @throws IndexOutOfBoundsException if index {@code k} is not within the
* sub-range {@code [0, a.length)}
*/
public static void select(double[] a, int k) {
IndexSupport.checkIndex(0, a.length, k);
doSelect(a, 0, a.length, k);
}
/**
* Partition the array such that indices {@code k} correspond to their correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param k Indices (may be destructively modified).
* @throws IndexOutOfBoundsException if any index {@code k} is not within the
* sub-range {@code [0, a.length)}
*/
public static void select(double[] a, int[] k) {
IndexSupport.checkIndices(0, a.length, k);
doSelect(a, 0, a.length, k);
}
/**
* Partition the array such that index {@code k} corresponds to its correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param fromIndex Index of the first element (inclusive).
* @param toIndex Index of the last element (exclusive).
* @param k Index.
* @throws IndexOutOfBoundsException if the sub-range {@code [fromIndex, toIndex)} is out of
* bounds of range {@code [0, a.length)}; or if index {@code k} is not within the
* sub-range {@code [fromIndex, toIndex)}
*/
public static void select(double[] a, int fromIndex, int toIndex, int k) {
IndexSupport.checkFromToIndex(fromIndex, toIndex, a.length);
IndexSupport.checkIndex(fromIndex, toIndex, k);
doSelect(a, fromIndex, toIndex, k);
}
/**
* Partition the array such that indices {@code k} correspond to their correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param fromIndex Index of the first element (inclusive).
* @param toIndex Index of the last element (exclusive).
* @param k Indices (may be destructively modified).
* @throws IndexOutOfBoundsException if the sub-range {@code [fromIndex, toIndex)} is out of
* bounds of range {@code [0, a.length)}; or if any index {@code k} is not within the
* sub-range {@code [fromIndex, toIndex)}
*/
public static void select(double[] a, int fromIndex, int toIndex, int[] k) {
IndexSupport.checkFromToIndex(fromIndex, toIndex, a.length);
IndexSupport.checkIndices(fromIndex, toIndex, k);
doSelect(a, fromIndex, toIndex, k);
}
/**
* Partition the array such that index {@code k} corresponds to its correctly
* sorted value in the equivalent fully sorted array.
*
* <p>This method pre/post-processes the data and indices to respect the ordering
* imposed by {@link Double#compare(double, double)}.
*
* @param fromIndex Index of the first element (inclusive).
* @param toIndex Index of the last element (exclusive).
* @param a Values.
* @param k Index.
*/
private static void doSelect(double[] a, int fromIndex, int toIndex, int k) {
if (toIndex - fromIndex <= 1) {
return;
}
// Sort NaN / count signed zeros.
// Caution: This loop contributes significantly to the runtime.
int cn = 0;
int end = toIndex;
for (int i = toIndex; --i >= fromIndex;) {
final double v = a[i];
// Count negative zeros using a sign bit check
if (Double.doubleToRawLongBits(v) == Long.MIN_VALUE) {
cn++;
// Change to positive zero.
// Data must be repaired after selection.
a[i] = 0.0;
} else if (v != v) {
// Move NaN to end
a[i] = a[--end];
a[end] = v;
}
}
// Partition
if (end - fromIndex > 1 && k < end) {
QuickSelect.select(a, fromIndex, end - 1, k);
}
// Restore signed zeros
if (cn != 0) {
// Use partition index below zero to fast-forward to zero as much as possible
for (int j = a[k] < 0 ? k : -1;;) {
if (a[++j] == 0) {
a[j] = -0.0;
if (--cn == 0) {
break;
}
}
}
}
}
/**
* Partition the array such that indices {@code k} correspond to their correctly
* sorted value in the equivalent fully sorted array.
*
* <p>This method pre/post-processes the data and indices to respect the ordering
* imposed by {@link Double#compare(double, double)}.
*
* @param fromIndex Index of the first element (inclusive).
* @param toIndex Index of the last element (exclusive).
* @param a Values.
* @param k Indices (may be destructively modified).
*/
private static void doSelect(double[] a, int fromIndex, int toIndex, int[] k) {
if (k.length == 0 || toIndex - fromIndex <= 1) {
return;
}
// Sort NaN / count signed zeros.
// Caution: This loop contributes significantly to the runtime for single indices.
int cn = 0;
int end = toIndex;
for (int i = toIndex; --i >= fromIndex;) {
final double v = a[i];
// Count negative zeros using a sign bit check
if (Double.doubleToRawLongBits(v) == Long.MIN_VALUE) {
cn++;
// Change to positive zero.
// Data must be repaired after selection.
a[i] = 0.0;
} else if (v != v) {
// Move NaN to end
a[i] = a[--end];
a[end] = v;
}
}
// Partition
int n = 0;
if (end - fromIndex > 1) {
n = k.length;
// Filter indices invalidated by NaN check
if (end < toIndex) {
for (int i = n; --i >= 0;) {
final int index = k[i];
if (index >= end) {
// Move to end
k[i] = k[--n];
k[n] = index;
}
}
}
// Return n, the count of used indices in k.
// Use this to post-process zeros.
n = QuickSelect.select(a, fromIndex, end - 1, k, n);
}
// Restore signed zeros
if (cn != 0) {
// Use partition indices below zero to fast-forward to zero as much as possible
int j = -1;
if (n < 0) {
// Binary search on -n sorted indices: hi = (-n) - 1
int lo = 0;
int hi = ~n;
while (lo <= hi) {
final int mid = (lo + hi) >>> 1;
if (a[k[mid]] < 0) {
j = mid;
lo = mid + 1;
} else {
hi = mid - 1;
}
}
} else {
// Unsorted, process all indices
for (int i = n; --i >= 0;) {
if (a[k[i]] < 0) {
j = k[i];
}
}
}
for (;;) {
if (a[++j] == 0) {
a[j] = -0.0;
if (--cn == 0) {
break;
}
}
}
}
}
/**
* Partition the array such that index {@code k} corresponds to its correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param k Index.
* @throws IndexOutOfBoundsException if index {@code k} is not within the
* sub-range {@code [0, a.length)}
*/
public static void select(int[] a, int k) {
IndexSupport.checkIndex(0, a.length, k);
if (a.length <= 1) {
return;
}
QuickSelect.select(a, 0, a.length - 1, k);
}
/**
* Partition the array such that indices {@code k} correspond to their correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param k Indices (may be destructively modified).
* @throws IndexOutOfBoundsException if any index {@code k} is not within the
* sub-range {@code [0, a.length)}
*/
public static void select(int[] a, int[] k) {
IndexSupport.checkIndices(0, a.length, k);
if (k.length == 0 || a.length <= 1) {
return;
}
QuickSelect.select(a, 0, a.length - 1, k, k.length);
}
/**
* Partition the array such that index {@code k} corresponds to its correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param fromIndex Index of the first element (inclusive).
* @param toIndex Index of the last element (exclusive).
* @param k Index.
* @throws IndexOutOfBoundsException if the sub-range {@code [fromIndex, toIndex)} is out of
* bounds of range {@code [0, a.length)}; or if index {@code k} is not within the
* sub-range {@code [fromIndex, toIndex)}
*/
public static void select(int[] a, int fromIndex, int toIndex, int k) {
IndexSupport.checkFromToIndex(fromIndex, toIndex, a.length);
IndexSupport.checkIndex(fromIndex, toIndex, k);
if (toIndex - fromIndex <= 1) {
return;
}
QuickSelect.select(a, fromIndex, toIndex - 1, k);
}
/**
* Partition the array such that indices {@code k} correspond to their correctly
* sorted value in the equivalent fully sorted array.
*
* @param a Values.
* @param fromIndex Index of the first element (inclusive).
* @param toIndex Index of the last element (exclusive).
* @param k Indices (may be destructively modified).
* @throws IndexOutOfBoundsException if the sub-range {@code [fromIndex, toIndex)} is out of
* bounds of range {@code [0, a.length)}; or if any index {@code k} is not within the
* sub-range {@code [fromIndex, toIndex)}
*/
public static void select(int[] a, int fromIndex, int toIndex, int[] k) {
IndexSupport.checkFromToIndex(fromIndex, toIndex, a.length);
IndexSupport.checkIndices(fromIndex, toIndex, k);
if (k.length == 0 || toIndex - fromIndex <= 1) {
return;
}
QuickSelect.select(a, fromIndex, toIndex - 1, k, k.length);
}
}