Category : Files from Magazines
Archive   : DDJ1089.ZIP
Filename : LETTERS.LST

_LETTERS TO THE EDITOR_

[LISTING ONE]

/***********************************************************************/
/* Determinant Calculator Using Near Heap and Upper Triangular Matrix */

double det(arg,n) /* arg = array name, n = # elements */
char *arg;
int n;
{
register int i,t,k,p;
double **a, ret, x; /* Array, Return value of Det, Temp variable */
char **sdim2(); /* defined in Listing Three */
int * row; /* row pointer for pivoting */

/* dynamically create 2 dim "array " array a from arg */
a = (double **)sdim2(arg,n,n,sizeof(double));

row = (int *)malloc(n*sizeof(int)); /* row pointer for pivoting allocation */
if (row == (int *) NULL) {
fprintf(stderr,"No heap space for Row Pointers for Pivoting ");
exit(1);
}

/* creating upper triangular matrix with test for 0 values on diagonal */

/* first initialize pointer vector */
for (i = 0; i < n ; i++)
row[i] = i ;

/* find pivot elements */
for (p = 0; p < n - 1; p++) {
for (k = p + 1;k < n ;k++) {
if ( fabs(a[row[k]][p]) > fabs(a[row[p]][p]) ) { /* switch index */
t = row[p];
row[p] = row[k];
row[k] = t;
}
}
/* In usual application this would be an error message that the
* matrix is singular and the program would exit here */
if ( a[row[p]][p] == 0 ) /* Determinant is 0 */
break; /* No need to continue on */

for (k = p+1;k < n; k++) {
x = a[row[k]][p]/a[row[p]][p];
for ( i = p + 1; i < n ;i++) /* do Gauss Jordan elimination */
a[row[k]][i] = a[row[k]][i] - x * a[row[p]][i];
}
}

/* if ( a[row[n-1]][n-1] == 0 ) Determinant is 0 - This would in
* normal application be a message Matrix is singular and an exit */

/* value of determinant is product of diagonals of upper triangular matrix */
for (ret = 1,i = 0;i < n; i++)
ret *= a[row[i]][i]; /* if any of diagonals are 0 Det = 0 */

free(row);
free(a);

return ret;
}