C vs. FORTRAN

[Richard A. O'Keefe]

In article <20480 at beta.lanl.gov> jlg at beta.lanl.gov (Jim Giles) writes:
>In article <143 at quintus.UUCP>, ok at quintus.uucp (Richard A. O'Keefe) writes:
>> Not to knock either C or Fortran, but this may be a case where C is better
>> to use than Fortran.

>The exponentiation operator isn't bad just because some people use it badly.

I didn't say that it was bad.  What I do say is that it is _MISLEADING_.
In C, for example, most of the operations for which there is special
syntax correspond to machine operations of similar cost.

What with prototypes and all, an ANSI C compiler is fully entitled to treat
	#include <math.h>
	double x;
	...
	x = pow(x, 3);
	x = pow(x, 3.5);
just the same as Fortran 77 would treat
	DOUBLE PRECISION X
	...
	X = X**3
	X = X**3.5
It is already the case on the UNIX systems where I've tried it that
	pow(x, (double)3)
uses the iterative method used for X**3.  (I personally don't think that
is a good thing, but that's another story.)

>X**3 is more readable than X*X*X (particularly is X is an expression).
pow(x,3) isn't _that_ hard to read, either.
>X**3.5 is more readable than EXP(3.5*LOG(X)).
True, but again, pow(x,3.5) isn't terribly hard to read, and in any
case the two expressions X**3.5 and EXP(3.5*LOG(X)) don't mean exactly
the same thing (may have different over/underflow conditions, yield
different results, &c).

>If you always do pow(X,Y), the compiler has no choices to make.
Misleading and soon to be wrong.  Misleading, because the run-time library
_does_ have some choices to make, and current libraries make them.  Soon
to be wrong, because ANSI C compilers are allowed to detect the standard
library functions and do special things with them.

>If the potential for misuse doesn't even result in incorrect code
>(as in this case), then there is no reason not to include the feature.

Well, ANSI C _does_ include the feature.  My point was simply that the
notation "pow(x, y)" _looks_ as though it might be costly, while the
notation "x ** y" _looks_ as though it is cheap and simple, and that
the first perception is the correct one, so that ANSI C's notation may
be better.  In fact it is _not_ the case in this example that the
misuse did not result in incorrect code.  Consider the two fragments:

C	Assume that N.ge.2		/* Assume that N >= 2 */
C	Fragment 1			/* Fragment 1 */
	T = A(1)			t = a[0];
	DO 10 I = 2,N			for (i = 1; i < n; i++) {
	   T = T + A(I) * X**(I-1)	    t += a[i]*pow(x, (double)i);
10	CONTINUE			}
C	Fragment 2			/* Fragment 2 */
	T = A(N)			t = a[N-1];
	DO 20 I = N-1,1,-1		for (i = n-2; i >= 0l i--) {
	    T = T*X + A(I)		    t = t*x + a[i];
20	CONTINUE			}

The two fragments are *not* equivalent.  It is easy to come up with an
example where fragment 2 correctly yields 1.0 and fragment 1 overflows.
{I am not asserting that fragment 2 is always the right thing to use!}

A good numeric programmer will use X**N __very__ seldom for other than
constant arguments.

What C _does_ lack that a numeric programmer might miss is a single-
precision version of pow().  Given
	float x, y, u, v;
	x = pow(u,v)/y;
it is not in general safe to evaluate pow(u,v) in single precision.
How much faster is float**float than double**double?
Sun-3/50 (-fsoft)	1.2
Sun-3/160 (-f68881)	2.5
Sequent (80387)		1.4
If this is typical (and I have no reason to expect that it is), the
lack of single precision elementary functions in C comes into the
"nuisance" category rather than the "major problem" one.  Since UNIX
systems often come with Fortran these days, and Fortran has single
precision elementary functions, it seems a pity not to let ANSI C
share them.