indexx() Числовые рецепты (C) Алгоритм сортировки индекса странным образом игнорирует первые два элемента

Question

Я пытаюсь использовать алгоритм indexx() из Numeric Recipes (NR) в C и обнаружил очень странную ошибку.

(NR общедоступен здесь: http://www2.units.it/ipl/students_area/imm2/files/Numerical_Recipes.pdf стр. 338, раздел 8.4)

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

Ниже приведен минимальный рабочий пример, показывающий, что алгоритм игнорирует первые два элемента. В выходном массиве первые два элемента всегда равны 0, за которым следует длина массива (в данном примере 9). Остальные элементы отображаются правильно отсортированными.

О, и я пытался спросить на форумах NR, но долго ждал активации моей учетной записи... Большое спасибо заранее!

[Отредактированные имена массивов]

#include "nr_c/nr.h"

#include <stdio.h>
#include <stdlib.h>


int main()
{
    float unsorted[9] = {4., 5., 2., 6., 3., 8., 1., 9., 7.};
    long unsigned int sort[9];

    printf("Unsorted input array:\n");
    for (int i=0; i<9; i++) {
        printf("%.1f  ", unsorted[i]);
    }
    printf("\n\n");

    indexx(9, unsorted, sort);

    printf("Indexx output array:\n");
    for (int i=0; i<9; i++) {
        printf("%d    ", sort[i]);
    }
    printf("\n\n");

    printf("Should-be-sorted array:\n");
    for (int i=0; i<9; i++) {
        printf("%.1f  ", unsorted[sort[i]]);
    }
    printf("\n\n");

    return 0;
}

Выход:

Unsorted input array:
4.0  5.0  2.0  6.0  3.0  8.0  1.0  9.0  7.0  

Indexx output array:
0    9    6    2    4    1    3    8    5    

Should-be-sorted array:
4.0  0.0  1.0  2.0  3.0  5.0  6.0  7.0  8.0 

Answer 1

В коде "Числовые рецепты" используется индексирование на основе 1. (Из-за своего происхождения на Фортране, первая версия была написана на Фортране, и Фортран использует индексирование на основе 1 для массивов и матриц. Версия C, вероятно, была основана на выводе f2c...) В исходном коде в вопросе indexx() Функция игнорирует первый элемент как unsorted[] а также sort[] массивы. Плюс: он обращается к одному элементу за последними элементами массива (что приводит к UB). В результате в индексном массиве присутствуют как 0, так и 9. (начальный 0 на самом деле неинициализированная память, но она не была затронута indexx() функция)

---

Если моя гипотеза верна, должно работать следующее:

---

#include "nr_c/nr.h"

#include <stdio.h>
#include <stdlib.h>


int main()
{
    float unsorted[9] = {4., 5., 2., 6., 3., 8., 1., 9., 7.};
    long unsigned int sort[9];

    printf("Unsorted input array:\n");
    for (int i=0; i<9; i++) {
        printf("%.1f  ", unsorted[i]);
    }
    printf("\n\n");

    indexx(9, unsorted-1, sort-1); // <<-- HERE

    printf("Indexx output array:\n");
    for (int i=0; i<9; i++) {
        printf("%d    ", sort[i]);
    }
    printf("\n\n");

    printf("Should-be-sorted array:\n");
    for (int i=0; i<9; i++) {
        printf("%.1f  ", unsorted[sort[i]-1]); // <<-- AND HERE
    }
    printf("\n\n");

    return 0;
}

---

numerical recipes in C код предполагает индексирование на основе 1: массив размером N имеет индексы 1..N, Это было сделано, чтобы не перепутать математиков. (в результате, целое поколение программистов было сбито с толку) Функции размещения используются как обертки malloc(), которые обманывают, возвращая указатель на -1член выделенного пространства. Для float вектор это может быть как:

---

#include <stdlib.h>

float * float_vector(unsigned size)
{
float * array;
array = calloc( size, sizeof *array);
if (!array) return NULL;
return array -1;
}