[comp.lang.c] Numerical Recipes is non-portable code

rob@kaa.eng.ohio-state.edu (Rob Carriere) (09/02/88)

In some article in a chain so long ``rn'' chokes on it, gwyn@smoke.ARPA 
(Doug Gwyn) writes:
> [ that I write that ]:
>>Trivial refutation time!  Surely it is obvious that ``numerical
>>problems'' forms a (large) superset of ``array/matrix code'' as far as
>>numerical analysis is concerned?
>Trivial indeed!  If the code does not involve arrays/matrices,
>the issue of 0-based or 1-based indexing doesn't even arise.

Well... I was going to apologize for the tone of my last posting, but
you're sure not making it any easier.  Guess I still will though, I
short-changed you yesterday, so I'll overpay you today.  So, please
consider that my last posting was affected by 3 days high on emotions
and low on sleep, making my temper about a mile shorter than it should
have been, and please accept my apologies for the tone of that
posting.  Peace now?

As for your refutation of my refutation, it should have been clear
from the following paragraph in the posting that I meant code other
than the generic invert-matrix, do-svd etc stuff.  Most of those can
indeed be written base 1 or base 0 without any difference.  There are
however problems (of a more specialized nature, and the arrays do not
necessarily model vectors or matrices) where the base does make a
difference.  I should know, I write the stuff!

Several people have suggested that for these cases one either carries
the address arithmatic through the code, or accesses the arrays
through a function that does the arithmatic for you.  The first
approach has the disadvantage that even using special lay-out tricks
the readability of the resulting code drops to just about zero (again,
this is experience, not just unwillingness to try).  Both approaches
have the disadvantage that they introduce inefficiencies into the
code.  Since I am talking about things that eventually are going to be
somewhere in a real-time system, efficiency is worth more than ``gee,
it'd be neat if it were a little faster''; it is access this array
5000 times per second or the voice link dies/ the plant blows up/ the
missile flies the other way/ fill in your favorite disaster here.

Rob Carriere
It's long, but I hope it clears the dust a little. *



*
Yup, it's a broom!

gwyn@smoke.ARPA (Doug Gwyn ) (09/02/88)

In article <557@accelerator.eng.ohio-state.edu> rob@kaa.eng.ohio-state.edu (Rob Carriere) writes:
>Peace now?

Sure; why not?

>Several people have suggested that for these cases one either carries
>the address arithmatic through the code, or accesses the arrays
>through a function that does the arithmatic for you.

Or you can allocate an extra element and go ahead and use 1-based
indexing.  Nobody is forcing you to use the [0] element of an array.
It is seldom that the small amount of wasted space would matter.

news@ism780c.isc.com (News system) (09/03/88)

In article <557@accelerator.eng.ohio-state.edu> rob@kaa.eng.ohio-state.edu (Rob Carriere) writes:
>                                                            There are
>however problems (of a more specialized nature, and the arrays do not
>necessarily model vectors or matrices) where the base does make a
>difference.  I should know, I write the stuff!
>
>
>Rob Carriere

Could you give a small example, perhaps I don't understand. with a
declaration (using an extended version of C):

     int x(5:49),y(0:44);
     int i,k;
     i=6;
     k=x(i);        /* easy to write access to second item */
     k=y(i-5);      /* more trouble to write access to second item */

However, with a reasonable compiler, both assignments to k will generate
identical code.  Your compiler may vary :-).  But if you are writing time
critical code with a faulty compiler, I guess you do what you can.

     Marv Rubinstein

rob@kaa.eng.ohio-state.edu (Rob Carriere) (09/05/88)

In article <15029@ism780c.isc.com> marv@ism780.UUCP (Marvin Rubenstein) writes:
>In article <557@accelerator.eng.ohio-state.edu> rob@kaa.eng.ohio-state.edu 
>(Rob Carriere) writes:
>> [ still on non-0 base arrays ]
>Could you give a small example, perhaps I don't understand. with a
>declaration (using an extended version of C):
>
>     int x(5:49),y(0:44);
>     int i,k;
>     i=6;
>     k=x(i);        /* easy to write access to second item */
>     k=y(i-5);      /* more trouble to write access to second item */
>
>However, with a reasonable compiler, both assignments to k will generate
>identical code.  Your compiler may vary :-).  But if you are writing time
>critical code with a faulty compiler, I guess you do what you can.

You more or less have the idea (except that the array would be float,
I'm one of them there dirty numerical geeks :-).  Actually, a
reasonable (as opposed to a *very* good) compiler would probably
generate better code for k=x[i]; because it cannot do the optimization
I'm trying to manually make either: illegal pointer!  To add insult to
the poor compiler's injury, in most cases the code would actually
read:

void foo( int n )
{
    float *x;
    int i;
    float k;

    i = 6;
    x = vector( 5, n );
...
    k = x[i];
...
}
...
foo( 49 );

As a further and at least as important point, the code is generated
from some underlying math that insists that x starts at (say) 5.  So I
can fake it and confuse everyone with my constant -5 offsets in all
the indeces, or I can use a non-zero base offset.

If you are curious as to what perverse math would generate such stuff
(and it's *not* because they don't know to start counting at zero!), I
sent off a mail message with a couple of examples to Marty Fouts who
asked the same question.  It's a bit too long and too marginal to this
group to post, but I'd gladly send you a copy (be sure to give me a
*very* complete email address, our mail deamon spends far more time
earning its last than its first name!)

Rob Carriere

karl@haddock.ima.isc.com (Karl Heuer) (09/08/88)

rob@kaa.eng.ohio-state.edu (Rob Carriere) writes:
>marv@ism780.UUCP (Marvin Rubenstein) writes:
>>[In a hypothetical extended C, with arrays x[5:49] and y[0:44], and an
>>integer variable i, the accesses `x[i]' and `y[i-5]' should generate the
>>same code on a reasonable compiler.]
>
>Actually, a reasonable compiler would probably generate better code for
>k=x[i]; because it cannot do the optimization I'm trying to manually make
>either: illegal pointer!

If the architecture in question is "normal", then the compiler *is* free to do
that optimization.  If the hardware will fault on the mere mention of an
illegal pointer, then it can't optimize it -- but how do you think it's going
to implement z=x[i], then?  They'll still generate the same code, but in this
case it'll be equally bad rather than equally good.

Karl W. Z. Heuer (ima!haddock!karl or karl@haddock.isc.com), The Walking Lint

rob@raksha.eng.ohio-state.edu (Rob Carriere) (09/08/88)

In article <7142@haddock.ima.isc.com> karl@haddock.ima.isc.com (Karl Heuer) 
writes:
>rob@kaa.eng.ohio-state.edu (Rob Carriere) writes:
>>marv@ism780.UUCP (Marvin Rubenstein) writes:
>>>[In a hypothetical extended C, with arrays x[5:49] and y[0:44], and an
>>>integer variable i, the accesses `x[i]' and `y[i-5]' should generate the
>>>same code on a reasonable compiler.]
>>
>> [ disagree ]
>If the architecture in question is "normal", then the compiler *is* free to do
>that optimization.  If the hardware will fault on the mere mention of an
>illegal pointer, then it can't optimize it -- but how do you think it's going
>to implement z=x[i], then?  They'll still generate the same code, but in this
>case it'll be equally bad rather than equally good.

Actually, if the system has memory to burn, the compiler might
allocate the 5 extra int's (probably only doing so for small lower
bounds, but those do make up almost all the cases anyway).  This is
similar to the fix to the Numerical Recipes vector() routine that I
posted a while ago.

Rob Carriere
"There's nothing quite like a broken architecture"