1-based to 0-based array in a function

[ferrad]

Numerical Recipes has a heap sort function below. It sorts an array [1..n], however C arrays start at 0, and I would like to modify this function to deal with zero-based arrays.

void hpsort(long n, double ra[]) {
// Sorts an array ra[1..n] into ascending numerical order using the Heapsort algorithm. n is
// input; ra is replaced on output by its sorted rearrangement.

    unsigned long i,ir,j,l;
    double rra;
    if (n < 2) return;
    l = (n >> 1)+1;
    ir = n;
    for (;;) {
        if (l > 1) {
            rra = ra[--l];
        } else {
            rra = ra[ir];
            ra[ir] = ra[1];
            if (--ir == 1) {
                ra[1] = rra;
                break;
            }
        }
        i = l;
        j = l+l;
        while (j <= ir) {
            if (j < ir && ra[j] < ra[j+1]) j++;
            if (rra < ra[j]) {
                ra[i] = ra[j];
                i = j;
                j <<= 1;
            } else break;
        }
        ra[i] = rra;
    }
}

I tried modifying the first bit to:

        if (l > 0) {
            rra = ra[--l];
        } else {
            rra = ra[ir];
            ra[ir] = ra[0];
            if (--ir == 0) {
                ra[0] = rra;
                break;
            }
        }

However it just goes in an endless loop in the while loop. Can anyone spot the problem?

[tipaye]

rra = ra[--l];

here, you're actually modifying l, so the value of l will always be 0 after that line of code; hence infinite loop.

Try rra = ra[l-1]; instead.

I have to say I haven't checked what you're doing, just tried to spot the source of your infinite loop.

[JLBorges]

#include <iostream>

// view of an array with the element positions offset by BASE
template < typename T, std::size_t BASE = 1 > struct array_view {

    std::size_t size() const { return sz ; }

    T& operator[] ( std::size_t pos ) { return ptr[ pos - BASE ] ; } // invariant pos < ( size() + BASE )
    // const T& operator[] ( std::size_t pos ) const { return ptr[ pos - BASE ] ; }

    T* ptr ;
    std::size_t sz ;
};

template < typename T > array_view<T,1> make_one_based( T* ptr, std::size_t n ) { return { ptr, n } ; }
template < typename T, std::size_t N > array_view<T,1> make_one_based( T (&arr)[N] ) { return { arr, N } ; }


void hpsort( /*long n, double ra[] */ array_view<double,1> ra ) { // *** modifioed

    const long n = ra.size() ; // *** added

    // *** rest of the function is retained verbatim

    // Sorts an array ra[1..n] into ascending numerical order using the Heapsort algorithm. n is
    // input; ra is replaced on output by its sorted rearrangement.

    unsigned long i,ir,j,l;
    double rra;
    if (n < 2) return;
    l = (n >> 1)+1;
    ir = n;
    for (;;) {
        if (l > 1) {
            rra = ra[--l];
        } else {
            rra = ra[ir];
            ra[ir] = ra[1];
            if (--ir == 1) {
                ra[1] = rra;
                break;
            }
        }
        i = l;
        j = l+l;
        while (j <= ir) {
            if (j < ir && ra[j] < ra[j+1]) j++;
            if (rra < ra[j]) {
                ra[i] = ra[j];
                i = j;
                j <<= 1;
            } else break;
        }
        ra[i] = rra;
    }
}

int main()
{
    double a[] = { 34, 43, 98, 56, 57, 68, 97, 14, 98, 74, 64, 63, 55, 35, 24, 25, 45, 67 } ;
    for( double v : a ) std::cout << v << ' ' ;
    std::cout << '\n' ;

    auto array_view = make_one_based(a) ;
    hpsort( array_view ) ;
    for( double v : a ) std::cout << v << ' ' ;
    std::cout << '\n' ;
}

http://coliru.stacked-crooked.com/a/8fd2c4d06c8d1aba

[dhayden]

At line 7, insert:
--ra;
and include a long comment about how the algorithm is from Numerical Recipes and uses 1-based arrays.

[ferrad]

Thanks dhayden, works like a charm.

[JLBorges]

> At line 7, insert: --ra;

This will engender undefined behaviour if ra points to the first element of an array.

[ferrad]

ra does initially point to the first element of an array. --ra just makes it point to the position one below, but this location will never be accessed as the code is designed to access element 1, not 0.

[JLBorges]

> ra does initially point to the first element of an array. --ra just makes it point to the position one below

--ra engenders undefined behaviour.


> but this location will never be accessed

It does not matter what is done, or not done later; even if ra is never used again, the entire program has undefined behaviour.

The gory details:

the expression --x is exactly equivalent to x -= 1, ...
All arithmetic conversion rules and pointer arithmetic rules defined for arithmetic operators apply.

http://en.cppreference.com/w/cpp/language/operator_incdec

The behavior of every builtin compound-assignment expression E1 op= E2 ... is exactly the same as the behavior of the expression E1 = E1 op E2, except that the expression E1 is evaluated only once and that it behaves as a single operation with respect to indeterminately-sequenced function calls

http://en.cppreference.com/w/cpp/language/operator_assignment

If the pointer P points to the ith element of an array, then the expressions P+n, n+P, and P-n are pointers of the same type that point to the i+nth, i+nth, and i-nth element of the same array, respectively. ...
Any other situations (that is, attempts to generate a pointer that isn't pointing at an element of the same array or one past the end) invoke undefined behavior.

http://en.cppreference.com/w/cpp/language/operator_arithmetic

[ferrad]

Sorry I'm totally lost in your response. It all works perfectly. --ra just lets ra point at the location just before it. Which is what I want. Thanks.

[dhayden]

JLBorges's point is that, technically, when ra initially points to the first element, the behavior of --ra is undefined. There are plenty of reasons for this: it might point outside the address space, or it might interfere with optimizations. The point is that compiler writers don't have to ensure that it works.

In practice, it might work just fine, but there are no guarantees.

[ferrad]

One of the many reasons I consider C/C++ as very dangerous odd languages. Nothing is as it seems.