Что не так с этой реализацией интерполяционного поиска?

Question

Это распространенная в C / C ++ реализация алгоритма интерполяционного поиска, встречающаяся в Интернете. Однако при использовании отсортированного массива из примерно 100000 целых чисел средняя переменная начинает генерировать отрицательные индексы массива, вызывая ошибку сегментации. В чем может быть проблема?

#include <stdlib.h>
#include <stdio.h>
#include <time.h>
int interpolationSearch(int sortedArray[], int toFind, int len) {
    // Returns index of toFind in sortedArray, or -1 if not found
    int low = 0;
    int high = len - 1;
    int mid;

    while (sortedArray[low] <= toFind && sortedArray[high] >= toFind) {
        mid = low + ((toFind - sortedArray[low]) * (high - low)) /
              (sortedArray[high] - sortedArray[low]);

        if (sortedArray[mid] < toFind) {
            low = mid + 1;
        } else if (sortedArray[mid] > toFind) {
            high = mid - 1;
        } else {
            return mid;
        }
    }

    if (sortedArray[low] == toFind)
        return low;
    else
        return -1; // Not found
}

int main(void) {
    srand(time(0));
    int arr[100000];
    for (int i=0; i<100000; i++) {
        arr[i] = rand()%100000;
    }

    int length = sizeof(arr)/sizeof(int);
    qsort(arr,length,sizeof(int),order);

    for (int j=0; j<10000; j++) {
        interpolationSearch(arr,rand()%100000,length);
    }
}

Answer 1

Подвыражение: ((toFind - sortedArray[low]) * (high - low))

... может легко вычислить что-то вроде: ((99999-0) * (99999-0)) == 99999^2

... что намного больше, чем 2 ^ 31 (==диапазон 32-разрядных целых чисел со знаком).

Как только оно превысит 2 ^ 31-1, целое число переполнится в отрицательные числа, отсюда и ваши отрицательные индексы.Если оно превысит 2 ^ 32 (что также возможно), то (скорее всего, технически не определено) вы потеряете старшие биты и в итоге получите фактически случайные смещения, как положительные, так и отрицательные.

Чтобы избежать всего этого, вам нужно тщательно выполнить математику, чтобы убедиться, что ни одно из ваших подвыражений не дает целочисленного переполнения.Обычно самый простой способ сделать это - преобразовать в число с плавающей запятой, диапазон которого на много порядков больше, чем 32-разрядные целые числа.

В конечном итоге такая интерполяция для двоичного поиска обычно не стоит - затраты на вычисление интерполяции обычно больше, чем несколько дополнительных итераций цикла, которые он «сохраняет».

Answer 2

Как объяснили другие ответы, вы пытаетесь вычислить выражение вида

A * B / C

но это не так, потому что A * B переполнен. Предложение пересмотреть выражение до

A * (B / C)

не будет работать, потому что обычно B меньше C, поэтому целочисленное деление усекается до нуля.

Предложение переключиться на число с плавающей точкой будет работать, но будет дорогостоящим. Но вы можете использовать фиксированную точку , преобразовав выражение в:

A * ((B * F) / C) / F

(где F - тщательно выбранная степень двух).

Answer 3

Проблема в выражении, которое вычисляет mid. Продукт может легко переполниться даже с 32-битными целыми числами. Тогда это становится отрицательным. Вероятно, было бы лучше выполнить разделение перед продуктом.

Изменение средних вычислений для использования 64-битных целых чисел (по крайней мере, для промежуточных вычислений) исправляет проблемы.

Ниже моя модифицированная версия (int64_t определено в <stdint.h>:

int interpolationSearch(int sortedArray[], int toFind, int len) {
    // Returns index of toFind in sortedArray, or -1 if not found
    int low = 0;
    int high = len - 1;
    int mid;

    int l = sortedArray[low];
    int h = sortedArray[high];

    while (l <= toFind && h >= toFind) {
        int64_t high_low = (high - low);
        int64_t toFind_l = (toFind - l);
        int64_t product = high_low*toFind_l;
        int64_t h_l = h-l;
        int64_t step = product / h_l;
        mid = low + step;

/*        mid = (low + high)/2;*/
        int m = sortedArray[mid];

        if (m < toFind) {
            l = sortedArray[low = mid + 1];
        } else if (m > toFind) {
            h = sortedArray[high = mid - 1];
        } else {
            return mid;
        }
    }

    if (sortedArray[low] == toFind)
        return low;
    else
        return -1; // Not found
}

Еще более простым решением было бы сделать дихотомический поиск вместо интерполяции, просто используя: mid = (low + high) / 2. даже если он сходится немного медленнее, чем интерполяция, он избегает нескольких операций, включая произведение и деление, тем самым делая внутренний цикл быстрее. Не уверен, что потенциальная более быстрая сходимость интерполяции компенсирует эту потерю простоты.

Я провел несколько тестов производительности. Источник моей тестовой программы включен в этот вопрос

Удивительно (для меня) использование float дает более эффективную программу, чем использование больших целых чисел. В моей системе бинарный поиск стал быстрее примерно для 1000 элементов в массиве. Для массивов размером 100000 интерполяционный поиск почти в два раза быстрее, чем простой двоичный поиск.