Linear Algebra and the C Language/a0n7


Install and compile this file in your working directory. 

/* ------------------------------------ */
/*  Save as:  c00g.c                    */
/* ------------------------------------ */
#include "v_a.h"
/* ------------------------------------ */
#define ARRAY   A3
/* ------------------------------------ */
#define RC     RC3
/* ------------------------------------ */
void fun(void)
{
double **A         = rsymmetric_mR(        i_mR(RC,RC),999);
double **AEVect    =     eigs_V_mR(A,      i_mR(RC,RC));
double **InvAEVect =  transpose_mR(AEVect, i_mR(RC,RC));
double **AEValue   =       eigs_mR(A,      i_mR(RC,C1)); 

double **T         =                       i_mR(RC,RC);
double **a         =                       i_mR(RC,RC);

double **b    [ARRAY];
double **r    [ARRAY];
double **br   [ARRAY];
double **EVabr[ARRAY];

int i;

  for(i=A0; i<ARRAY; i++)
     {
          b[i] =  c_c_mR(   AEVect,i+C1,              i_mR(RC,C1),C1);
          r[i] =  c_r_mR(InvAEVect,i+R1,              i_mR(R1,RC),R1);
      
         br[i] =  mul_mR(b[i],r[i],                   i_mR(RC,RC));
         
      EVabr[i] = smul_mR(1./AEValue[i+R1][C1], br[i], i_mR(RC,RC)); 
      }
      
  clrscrn();
  printf(" A:");
  p_mR(A, S10,P4,C6);

  printf(" AEVect:        Eingvectors of A");
  p_mR(AEVect, S10,P4,C6);
  
  printf(" AEValue:       Eignvalues of A");
  p_mR(AEValue, S10,P4,C6);  
  stop();

  clrscrn(); 
  printf(" inv(A):");
  inv_mR(A, T);
  pE_mR(T, S10,P4,C6);
  
  add_mR(EVabr[0], EVabr[1], T); 
  add_mR(     T, EVabr[2], a);  
  printf(" We can calculate the inverse of A\n"
         " by simply taking the inverse of  \n"
         " each of the eigenvalues.   \n\n\n\n"
         " 1/E1 b1r1 + 1/E2 b2r2 + 1/E3 b3r3 = inv(A)");
  pE_mR(a, S10,P4,C6);
   
  f_mR(A);
  f_mR(AEVect);
  f_mR(InvAEVect);
  f_mR(AEValue);
     
  f_mR(T);
  f_mR(a);
  
  for(i=A0; i<ARRAY; i++)
     {
       f_mR(    b[i]);
       f_mR(    r[i]);
       f_mR(   br[i]); 
       f_mR(EVabr[i]);          
     } 
}
/* ------------------------------------ */
int main(void)
{
time_t t;

  srand(time(&t));

do
{
 fun();

} while(stop_w());

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

A property of spectral decomposition:

     1/E1 * b1r1 + 1/E2 * b2r2 + 1/E3 * b3r3 = inv(A)
     
We can calculate the inverse of A by simply taking the inverse of each of the eigenvalues.

  b: The columns of the eigenvectors of A
  r: The rows of eigenvectors of the inverse of A
  E: The eignvalues of A  

  b1r1 is obtained by multiplying the first column of the eigenvector of A by the first row of the eigenvector of the inverse of A.

Screen output example: 

                                                                                       
 A:
 +426.0000  +704.0000  +935.0000 
 +704.0000   +66.0000   +22.0000 
 +935.0000   +22.0000  +823.0000 

 AEVect:        Eingvectors of A
   +0.6591    +0.6910    +0.2969 
   +0.2897    -0.5975    +0.7477 
   +0.6940    -0.4067    -0.5940 

 AEValue:       Eignvalues of A
+1720.0870 
 -733.1394 
 +328.0524 

 Press return to continue. 


 inv(A):
-1.3013e-04 +1.3508e-03 +1.1173e-04 
+1.3508e-03 +1.2657e-03 -1.5685e-03 
+1.1173e-04 -1.5685e-03 +1.1301e-03 

 We can calculate the inverse of A
 by simply taking the inverse of  
 each of the eigenvalues.   



 1/E1 b1r1 + 1/E2 b2r2 + 1/E3 b3r3 = inv(A)
-1.3013e-04 +1.3508e-03 +1.1173e-04 
+1.3508e-03 +1.2657e-03 -1.5685e-03 
+1.1173e-04 -1.5685e-03 +1.1301e-03 


 Press   return to continue
 Press X return to stop
This article is issued from Wikibooks. The text is available under Creative Commons Attribution-Share Alike 4.0 unless otherwise noted. Additional terms may apply for the media files.