Linear Algebra and the C Language/a0d7

/* ------------------------------------ */
/*  Save as :  c00c.c                   */
/* ------------------------------------ */
#include "v_a.h"
/* ------------------------------------ */
#define   RCA        RC4
/* ------------------------------------ */
/* ------------------------------------ */
double f(
double x)
{
        return(1.0/(x));
}
char  feq[] =  "1.0/(x)";
/* ------------------------------------ */
/* ------------------------------------ */
void fun(void)
{
double **A             = rsymmetric_mR(            i_mR(RCA,RCA),999.);
double **invA          =                           i_mR(RCA,RCA);

double **EigsVector    =     eigs_V_mR(A,          i_mR(RCA,RCA));   
double **T_EigsVector  =  transpose_mR(EigsVector, i_mR(RCA,RCA)); 

double **EigsValue     =                           i_mR(RCA,RCA);
double **invEigsValue  =                           i_mR(RCA,RCA);

double **T             =                           i_mR(RCA,RCA);

  clrscrn();
  printf(" A :");
  p_mR(A,S10,P2,C6);  
  
/* EigsValue = T_EigsVector * A * EigsVector */  
  mul_mR(T_EigsVector,A,T);
  mul_mR(T,EigsVector,EigsValue);   
    
  f_eigs_mR(f,EigsValue,invEigsValue);

/* invA = EigsVector * invEigsValue * T_EigsVector */ 
  mul_mR(EigsVector,invEigsValue,T);
  mul_mR(T,T_EigsVector,invA);
  
  printf(" invA = EigsVector * invEigsValue * T_EigsVector");  
  pE_mR(invA,S10,P4,C6); 
  
  printf(" Ide = A * invA "); 
  p_mR(mul_mR(A,invA,T),S10,P5,C10);     
   
  f_mR(A);
  f_mR(invA);
    
  f_mR(EigsVector);
  f_mR(T_EigsVector);
  
  f_mR(EigsValue);
  f_mR(invEigsValue);    
     
  f_mR(T);
}
/* ------------------------------------ */
int main(void)
{
time_t t;

  srand(time(&t));

do
{
  fun();

} while(stop_w());

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

                                                                                       
 A :
   +421.00    -342.00    +330.00    +597.00 
   -342.00    -623.00    -182.00    -388.00 
   +330.00    -182.00    +660.00    +216.00 
   +597.00    -388.00    +216.00    -182.00 

 invA = EigsVector * invEigsValue * T_EigsVector
+2.2760e-04 -8.1928e-04 -8.3237e-04 +1.5053e-03 
-8.1928e-04 -1.1110e-03 +1.4953e-04 -1.4138e-04 
-8.3237e-04 +1.4953e-04 +2.1395e-03 -5.0998e-04 
+1.5053e-03 -1.4138e-04 -5.0998e-04 -8.6062e-04 

 Ide = A * invA 
  +1.00000   -0.00000   +0.00000   -0.00000 
  +0.00000   +1.00000   +0.00000   +0.00000 
  -0.00000   -0.00000   +1.00000   -0.00000 
  -0.00000   -0.00000   +0.00000   +1.00000 


 Press   return to continue
 Press X return to stop