Linear Algebra and the C Language/a0gk

/* ------------------------------------ */
/*  Save as :  c00d.c                   */
/* ------------------------------------ */
#include "v_a.h"
/* ------------------------------------ */
#define  ARRAY  A2 
#define    RCA  R5 
#define    POW  3  
/* ------------------------------------ */
void fun(void)
{
double **A[ARRAY];
double **A_n[ARRAY];

double **P          =   r_mR(   i_mR(RCA,RCA),99); 
double **InvP       = inv_mR(P, i_mR(RCA,RCA)   );
double **T          =           i_mR(RCA,RCA);
double **InvP_A0n_P =           i_mR(RCA,RCA);
int i;

  for(i=A0; i<ARRAY; i++)
     {   
         A[i] = i_mR(RCA,RCA); 
       A_n[i] = i_mR(RCA,RCA);   
	   }
	    
  rsymmetric_mR(A[R0],99);
  
  mul_mR(InvP,A[A0],T); 
  mul_mR(T,P,A[A0+A1]);

   clrscrn();
   printf(" The two similar matrices :\n\n"
          " A[%d] = InvP A[%d] P   \n\n",A1,A0);
          
   for(i=A0; i<ARRAY; i++)
      {
	   printf(" A[%d] : ",i);      
       p_mR(A[i], S12,P2,C6); 
      }   
   stop();    
  
   clrscrn();
   
   for(i=A0; i<ARRAY; i++)
      {
       pow_mR(POW,A[i],A_n[i]);
	   }
	   
   mul_mR(InvP,A_n[A0],T);   
   mul_mR(T,P,InvP_A0n_P); 
      
   printf(" If   A[%d] = InvP A[%d] P   Then   "
          " A[%d]**%d = InvP A[%d]**%d P\n\n",
            A1,A0, A1,POW, A0,POW); 
            
   printf(" A[%d]**%d : ",A1,POW);      
   p_mR(A_n[1], S12,P2,C6); 
   
   printf(" InvP A[%d]**%d P : ",A0,POW);      
   p_mR(InvP_A0n_P, S12,P2,C6); 


  for(i=A0; i<ARRAY; i++)
     { 
	   f_mR(A[i]);
	   f_mR(A_n[i]);
	   }
	   
  f_mR(P);		 
  f_mR(InvP); 	   
  f_mR(T); 
  f_mR(InvP_A0n_P);       
}
/* ------------------------------------ */
int main(void)
{
time_t t;

  srand(time(&t));
  
  do{
        fun();
        
  }while(stop_w());

  return 0;
}
/* ------------------------------------ */
/* ------------------------------------ */

                                                                                       
 The two similar matrices :

 A[1] = InvP A[0] P   

 A[0] : 
      +32.00       +29.00       -93.00       +79.00        +5.00 
      +29.00       +24.00        -3.00       -17.00       +11.00 
      -93.00        -3.00       +90.00        -2.00       -74.00 
      +79.00       -17.00        -2.00       +48.00       +89.00 
       +5.00       +11.00       -74.00       +89.00       -60.00 

 A[1] : 
    -1069.53      +703.66     -2295.81     -1156.10      -224.37 
     -621.42      +373.86     -1322.60      -633.61      -118.25 
     +496.54      -204.15      +994.76      +479.69       +63.30 
     +335.79      -175.46      +732.37      +379.77       +64.87 
    -3701.98     +1796.98     -7475.12     -3572.85      -544.86 

 Press return to continue. 


 If   A[1] = InvP A[0] P   Then    A[1]**3 = InvP A[0]**3 P

 A[1]**3 : 
-36591496.94  +5453805.88 -57862479.38 -26022892.19    +33921.00 
-21684269.01  +2309842.55 -33840043.76 -15006567.92   +311882.34 
+17943996.69  +1576847.35 +26629252.44 +11263809.81  -1283894.72 
+12008119.64   -986211.15 +18711732.09  +8319472.35   -237431.71 
-130824153.47  +2592837.44 -198856855.35 -85986771.15  +5390651.60 

 InvP A[0]**3 P : 
-36591496.94  +5453805.88 -57862479.38 -26022892.19    +33921.00 
-21684269.01  +2309842.55 -33840043.76 -15006567.92   +311882.34 
+17943996.69  +1576847.35 +26629252.44 +11263809.81  -1283894.72 
+12008119.64   -986211.15 +18711732.09  +8319472.35   -237431.71 
-130824153.47  +2592837.44 -198856855.35 -85986771.15  +5390651.60 


 Press   return to continue
 Press X return to stop