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 **EigsValue_cosh2A_mZ(
double **EigsValue,
double **EigsValue_cosh2A
)
{
int  r;
int  c;

nb_Z cosh2A;

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

        if(((r-R1)*C2)==c-C1)
          {
           cosh2A = cosh_Z(i_Z(2.*EigsValue[r][c],2.*EigsValue[r][c+C1]));
           
           EigsValue_cosh2A[r][c   ] = cosh2A.r;
           EigsValue_cosh2A[r][c+C1] = cosh2A.i;
          }
              
 return(EigsValue_cosh2A);
}
/* ------------------------------------ */
/* ------------------------------------ */
double **EigsValue_coshAP2_pls_sinhAP2_mZ(
double **EigsValue,
double **EigsValue_coshAP2_pls_sinhAP2
)
{
int  r;
int  c;

nb_Z coshAP2_pls_sinhAP2;

nb_Z cosh;
nb_Z coshP2;

nb_Z sinh;
nb_Z sinhP2;

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

        if(((r-R1)*C2)==c-C1)
          {
           cosh   = cosh_Z(i_Z(EigsValue[r][c],EigsValue[r][c+C1]));
           coshP2 = mul_Z(cosh,cosh);			  
			  
           sinh   = sinh_Z(i_Z(EigsValue[r][c],EigsValue[r][c+C1]));
           sinhP2 = mul_Z(sinh,sinh);
           
           coshAP2_pls_sinhAP2 = add_Z(coshP2,sinhP2);
           
           EigsValue_coshAP2_pls_sinhAP2[r][c   ] = coshAP2_pls_sinhAP2.r;
           EigsValue_coshAP2_pls_sinhAP2[r][c+C1] = coshAP2_pls_sinhAP2.i;
          }
              
 return(EigsValue_coshAP2_pls_sinhAP2);
}
/* ------------------------------------ */
/* ------------------------------------ */
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 **EigsValue                     = i_mZ(RCA,RCA);
double **EigsValue_cosh2a              = i_mZ(RCA,RCA);
double **EigsValue_coshAP2_pls_sinhAP2 = i_mZ(RCA,RCA);

double **cosh2a                        = i_mZ(RCA,RCA);
double **coshAP2_pls_sinhAP2           = 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 :");
  eigs_V_mZ(A,V,FACTOR_E);
  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(" EigsValue = invV * A * V");
  mul_mZ(invV,A,T);
  mul_mZ(T,V,EigsValue); 
  pE_mZ(clean_eyes_mZ(EigsValue), S12,P4, S12,P4, C3);
  stop();

  clrscrn();  
  printf(" cosh(2*A)        :");
  EigsValue_cosh2A_mZ(EigsValue,EigsValue_cosh2a);
   
  mul_mZ(V,EigsValue_cosh2a,T);
  mul_mZ(T,invV,cosh2a); 
  pE_mZ(cosh2a, S12,P4, S8,P4, C3);

  printf(" cosh(A)**2+sinh(A)**2 :");
  EigsValue_coshAP2_pls_sinhAP2_mZ(EigsValue,EigsValue_coshAP2_pls_sinhAP2);
   
  mul_mZ(V,EigsValue_coshAP2_pls_sinhAP2,T);
  mul_mZ(T,invV,coshAP2_pls_sinhAP2); 
  pE_mZ(coshAP2_pls_sinhAP2, S12,P4, S8,P4, C3);  
  
  f_mZ(A);
  
  f_mZ(V);  
  f_mZ(invV);  
  f_mZ(T);  
  
  f_mZ(EigsValue);
  
  f_mZ(EigsValue_cosh2a);   
  f_mZ(cosh2a); 
  
  f_mZ(EigsValue_coshAP2_pls_sinhAP2);   
  f_mZ(coshAP2_pls_sinhAP2);     
}
/* ------------------------------------ */
int main(void)
{
time_t t;

  srand(time(&t));

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

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


Exemple de sortie écran :
 A :
  +0.0631 -0.0058i   +0.0825 +0.0468i   +0.0862 -0.0453i 
  -0.0571 +0.0469i   -0.0689 +0.0290i   +0.0960 -0.0539i 
  +0.0488 +0.0110i   -0.0587 +0.0919i   +0.0777 +0.0401i 

 V :
 -1.9061e-01 -7.1959e-01i  +7.3448e-01 -4.0806e-01i  -1.8452e-01 -2.4072e-01i 
 +5.1565e-01 -2.3206e-01i  +1.2721e-01 +1.8056e-01i  +8.7213e-02 +7.9767e-01i 
 +3.5512e-01 +1.9713e-17i  +4.9522e-01 +0.0000e+00i  +5.1392e-01 -1.9746e-18i 

 inv(V) ... Some time the matrix is not invertible :
 -4.9168e-01 +3.1511e-01i  +6.8143e-01 +6.1835e-01i  +5.1999e-01 -1.2798e+00i 
 +8.9119e-01 -1.0436e-01i  -4.0095e-01 +3.7864e-01i  +1.0246e+00 +9.3803e-01i 
 -5.1901e-01 -1.1717e-01i  -8.4504e-02 -7.9213e-01i  +5.9920e-01 -1.9592e-02i 

 EigsValue = invV * A * V
 +4.8615e-02 +1.0711e-01i  +0.0000e+00 +0.0000e+00i  +0.0000e+00 +0.0000e+00i 
 +0.0000e+00 +0.0000e+00i  +1.1056e-01 +1.8408e-02i  +0.0000e+00 +0.0000e+00i 
 +0.0000e+00 +0.0000e+00i  +0.0000e+00 +0.0000e+00i  -8.7270e-02 -6.2222e-02i 

 Press return to continue. 


 cosh(2*A)        :
 +1.0036e+00-1.6439e-03i  -5.0415e-03+2.4424e-02i  +4.8360e-02-6.5205e-03i 
 +9.1007e-03-6.2711e-03i  +9.9261e-01+1.8343e-02i  +3.3921e-03+2.5588e-02i 
 +1.1150e-02-9.5216e-03i  -6.6912e-03-1.1782e-04i  +1.0168e+00+3.3984e-02i 

 cosh(A)**2+sinh(A)**2 :
 +1.0036e+00-1.6439e-03i  -5.0415e-03+2.4424e-02i  +4.8360e-02-6.5205e-03i 
 +9.1007e-03-6.2711e-03i  +9.9261e-01+1.8343e-02i  +3.3921e-03+2.5588e-02i 
 +1.1150e-02-9.5216e-03i  -6.6912e-03-1.1782e-04i  +1.0168e+00+3.3984e-02i 


 Press   return to continue
 Press X return to stop