Conformant Arrays in C

[David Hough]

I have been wondering about a problem that's been discussed,
often incorrectly, in comp.lang.c.  The problem is easy enough
to understand if you've ever written linear algebra software,
but the correct "in the spirit of C" solution wasn't clear
to me.

Fortunately Richard O'Keefe has volunteered some ideas.
I'd appreciate feedback on them.  If you're simply opposed
to making C suitable for scientific programming, the net has
already heard from you, so you needn't repeat your views,
but if you have ideas for better syntax or better ways to
achieve the goal, please speak up.  The following is one of my
comments on the January 1988 X3J11 draft for ANSI C,
followed by the proposed solution.


Comment #24, Section 3.5.4.2:   fix arrays

     I know no C translation that's as clear as the  follow-
ing Fortran code:

        SUBROUTINE MATMUL(X,LX,Y,LY,Z,LZ,NX,NY,NZ)
        REAL X(LX,*),Y(LY,*),Z(LZ,*)
        DO 1 I=1,NX
                DO 2 J=1,NZ
                SUM=0
                        DO 3 K=1,NY
 3                      SUM=SUM+X(I,K)*Y(K,J)
 2              Z(I,J)=SUM
 1      CONTINUE
        END

Code like this is at the heart of most of the major portable
Fortran  libraries  of  mathematical software developed over
the last twenty years. The declared leading dimensions of X,
Y, and Z are not known until runtime.

     The Draft, like traditional C, disallows the equivalent

        void matmul(x,lx,cx,...)
        int lx, cx;
        double x[lx][cx] ;

unless lx and cx are  known  at  compile  time.   Equivalent
functionality  can  be  obtained  by  treating all arrays as
one-dimensional and doing the  subscripting  "by  hand",  so
that  it's harder to get right and harder to optimize.  Fix-
ing C's handling of arrays would be a far worthier task  for
the  X3J11 experts than (for instance) specifying details of
the math library.  Maybe there's some good fundamental  rea-
son  why  variables  are disallowed as array dimensions; C's
treatment of arrays has always been confusing.

     Note that the section numberings in the Draft  and  Ra-
tionale are out of synch in the declarations subchapter.

     Recommendation: Fix, by language design  if  necessary,
the  treatment  of  arrays  in  C  to  properly  accommodate
multiply-dimensioned arrays whose  dimensions  vary  at  run
time.   The  goal is not to duplicate Fortran's features ex-
actly, but rather to insure that portable linear algebra li-
braries are as easy to create in C as in Fortran.

Array proposal
----- --------

     Richard O'Keefe has kindly provided the following  out-
line of how an array facility might work in C:

     The  smallest  possible  change  would  be  "conformant
arrays" much as in ISO Pascal or in Turing, where the param-
eter declaration was something like

        void matmul(double a[p][q], b[q][r], c[p][r])
            {
                int i, j, k;
                double t;
                for (i = p; --i >= 0; )
                    for (j = r; --j >= 0; ) {
                        t = 0.0;
                        for (k = q; --k >= 0; ) t *= a[i][j]*b[j][k];
                        c[i][j] = t;
                    }
            }

If p, q, and r are #defined to be constant expressions, this
is already legal, so we need one more thing to indicate that
p,q,r are being defined here, not used.  Consider  the  fol-
lowing:

        declarator:     ...
                  |     declarator '[' subscript_spec ']'
                  ;

        subscript_spec: '?' identifier
                      | constant_expression
                      | /* empty */
                      ;

where  /*  empty  */  is   only   allowed   as   the   first
subscript_spec,  and  ?id  is  only  allowed  in  a function
header.  To avoid having to specify what happens if  ?x  ap-
pears several times but the arguments don't agree with that,
make it illegal, so the first suggestion would  have  to  be
written

        void matmul(double a[?ar][?ac], /* ar >= p, ac >= q */
                    double b[?br][?bc], /* br >= q, bc >= r */
                    double c[?cr][?cc], /* cr >= p, cc >= r */
                    int p,              /* 0 <= p <= min(ar,cr) */
                    int q,              /* 0 <= q <= min(ac,br) */
                    int r)              /* 0 <= r <= min(bc,cc) */
            {
                int i, j, k;
                double t;
                for (i = p; --i >= 0; )
                    for (j = r; --j >= 0; ) {
                        t = 0.0;
                        for (k = q; --k >= 0; ) t *= a[i][j]*b[j][k];
                        c[i][j] = t;
                    }
            }

A conformant array may be passed as a parameter to  a  func-
tion  provided  the function's prototype was visible to con-
firm that a conformant array parameter was intended.

     The simplest way of treating sizeof is to rule it  out:
if  the  description of 'a' involves any ?s, you can't apply
sizeof to it.  So given a parameter

        float fred[?f1][?f2][10];

sizeof fred and sizeof fred[1] would be illegal, but  sizeof
fred[1][2] and sizeof fred[1][2][3] would be legal.


David Hough

ARPA: dhough at sun.com
UUCP: {ucbvax,decvax,decwrl,seismo}!sun!dhough