/* DetProduct.c */
/* determinant of a product = product of determinants */

#include <stdio.h>

#define NN 220
#define INF 0XEFFFFF

void MatInput(int*);
void MatDisplay(int*);
int MatProd(int*,int*,int*);
int Det(int*);

void main(void)
{
 int M,N,   DVal[4],       Mat1[NN], Mat2[NN], Mat3[NN];

LP:
  printf("input size of (square) matrix nxn or 0x0 to terminate ");
   scanf("%dx%d",&M,&N); getchar();
  if (M != N) {puts("ERROR: matrix must be square"); goto LP;}
 if (N>12) {puts("size must be less than 13x13"); goto LP;}

 Mat1[0] = Mat2[0] = Mat3[0] = N*101;

 puts("\nInput first matrix");  MatInput(Mat1);  DVal[1]=Det(Mat1);
     MatDisplay(Mat1);   printf("Determinant value = %d\n",DVal[1]);
     
 puts("\nInput second matrix"); MatInput(Mat2);   DVal[2]=Det(Mat2);
    MatDisplay(Mat2);   printf("Determinant value = %d\n",DVal[2]);

 puts("\nCreate product of matrices");
 MatProd(Mat1,Mat2,Mat3);   DVal[3]=Det(Mat3);
 MatDisplay(Mat3); printf("Determinant value of product matrix = %d\n",DVal[3]);
 puts("\nProduct of determinant values of first two matrices:");
  printf("%d x %d = %d\n",DVal[1],DVal[2],DVal[1]*DVal[2]);

  puts("\n");
 goto LP;
}

/******************************************************************************/
/*****************************************************************************/

/********************** MatProd C = Mat A x Mat B *************************/

int MatProd(int *A, int *B, int *C)
{
 int I,J,K,M1,N1,M2,N2,Sum,   *PP, *P, *Q, *R;  

 M1 = A[0]/100; N1 = A[0]%100;  M2 = B[0]/100; N2 = B[0]%100;

 if (N1 != M2)
  {
   puts("ERROR:");
   puts("First matrix not conformable to second matrix for multiplication");
   return 0;
  }

  PP= &A[1];    C[0] = 100*M1 + N2;  R = &C[1];

 for (I=1; I<=M1; I++)
 {
  for (J=1; J<=N2; J++)
   {
    P=PP; Q=&B[J]; Sum=0;
    for (K=1; K<=N1; K++)
     {
      Sum = Sum + *P**Q;  /*  printf("I=%d J=%d *P=%d *Q=%d\n",I,J,*P,*Q); */
      P++; Q = Q + N2;
     }
   *R = Sum;  R++;          /*    printf("I=%d J=%d Sum = %d\n\n",I,J,Sum); */
   } /* end J */

  PP = PP + N1;

 } /* end I */
 *R = INF;

 return C[0];
}

/*****************************************************************************/

void MatInput(int *P)
{
 int I,J,K,M,N;

 M = N = *P%100;

 P++;
 for (I=1; I<=M; I++)
  {
   printf("Prepare to enter %d numbers in row %d\n",N,I);
   for (J=1; J<=N; J++)
    {
     if (J==1) puts("first number");
     if ((J>1) && (J<N)) puts("next number");
     if (J==N) puts("last number");
     scanf("%d",&K); getchar();
     *P=K; P++;
    }
  }

 *P=INF;
 return;
}

/*************************** Print A as a Matrix *****************************/

void MatDisplay(int *P)
 {
  int I,J,K,M,N,Max=0,T;    char S1[] = "%1d%2d%3d%4d%5d%6d%7d%8d%9d",S2[3];

  M=*P/100;  N=*P%100;  P++;
  for (I=1; I<=M; I++)
   {
    for (J=1; J<=N; J++)
     {
      if (*P>=0) K=*P; else K=-*P;
      if (K>Max) Max=K;
    }
  }

 T=1; K=1;
 while (T<=Max) {T=T*10; K++;} 
 for (I=0; I<3; I++) S2[I] = S1[3*K+I];  S2[3]='\0';


 for (I=1; I<=M; I++)
   {
    for (J=1; J<=N; J++) {printf(S2,*P); P++;}
    puts("");
  }
 return;
}

/**************************************************************************/

int Det(int *A)
   /*  NN has been defined globally so that NN-1 is the maximum length of
         any permutation (to prevent rollover of counters).
       Receive integer N, 2<=N<NN, and a pointer A to an array of length NN+1;
       The first permutation is odd; after that the odd-even alternate */
{
 int CNT,H,I,J,K,N,PROD,COPY;    int ODD=1,SUM=0;
 int R[]={0,0,0,0,0,0,0,0,0,0,0,0,0,0},
     D[]={1,1,1,1,1,1,1,1,1,1,1,1,1,1},
  PERM[]={0,1,2,3,4,5,6,7,8,9,10,11,12,13};

 N=A[0]%100;

 if (N<2 || N>NN)
  {printf("Input %d  not between limits  2 and %d\n",N,NN); return 0;}
 CNT=1; for (I=1; I<=N; I++) CNT*=I;   /* CNT = N! */

 while (CNT)
  {              /* form a sum of products of N! terms from matrix A */
   K=0; I=N;
  INDEX:      /* but first, construct a new PERMutation from the old one */
     R[I] = J = R[I] + D[I];
     if (J==I) {D[I]=-1; goto LOOP;}
     if (J) goto TRANSPOSE;
     D[I]=1; K++;
  LOOP:
     if (I>2) {I--; goto INDEX;}
     J=1;
  TRANSPOSE:
     J+=K; COPY=PERM[J]; PERM[J]=PERM[J+1]; PERM[J+1]=COPY;
     /* note: odd and even permutations alternate during the generation */

   /* Using permutation PERM form a product PROD of entries from the matrix A */
     PROD=1;
     for (H=1; H<=N; H++) PROD*= A[(H-1)*N + PERM[H]];   /* row H */
    /* PROD = A[1,PERM[I]]* **** * A[N,PERM[N]] */

/*If permutation even, PROD is added to SUM; if odd then -PROD is added */
     if(ODD) SUM-=PROD; else SUM+=PROD;
     ODD=1-ODD;
     CNT--;
  } /* end of while CNT */   
 return SUM;
}

