Aller au contenu

Mathc complexes/a345

Un livre de Wikilivres.


Décomposition d'une matrice en éléments propres

Installer et compiler ces fichiers dans votre répertoire de travail.


c00d.c
/* ------------------------------------ */
/*  save as :   c00d.c                  */
/* ------------------------------------ */
#include "w_a.h"
/* ------------------------------------ */
#define FACTOR_E        +1.E-4  
/* ------------------------------------ */        
/* ------------------------------------ */
double f(
double x)
{  
        return(1./x);
}
char  feq[] = "1./x";
/* ------------------------------------ */
/* ------------------------------------ */
void fun(int rc)
{
double **A      = rcsymmetric_mZ(            i_mZ(rc,rc),9);
double **b      =           r_mZ(            i_mZ(rc,C1),99);

double **A_T    =  ctranspose_mZ(A,           i_mZ(rc,rc));
double **U      =         X_U_mZ(A_T,         i_mZ(rc,rc),FACTOR_E);
double **U_T    =  ctranspose_mZ(U,           i_mZ(rc,rc));

double **U_TA   =         mul_mZ(U_T, A,      i_mZ(rc,rc));
double **X      =         mul_mZ(U_TA, U,     i_mZ(rc,rc));
            
double **invX   =     f_eigs_mZ(f,X,          i_mZ(rc,rc));
double **UinvX  =        mul_mZ(U,invX,       i_mZ(rc,rc));
double **inv    =        mul_mZ(UinvX,U_T,    i_mZ(rc,rc));

double **x      =        mul_mZ(inv, b,       i_mZ(rc,C1));

  clrscrn();  
  printf(" A :");
  p_mZ(A, S9,P2, S8,P2, C4);
  printf(" b :");
  p_mZ(b, S8,P2, S6,P2, C4); 
  printf(" inv :");
  pE_mZ(inv, S12,P4, S12,P4, C4);   
  stop();  
 
  clrscrn();       
  printf(" x = inv b  ");   
  pE_mZ(x, S12,P3, S12,P3, C4);  

  printf(" B :");
  p_mZ(b, S8,P2, S8,P2, C4);
  printf("  Ax :");
  p_mZ(mul_mZ(A,x,b), S8,P2, S8,P2, C4);
     
   f_mZ(A);
   f_mZ(A_T);

   f_mZ(b);
   f_mZ(x);
     
   f_mZ(U);
   f_mZ(U_T);
   f_mZ(U_TA);
   f_mZ(X);
   
   f_mZ(invX);
   f_mZ(UinvX);
   f_mZ(inv);  
}
/* ------------------------------------ */
int main(void)
{
time_t t;

  srand(time(&t));

do
{
    fun(R3);
    
} while(stop_w());

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


Exemple de sortie écran :
 A :
  +275.00   +0.00i    -21.00 -129.00i   +107.00   -5.00i 
   -21.00 +129.00i   +106.00   +0.00i    -10.00  +87.00i 
  +107.00   +5.00i    -10.00  -87.00i    +74.00   +0.00i 

 b :
   +5.00+78.00i 
  +78.00 +6.00i 
  -50.00+18.00i 

 inv :
 +8.8183e-03 +0.0000e+00i  +2.4691e-03 +1.4462e-02i  +4.5855e-03 -3.5273e-04i 
 +2.4691e-03 -1.4462e-02i  +4.4727e-01 +7.0473e-18i  +5.7848e-02 -5.0476e-01i 
 +4.5855e-03 +3.5273e-04i  +5.7848e-02 +5.0476e-01i  +6.0811e-01 +0.0000e+00i 

 Press return to continue. 


 x = inv b  
  -7.302e-02  +1.931e+00i 
  +4.222e+01  +2.908e+01i 
  -2.893e+01  +5.102e+01i 

 B :
   +5.00  +78.00i 
  +78.00   +6.00i 
  -50.00  +18.00i 

  Ax :
   +5.00  +78.00i 
  +78.00   +6.00i 
  -50.00  +18.00i 


 Press   return to continue
 Press X return to stop