Should I convert FORTRAN code to C? (Really, multidim. arrays in C)

[Peter Shenkin]

In article<10032 at ames.arc.nasa.gov> eugene at pioneer.UUCP (Eugene N. Miya) writes:
>
>I have been asked to post this.  Please send mail to him, not me.
[ note:  "him" is:
    Matt Wette
    mwette%gauss at hub.ucsb.edu
]
>                                                     ...I see two
>issues in making transition to something better.  The first is the
>choosing of a suitable replacement language.  The second is that of
>translating Fortran Dusty Decks.
>
>1.  Choice of Language
>
>                           Multidimensioned arrays in C are hard to
>implement well and in a way that will be portable.  One person's idea
>of how to implement multidimesioned arrays will differ from another's.
>Now try to write code that will handle both!  YECH!

Just on this issue, I wonder what other numerikers in the audience think
of Press et al's solution to this problem, given in "Numerical Recipes in
C", Wiley, 1988.  To make a long story short, they supply a routine:

    float **matrix( int nrl, int nrh, int ncl, int nch )

which uses dope vectors to dynamically allocate a 2-dimensional array of floats
with range [nrl...nrh][ncl...nch].  They also supply similar routines
for double-precision and integer matrices (called dmatrix() and imatrix() ),
and the corresponding "free" routines.  They do a similar thing for
vectors, and obviously an extension of this approach using varargs
-- or whatever the hell it's called in the draft ANSI standard :-) --
could easily be accomplished for arbitrary dimensionality.

This is, in Matt's words, just "one person's idea of how to implement
multidimensioned arrays", but it seems like a reasonable one.  The reason
for raising this issue here on the net is simply that if we who do write
numerical code in C (for better or for worse!) were to standardize on
this practice, we could better understand each others' code (or at least
this part of it!)

For the record, this generalizes to the following:

	float **array( int Ndim, 
	  int n1l, int n1h, int n2l, int n2h, ..., int nNl, int nNh )
	  /* allocate an N-dimensional array of floats with range
	   *   [n1l, n1h][n2l,n2h]...[nNl,nNh]
	   */
	void free_array( float **ar, int Ndim,
	  int n1l, int n1h, int n2l, int n2h, ..., int nNl, int nNh )
	  /* free array allocated with array() */

Then we also would have the following routines:
	int **iarray()     /* allocate an N-dimensional array of ints....  */
	double **darray()  /* allocate an N-dimensional array of doubles....  */
and the corresponding free routines.

Of course, we could then have:

#define vector( NL,NH )              array( 1, (NL), (NH) )
#define matrix( NRL,NRH,NCL,NCH )    array( 2, (NRL), (NRH), (NCL), (NCH) )

and perhaps even (for "real" C programmers):

#define vector0( N )              array( 1, 0, (N) )
#define matrix0( NR,NC )          array( 2, 0, (NR), 0, (NC) )

...and so on.  All presumably in some header file with a name like "arrays.h".

Any thought on this?

("Numerical analysts of the world, unite....")
-- 
*******************************************************************************
Peter S. Shenkin,    Department of Biological Sciences,    Columbia University,
New York, NY   10027         Tel: (212) 280-5517 (work);  (212) 829-5363 (home)
shenkin at cubsun.bio.columbia.edu    shenkin%cubsun.bio.columbia.edu at cuvma.BITNET