Mathc complexes/h05c
Bibliothèque
Installer ce fichier dans votre répertoire de travail.
/* Save as : wpoctave.h */
void p_Octave_mZ(
double **A,
char name[],
int Pr,
int Pi
 int r;
 int c;
 printf(" %s=[\n",name);
 for(r=R1; r<A[R_SIZE][C0]; r++)
 for(c=C1; c<A[C_SIZE][C0]; c+=C2)
 printf("%+.*f%+.*f*i,",Pr,A[r][c], Pi,A[r][c+C1]);
 if((r+1)<A[R_SIZE][C0]) printf("\b;\n");
 else printf("\b]\n\n");
void pE_Octave_mZ(
double **A,
char name[],
int Pr,
int Pi
 int r;
 int c;
 printf(" %s=[\n",name);
 for(r=R1; r<A[R_SIZE][C0]; r++)
 for(c=C1; c<A[C_SIZE][C0]; c+=C2)
 printf("%+.*e%+.*e*i,",Pr,A[r][c], Pi,A[r][c+C1]);
 if((r+1)<A[R_SIZE][C0]) printf("\b;\n");
 else printf("\b]\n\n");
void pAi_Octave_mZ(
double **A,
char name[],
int i,
int Pr,
int Pi
 int r;
 int c;
 printf(" %s%d=[\n",name,i);
 for(r=R1; r<A[R_SIZE][C0]; r++)
 for(c=C1; c<A[C_SIZE][C0]; c+=C2)
 printf("%+.*f%+.*f*i,",Pr,A[r][c], Pi,A[r][c+C1]);
 if((r+1)<A[R_SIZE][C0]) printf("\b;\n");
 else printf("\b]\n\n");
/* ------------------------------------ */
/* ------------------------------------ */
void p_Octave_mRZ(
double **A,
char name[],
int Pr
 int r;
 int c;
 printf(" %s=[\n",name);
 for(r=R1; r<A[R_SIZE][C0]; r++)
 for(c=C1; c<A[C_SIZE][C0]; c+=C2)
 printf("%+.*f,",Pr,A[r][c]);
 if((r+1)<A[R_SIZE][C0]) printf("\b;\n");
 else printf("\b]\n\n");
void pE_Octave_mRZ(
double **A,
char name[],
int Pr
 int r;
 int c;
 printf(" %s=[\n",name);
 for(r=R1; r<A[R_SIZE][C0]; r++)
 for(c=C1; c<A[C_SIZE][C0]; c+=C2)
 printf("%+.*e,",Pr,A[r][c]);
 if((r+1)<A[R_SIZE][C0]) printf("\b;\n");
 else printf("\b]\n\n");
/* ------------------------------------ */
/* ------------------------------------ */
 Déclaration des fichiers h.

Mathc complexes/h05d
Bibliothèque
Installer ce fichier dans votre répertoire de travail.
/* Save as : wddotu.h */
/* ------- u.v = conj(v)^t u ---------- */
nb_Z dotuv_Z(
double **u,
double **v
double **vt = i_mZ(R1, rsize_Z(v));
double **vtu = i_mZ(R1, C1);
nb_Z udotv;
 ctranspose_mZ(v,vt);
 mul_mZ(vt,u,vtu);
 udotv = i_Z( vtu[R1][C1],vtu[R1][C1+C1]);
 f_mZ(vt);
 f_mZ(vtu);
 return (udotv);
/* ------ ||u|| = (u.u)^(1/2) ------- */
double normuv_Z(
double **u
double x = Re_Z(dotuv_Z(u,u));
 return (pow(x,(1./2.)));
/* -- d(u,v) = ||u-v|| -- */
double distuv_Z(
double **u,
double **v
double **umnsv = i_mZ(rsize_Z(v),csize_Z(v));
double d;
 sub_mZ(u,v,umnsv);
 d = normuv_Z(umnsv);
 f_mZ(umnsv);
 return (d);
/* ------- (<u,v>/||v||^2).v ---------- */
double **projuv_mZ(
double **u,
double **v,
double **projuv
nb_Z z = div_Z(dotuv_Z(u,v),dotuv_Z(v,v)); 
 return(zmul_mZ(z,v,projuv)); /* (<u,v>/||v||^2).v */
/* ------- u.v = conj(v)^t u ---------- */
double **orthuv_mZ(
double **U,
double **Orth
double **Vect_U = i_mZ(rsize_Z(U),C1);
double **Vect_T = i_mZ(rsize_Z(U),C1);
double **Vect_Orth = i_mZ(rsize_Z(U),C1);
double **Vect_projuv = i_mZ(rsize_Z(U),C1);
int c_U; 
int c_Orth;
 c_c_mZ(U,C1, Orth,C1);
 for( c_U = C2; c_U <= csize_Z(U); c_U++)
 c_c_mZ(U,c_U, Vect_U,C1);
 c_c_mZ(U,c_U, Vect_T,C1);
 for( c_Orth = C1; c_Orth < c_U; c_Orth++)
 c_c_mZ(Orth,c_Orth, Vect_Orth,C1); 
 projuv_mZ(Vect_U,Vect_Orth,Vect_projuv); 
 sub_mZ(Vect_T,Vect_projuv,Vect_Orth);
 c_mZ(Vect_Orth,Vect_T); 
 c_c_mZ(Vect_T,C1, Orth,c_Orth);

f_mZ(Vect_U); 
 f_mZ(Vect_T); 
 f_mZ(Vect_Orth); 
 f_mZ(Vect_projuv); 
 return(Orth);
/* ------- u.v = conj(v)^t u ---------- */
double **Normalizeuv_mZ(
double **Orth
double **u = i_mZ(rsize_Z(Orth),C1);
double **Nu = i_mZ(rsize_Z(Orth),C1); 
int c_Orth = csize_Z(Orth);
int c;
 for( c = C1; c <= c_Orth; c++)
 c_c_mZ(Orth,c, u,C1); 
 smul_mZ(1./normuv_Z(u),u,Nu);
 c_c_mZ(Nu,C1, Orth,c); 
 f_mZ(u); 
 f_mZ(Nu); 
 return(Orth);
/* ------------------------------------ */
/* ------- u.v = conj(v)^t u ---------- */
double **Normalizeuv2_mZ(
double **Orth,
double **Norm
double **u = i_mZ(rsize_Z(Orth),C1);
double **Nu = i_mZ(rsize_Z(Orth),C1); 
int c_Orth = csize_Z(Orth);
int c;
 for( c = C1; c <= c_Orth; c++)
 c_c_mZ(Orth,c, u,C1); 
 smul_mZ(1./normuv_Z(u),u,Nu);
 c_c_mZ(Nu,C1, Norm,c); 
 f_mZ(u); 
 f_mZ(Nu); 
 return(Norm);
/* ------------------------------------ */
/* ------- u.v = conj(v)^t u ---------- */
void QRuv_mZ(
double **U,
double **Q,
double **R)
double **Q_T = i_mZ(csize_Z(Q),rsize_Z(Q));
 orthuv_mZ(U,Q);
 Normalizeuv_mZ(Q); 
 ctranspose_mZ(Q,Q_T);
 mul_mZ(Q_T,U, R);
 f_mZ(Q_T); 
/* ------------------------------------ */
double **eigsuv_mZ(
double **A,
double **EigsValue
int r = rsize_Z(A);
double **T = i_mZ(r,r);
double **Q = i_mZ(r,r);
double **R = i_mZ(r,r);
int i;
int re, ce;
 c_mZ(A,T);
 for(i=0;i<1000;i++)
 QRuv_mZ(T,Q,R); 
 mul_mZ(R,Q,T);

for(re=R1,ce=C1;re<=r;re++,ce+=C2)
 EigsValue[re][C1] = T[re][ce];
 EigsValue[re][C2] = T[re][ce+C1];

f_mZ(T);
 f_mZ(Q);
 f_mZ(R); 
 return(EigsValue); 
/* ------------------------------------ */
double **svduv_mZ(
double **A,
double **SvdValue)
int rA = rsize_Z(A);
int cA = csize_Z(A);
double **A_T = i_mZ(cA,rA);
double **S = i_mZ(cA,cA);
double **SvdValue_2 = i_mZ(cA,C1);
 ctranspose_mZ(A,A_T);
 mul_mZ(A_T,A, S); 
 eigsuv_mZ(S,SvdValue_2);
 sqrt_svd_mZ(SvdValue_2,SvdValue);
 f_mZ(A_T);
 f_mZ(S);
 f_mZ(SvdValue_2);
 return(SvdValue); 
/* ------------------------------------ */
 Déclaration des fichiers h.

Mathc complexes/c05a
Installer et compiler ces fichiers dans votre répertoire de travail.
/* ------------------------------------ */
/* Save as : c01a.c */
void fun(int r, int rn)
double **A = r_mZ(i_mZ(r+rn,r),9.);
double **Q = i_mZ(r+rn,r);
double **R = i_mZ(r,r);
 clrscrn();
 printf(" Copy/Past into the octave windows \n\n\n");
 p_Octave_mZ(A,"a",P0,P0); 
 printf(" [Q, R] = qr (a,0) \n\n");
 printf(" A :");
 p_mZ(A, S3,P0, S3,P0, C10);
 stop();
 clrscrn(); 
 QRuv_mZ(A,Q,R); 
 printf(" Q :");
 p_mZ(Q, S10,P4, S10,P4, C2); 
 printf(" R :");
 p_mZ(R, S10,P4, S10,P4, C2); 
 stop(); 
 clrscrn();
 printf(" A :");
 p_mZ(A, S3,P0, S3,P0, C10); 
 printf(" A = Q * R :"); 
 mul_mZ(Q,R, A);
 p_mZ(A, S3,P0, S3,P0, C10);

f_mZ(A);
 f_mZ(Q);
 f_mZ(R);
/* ------------------------------------ */
int main(void)
time_t t;
 srand(time(&t));
do
 fun(rp_I(R3)+R1,rp_I(R3));
} while(stop_w());
 return 0;
/* ------------------------------------ */
 Vérifier avec octave les calculs.

Cette fonction n'est correcte que pour les matrices ayant un nombre de lignes égales ou supérieures aux nombres de colonnes. Utiliser pour calculer les valeurs propres et les valeurs singulière, elle fait bien son travail puisque dans les deux cas elle travaille sur des matrices carrées.