----- Lecture Notes: Storing Matrices in 2-D Arrays 
  ==================================================
   Illustrates the difference between an array and
   a matrix, and how to store a matrix in an array.
   -------------------------------------------------
   
          A matrix is a two-dimensional table of data; a
   2-D array is a set of contiguous memory storage locations, 
   each indexed by a pair of integers. Matrices are stored
   inside arrays; usually, a matrix will be stored in
   in the "upper left-hand corner" of the array, which
   means each row of the matrix begins at the beginning
   of a row of the array. 

The following picture gives an accurate visualization
of the relationship of the elements of the matrix to
the array containing them:

                -----------------------------
               |               |             |
               |               |             |
               |    MATRIX     |             | 
               |               |             |
               |---------------              |
               |                 ARRAY       |
               |                             |
                -----------------------------

              Figure 1:  A matrix inside an array.

Notice that there are 4 important sizes here.  We name them as follows.
    const int N_ROWS_A  =  the number of rows in the array.  
    const int N_COLS_A  =  the number of columns in the array.
    int nrows           =  the number of rows in the matrix.
    int ncols           =  the number of columns in the matrix.

The programmer must be careful to enforce the constraints 
nrows <= N_ROWS_A and ncols <= N_COLS_A.

Because the dimensions of the array, N_ROWS_A and N_COLS_A, must 
be specified at compile time, they are stored as named constants.
Because the same array may store different matrices of different
sizes at different points in the program, the matrix dimensions,
nrows and ncols, are stored as variables. 

Storing the matrix in only part of the array is of great practical
importance. Otherwise, programs would have to constantly reset
storage parameters at run-time to make sure every user always filled
every array. This would be impractical and very tedious to code. 

      Care must be taken to ensure that the matrix is properly 
   stored inside the upper left-hand corner of the array.
   Without extra row delimiters in the initialization,
   the compiler simply fills the array one element at 
   time, without grouping the rows as they appear in 
   the matrix.  See Example init2.cpp.

Although Figure 1 correctly describes most partially filled 2-D arrays,
it does not fully express how a 2-D array is actually stored in memory.
In fact, a 2-D array in C or C++ is stored as a one-dimensionsal array
of one-dimensional arrays, each of which represents a "row" of the
2-dimensional table.  Here is a more realistic picture of the storage 
of a matrix in an array.  Elements denoted by 'm' indicate matrix entries, 
and elements denoted by 'A' indicate unused array space.


|-----Array row 0 -----||------Array row 1-----||-----Array row 2------|| ..

|----------------------||----------------------||----------------------|| ..
 m m m m m A A A A A A   m m m m m A A A A A A   m m m m m A A A A A A   
|----------------------||----------------------||----------------------|| ..


|--- ----|              |--- ----|              |--- ----|             
    |                       |                       |
  Matrix                  Matrix                  Matrix
   row 0                   row 1                   row 2


      Figure 2:  A 2-D array and its matrix, as stored in memory.


For example, suppose you are storing a 3x5 matrix inside a 7x11 array.
Suppose you want the computer to access matrix element a[ 2 ][ 3 ]. 
In order to do this, the computer 

1) Starts at the first array element, a[ 0 ][ 0 ].
2) Skips over the first 22 elements to get to row 2 (the third row).
3) Skips over another 3 elements to get to the fourth element in row 2,
   which is a[ 2 ][ 3 ].

Thus, to access a specific array element a[i][j], the compiler
must "know" how many elements are in each ROW of the ARRAY. Note, however, 
that the number of elements in each column of the array is irrelevant 
to this task. This explains why when 2-D arrays are used as function 
parameters, the COLUMN dimension of the ARRAY  
(i.e., the number of elements in each row of the array) 
must be numerically specified at compile time. 
The row dimension of the array is ignored by the compiler.  
It should therefore be omitted from the parameter list. 

  =================================================================

For example, here are the prototype, sample call, and header for a 
function to write an integer array to the screen (write_iats()):

	------------------------------------------------
Prototype:
	void write_iats( int [  ][ N_COLS_A ], int, int );   

Sample call:
	write_iats( array1, N_ROWS_M, N_COLS_M );

Header:
	void write_iats( int a[ ][ N_COLS_A ], int n_rows_m, int n_cols_m)
	------------------------------------------------

The function must know the column dimension of the array and both
row and column dimension of the matrix, but the row dimension of 
the array is irrelevant.

  =================================================================

In C and C++, 2-D arrays are "stored by row" as one-dimensional "arrays
of arrays". To access the elements of the array sequentially in memory,
you should traverse the ROWS of the matrix one at a time. This is far
more efficient than walking down the matrix "column by column".
It is important to remember this when writing loops that access the 
elements of the matrix one by one, because such loops are easily written
in either row-by-row OR column-by-column orientation.

For example,
consider the following function (see Example linalg.cpp) for adding a scalar
multiple s of one real matrix A to another matrix B ( B <-- B + s*A ).


const int N_COLS_A;      // These named constants should be global
const int N_COLS_B;      

//=============== Begin dmaxpy() ==============================

void dmaxpy( double A[ ][ N_COLS_A ], double B[ ][ N_COLS_B ], 
             double s, int nrows, int ncols )
{

// ---------------- Replace B with B + s*A ----------------
//                  using ROW-ORIENTED code

  int i;                 // i = row index
  int j;                 // j = column index

  for (i = 0; i < nrows; ++i)           //  For each fixed row i ...
    {
       for (j = 0; j < ncols; ++j)      //  Walk across row i,
          B[i][j] += s * A[i][j] ;      //  updating B[i][j] at each step.
    }
}
//=============== End dmaxpy() ==============================


The **column-oriented** equivalent for the dmaxpy() function body is: 

{
  int i;                 // i = row index
  int j;                 // j = column index

  for (j = 0; j < ncols; ++j)         //  For each fixed column j ...
    {
       for (i = 0; i < nrows; ++i)    //  Walk down column j,
          B[i][j] += s * A[i][j] ;    //  updating B[i][j] at each step.
    }                                 //   (inefficient in C and C++).
}
	(The function header is unchanged.)

 While this code will perform the desired computation in a
fashion numerically equivalent to the row-oriented version, 
the mechanical overhead associated with "jumping around" in
the array will make its execution considerably slower. While you
may not notice the difference for a single addition of 2 small matrices,
you might well notice significantly diminished performance if your code 
requires several thousand calls to such a function, or a few calls 
with very large matrices. In practice, such circumstances are not 
unusual.

Summary: in C and C++, write loops accessing arrays in
	 ROW-ORIENTED style.

  =================================================================

   Finally, you should know that in some computer languages
(notably, FORTRAN), 2-D arrays are stored by column rather than
by row. In this case, it is the row dimension that must be numerically 
specified in the function that receives the array in its parameter list.
Moreover, arrays stored by column should be accessed by column
for best performance. Thus, loops accessing arrays in these 
languages should be written in column-oriented style.

  =================================================================