indexx() 数值食谱 (C) 奇怪地忽略前两个元素的索引排序算法

Question

我正在尝试在 C 中使用 Numerical Recipes (NR) 中的 indexx() 算法，但发现了一个非常奇怪的错误。

（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

数值食谱 C 代码使用从 1 开始的索引。（因为它的 FORTRAN 起源，第一个版本是用 FORTRAN 编写的，fortran 对数组和矩阵使用基于 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() 的包装器，它通过返回指向已分配 space 的第 -1 个成员的指针来作弊.对于 float 向量，这可能是这样的：

#include <stdlib.h>

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

indexx() 数值食谱 (C) 奇怪地忽略前两个元素的索引排序算法

indexx() Numerical Recipes (C) index sorting algorithm weirdly ignoring first two elements

c

sorting

indexing

numerical-recipes