Fonctions matricielles dans C. Matrices Non symétriques conjugués


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


c00a.c
/* ------------------------------------ */
/*  Save as :   c00a.c                  */
/* ------------------------------------ */
#include "w_a.h"
/* ------------------------------------ */
#define FACTOR_E    +1.E-0         
#define RCA          RC3  
/* ------------------------------------ */
/* ------------------------------------ */
double **EV_sin2A_mZ(
double **EV,
double **EV_sin2A
)
{
int  r;
int  c;
nb_Z sin2A;

  for (   r=R1; r<EV[R_SIZE][C0]; r++)
    for ( c=C1; c<EV[C_SIZE][C0]; c+=C2)

        if(((r-R1)*C2)==c-C1)
          {
           sin2A = sin_Z(i_Z(2.*EV[r][c],2.*EV[r][c+C1]));
           
           EV_sin2A[r][c   ] = sin2A.r;
           EV_sin2A[r][c+C1] = sin2A.i;
          }
              
 return(EV_sin2A);
}
/* ------------------------------------ */
/* ------------------------------------ */
double **EV_TwocosAsinA_mZ(
double **EV,
double **EV_TwocosAsinA
)
{
int  r;
int  c;

nb_Z cosa;
nb_Z sina;

nb_Z TwocosAsinA;

  for (   r=R1; r<EV[R_SIZE][C0]; r++)
    for ( c=C1; c<EV[C_SIZE][C0]; c+=C2)

        if(((r-R1)*C2)==c-C1)
          {
           cosa = cos_Z(i_Z( EV[r][c], EV[r][c+C1]));
           sina = sin_Z(i_Z( EV[r][c], EV[r][c+C1]));
             
           TwocosAsinA = smul_Z(2., mul_Z(cosa,sina) ); 
           
           EV_TwocosAsinA[r][c   ] = TwocosAsinA.r;
           EV_TwocosAsinA[r][c+C1] = TwocosAsinA.i;     
          }
              
 return(EV_TwocosAsinA);
}
/* ------------------------------------ */    
/* ------------------------------------ */
void fun(void)
{                               
double **A              =      i_mZ(RCA,RCA);

double **V              =      i_mZ(RCA,RCA);
double **invV           =      i_mZ(RCA,RCA);
double **T              =      i_mZ(RCA,RCA);

double **EV             =      i_mZ(RCA,RCA);
double **EV_sin2A       =      i_mZ(RCA,RCA);
double **EV_TwocosAsinA =      i_mZ(RCA,RCA);

double **sin2A          =      i_mZ(RCA,RCA);
double **TwocosAsinA    =      i_mZ(RCA,RCA);

   do
  {
       rE_mZ(A,999,+1.E-4);  
   eigs_V_mZ(A,V,FACTOR_E);
  }while(!det_Z(V).r);
  
  clrscrn();
  printf(" A :");
  p_mZ(A, S9,P4, S8,P4, C3);

  printf(" V :");
  pE_mZ(V, S12,P4, S12,P4, C3);
  
  printf(" inv(V) ... Some time the matrix is not invertible :");
  inv_mZ(V,invV);
  pE_mZ(invV, S12,P4, S12,P4, C3);
  
  printf(" EV = invV * A * V");
  mul_mZ(invV,A,T);
  mul_mZ(T,V,EV); 
  pE_mZ(clean_eyes_mZ(EV), S12,P4, S12,P4, C3);
  stop();

  clrscrn();  
  printf(" sin(2*A)        :");
  EV_sin2A_mZ(EV,EV_sin2A);
   
  mul_mZ(V,EV_sin2A,T);
  mul_mZ(T,invV,sin2A); 
  pE_mZ(sin2A, S12,P4, S8,P4, C3);

  printf(" 2 cos(A)sin(A) :");
  EV_TwocosAsinA_mZ(EV,EV_TwocosAsinA);
   
  mul_mZ(V,EV_TwocosAsinA,T);
  mul_mZ(T,invV,TwocosAsinA); 
  pE_mZ(TwocosAsinA, S12,P4, S8,P4, C3);  
  
  f_mZ(A);
  
  f_mZ(V);  
  f_mZ(invV);  
  f_mZ(T);  
  
  f_mZ(EV);
  
  f_mZ(EV_sin2A);   
  f_mZ(sin2A); 
  
  f_mZ(EV_TwocosAsinA);   
  f_mZ(TwocosAsinA);     
}
/* ------------------------------------ */
int main(void)
{
time_t t;

  srand(time(&t));

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

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


Exemple de sortie écran :
 A :
  +0.0408 +0.0203i   -0.0127 +0.0039i   +0.0447 -0.0388i 
  +0.0360 -0.0293i   -0.0202 +0.0255i   -0.0349 +0.0807i 
  +0.0318 +0.0755i   -0.0161 -0.0050i   -0.0040 -0.0864i 

 V :
 +7.8458e-01 -2.5405e-01i  -7.4294e-02 +3.6063e-01i  +3.0020e-02 -1.5326e-01i 
 -1.2918e-02 +1.2384e-01i  -6.5551e-01 -2.7760e-01i  +6.1197e-01 -7.7270e-01i 
 +5.5172e-01 +1.1485e-17i  +5.9805e-01 -4.1308e-17i  +6.3525e-02 -3.5264e-18i 

 inv(V) ... Some time the matrix is not invertible :
 +8.3162e-01 +6.8435e-01i  -2.0647e-01 -2.9465e-02i  +3.0334e-01 -5.4459e-01i 
 -8.0663e-01 -5.6198e-01i  +1.1882e-01 -6.5890e-02i  +1.3938e+00 +3.9953e-01i 
 +3.7130e-01 -6.5289e-01i  +6.7456e-01 +8.7622e-01i  -1.4909e-02 +9.6851e-01i 

 EigsValue = invV * A * V
 +7.7485e-02 +2.8254e-03i  +0.0000e+00 +0.0000e+00i  +0.0000e+00 +0.0000e+00i 
 +0.0000e+00 +0.0000e+00i  -3.8152e-02 -6.3650e-02i  +0.0000e+00 +0.0000e+00i 
 +0.0000e+00 +0.0000e+00i  +0.0000e+00 +0.0000e+00i  -2.2734e-02 +2.0224e-02i 

 Press return to continue. 


 sin(2*A)        :
 +8.0976e-02+4.0490e-02i  -2.5309e-02+7.7647e-03i  +8.9424e-02-7.7565e-02i 
 +7.1831e-02-5.8941e-02i  -4.0377e-02+5.0973e-02i  -6.9379e-02+1.6166e-01i 
 +6.3612e-02+1.5091e-01i  -3.2168e-02-9.9614e-03i  -8.5821e-03-1.7273e-01i 

 2 cos(A)sin(A) :
 +8.0976e-02+4.0490e-02i  -2.5309e-02+7.7647e-03i  +8.9424e-02-7.7565e-02i 
 +7.1831e-02-5.8941e-02i  -4.0377e-02+5.0973e-02i  -6.9379e-02+1.6166e-01i 
 +6.3612e-02+1.5091e-01i  -3.2168e-02-9.9614e-03i  -8.5821e-03-1.7273e-01i 


 Press   return to continue
 Press X return to stop