Mathc matrices/h20a
Application
Installer ce fichier dans votre répertoire de travail.
/* Save as : dot_au.h */
/* ------- <u.v> = v^t A^t A u ------- */
double dot_Au_R(
double **A,
double **u,
double **v
double **At = i_mR(csize_R(A), rsize_R(A));
double **vt = i_mR(R1, rsize_R(v));
double **T_vtAt = i_mR(R1, rsize_R(v));
double **T_vtAtA = i_mR(R1, rsize_R(v));
double **T_vtAtAu = i_mR(R1, C1);
double dotAuAv;
 transpose_mR(A,At);
 transpose_mR(v,vt);
 mul_mR(vt,At,T_vtAt);
 mul_mR(T_vtAt,A,T_vtAtA);
 mul_mR(T_vtAtA,u,T_vtAtAu);
 dotAuAv = T_vtAtAu[R1][C1];
 f_mR(At);
 f_mR(vt);
 f_mR(T_vtAt);
 f_mR(T_vtAtA);
 f_mR(T_vtAtAu);
 return (dotAuAv);
/* ------------------------------------ */
/* ------ ||u|| = <u.u>^(1/2) ------- */
double norm_Au_R(
double **A,
double **u
 return ( pow( dot_Au_R(A,u,u),(1./2.)) );
/* ------------------------------------ */
/* ------ d(u,v) = ||u-v||^(1/2) ------- */
double dist_Au_R(
double **A,
double **u,
double **v
double **umnsv = i_mR(rsize_R(v),csize_R(v));
double d;
 sub_mR(u,v,umnsv);
/* d(u,v) = ||u-v||^(1/2) */
 d = norm_Au_R(A,umnsv);
 f_mR(umnsv);
 return (d);
/* ------------------------------------ */
/* ------- <u.v> = v^t A^t A u ------- */
double **proj_Au_mR(
double **A,
double **u,
double **v,
double **projuv
 /* <u,v> / ||v||^2 */
double s = dot_Au_R(A,u,v) / dot_Au_R(A,v,v);
/* (<u,v>/||v||^2) v */
 smul_mR(s,v,projuv);
 return(projuv); 
/* ------------------------------------ */
/* ------- <u.v> = v^t A^t A u ------- */
double **orth_Au_mR(
double **A,
double **U,
double **Orth
double **Vect_U = i_mR(rsize_R(U),C1);
double **Vect_T = i_mR(rsize_R(U),C1);
double **Vect_Orth = i_mR(rsize_R(U),C1);
double **Vect_projuv = i_mR(rsize_R(U),C1);
int c_U; 
int c_Orth;
 c_c_mR(U,C1, Orth,C1);
 for( c_U = C2; c_U < U[C_SIZE][C0]; c_U++)
 c_c_mR(U,c_U, Vect_U,C1);
 c_c_mR(U,c_U, Vect_T,C1);
 for( c_Orth = C1; c_Orth < c_U; c_Orth++)
 c_c_mR(Orth,c_Orth, Vect_Orth,C1); 
 proj_Au_mR(A,Vect_U,Vect_Orth,Vect_projuv); 
 sub_mR(Vect_T,Vect_projuv,Vect_Orth);
 c_mR(Vect_Orth,Vect_T); 
 c_c_mR(Vect_T,C1, Orth,c_Orth);

f_mR(Vect_U); 
 f_mR(Vect_T); 
 f_mR(Vect_Orth); 
 f_mR(Vect_projuv); 
 return(Orth);
/* ------------------------------------ */
/* ------- <u.v> = v^t A^t A u ------- */
double **Normalize_Au_mR(
double **A,
double **Orth
double **u = i_mR(rsize_R(Orth),C1);
double **Nu = i_mR(rsize_R(Orth),C1); 
int c_Orth = csize_R(Orth);
int c;
 for( c = C1; c <= c_Orth; c++)
 c_c_mR(Orth,c, u,C1); 
 smul_mR(1./norm_Au_R(A,u),u,Nu);
 c_c_mR(Nu,C1, Orth,c); 
 f_mR(u); 
 f_mR(Nu); 
 return(Orth);
/* ------------------------------------ */
/* ------- <u.v> = v^t A^t A u ------- */
double **mul_Au_mR(
double **A,
double **U,
double **V,
double **UV
double **u = i_mR(rsize_R(U),C1);
double **v = i_mR(rsize_R(U),C1);
int c1;
int c2;
 for ( c1=C1; c1<U[C_SIZE][C0]; c1++)
 for ( c2=C1; c2<U[C_SIZE][C0]; c2++) 
 c_c_mR(U,c1, u,C1);
 c_c_mR(V,c2, v,C1);
 UV[c1][c2] = dot_Au_R(A,u,v);
 f_mR(u);
 f_mR(v); 
 return(UV);
/* ------------------------------------ */
/* ---------- u.v = v^t u ------------- */
void QR_Au_mR(
double **A,
double **U,
double **Q,
double **R)
 orth_Au_mR(A,U,Q);
 Normalize_Au_mR(A,Q); 
 mul_Au_mR(A,Q,U,R); 
/* ------------------------------------ */
/* ------------------------------------ */
double **eigs_value_Au_mR(
double **A,
double **U,
double **EigsValue
int r = rsize_R(A);
double **T = i_mR(r,r);
double **Q = i_mR(r,r);
double **R = i_mR(r,r);
int i;
int rc;
 c_mR(U,T);
 for(i=0;i<LOOP_EIGSUV;i++)
 QR_Au_mR(A,T,Q,R); 
 mul_mR(R,Q,T);

f_mR(T);
 f_mR(Q);
 f_mR(R); 
 return(EigsValue); 
/* ------------------------------------ */
double **xeigs_value_Au_mR(
double **A,
double **U,
double **EigsValue
int r = rsize_R(A);
double **T = i_mR(r,r);
double **Q = i_mR(r,r);
double **R = i_mR(r,r);
int i;
int rc;
 c_mR(U,T);
do
 QR_Au_mR(A,T,Q,R); 
 mul_mR(R,Q,T);
 clrscrn();
 pE_mR(T,S5,P3,C6); 
} while(stop_w());
 getchar();
 for(rc=R1;rc<=r;rc++)
 EigsValue[rc][C1] = T[rc][rc];

f_mR(T);
 f_mR(Q);
 f_mR(R); 
 return(EigsValue); 
/* ------------------------------------ */
 Déclaration des fichiers h.

Mathc matrices/h20b
Application
Installer ce fichier dans votre répertoire de travail.
/* Save as : dot_diag.h */
/* ------------------------------------ */
double **f_sqrtdiag_mR(
double **Diag,
double **SqrtDiag
int r;
 dif_sizes_mR(Diag,SqrtDiag,"f_sqrtdiag_mR;","Diag SqrtDiag");
 for ( r=R1; r<Diag[R_SIZE][C0]; r++)
 SqrtDiag[r][r] = sqrt(Diag[r][r]);
 return(SqrtDiag);
/* ------------------------------------ */
/* ------- <u.v> = v^t SqrtA_T SqrtA u */
/* SqrtA_T == SqrtA */
double dot_Diag_R(
double **A,
double **u,
double **v
double **SqrtA = i_mR(csize_R(A), rsize_R(A)); 
double **vt = i_mR(R1, rsize_R(v));
double **T_vtAt = i_mR(R1, rsize_R(v));
double **T_vtAtA = i_mR(R1, rsize_R(v));
double **T_vtAtAu = i_mR(R1, C1);
double dotDiag;
 transpose_mR(v,vt);
 f_sqrtdiag_mR(A,SqrtA); 
 mul_mR(vt,SqrtA,T_vtAt);
 mul_mR(T_vtAt,SqrtA,T_vtAtA);
 mul_mR(T_vtAtA,u,T_vtAtAu);
 dotDiag = T_vtAtAu[R1][C1];
 f_mR(vt);
 f_mR(T_vtAt);
 f_mR(T_vtAtA);
 f_mR(T_vtAtAu);
 return (dotDiag);
/* ------------------------------------ */
/* ------ ||u|| = <u.u>^(1/2) ------- */
double norm_Diag_R(
double **A,
double **u
 return ( pow( dot_Diag_R(A,u,u),(1./2.)) );
/* ------------------------------------ */
/* ------ d(u,v) = ||u-v||^(1/2) ------- */
double dist_Diag_R(
double **A,
double **u,
double **v
double **umnsv = i_mR(rsize_R(v),csize_R(v));
double d;
 sub_mR(u,v,umnsv);
/* d(u,v) = ||u-v||^(1/2) */
 d = norm_Diag_R(A,umnsv);
 f_mR(umnsv);
 return (d);
/* ------------------------------------ */
/* ------- <u.v> = v^t SqrtA_T SqrtA u */
double **proj_Diag_mR(
double **A,
double **u,
double **v,
double **projuv
 /* <u,v> / ||v||^2 */
double s = dot_Diag_R(A,u,v) / dot_Diag_R(A,v,v);
/* (<u,v>/||v||^2) v */
 smul_mR(s,v,projuv);
 return(projuv); 
/* ------------------------------------ */
/* ------- <u.v> = v^t SqrtA_T SqrtA u */
double **orth_Diag_mR(
double **A,
double **U,
double **Orth
double **Vect_U = i_mR(rsize_R(U),C1);
double **Vect_T = i_mR(rsize_R(U),C1);
double **Vect_Orth = i_mR(rsize_R(U),C1);
double **Vect_projuv = i_mR(rsize_R(U),C1);
int c_U; 
int c_Orth;
 c_c_mR(U,C1, Orth,C1);
 for( c_U = C2; c_U < U[C_SIZE][C0]; c_U++)
 c_c_mR(U,c_U, Vect_U,C1);
 c_c_mR(U,c_U, Vect_T,C1);
 for( c_Orth = C1; c_Orth < c_U; c_Orth++)
 c_c_mR(Orth,c_Orth, Vect_Orth,C1); 
 proj_Diag_mR(A,Vect_U,Vect_Orth,Vect_projuv); 
 sub_mR(Vect_T,Vect_projuv,Vect_Orth);
 c_mR(Vect_Orth,Vect_T); 
 c_c_mR(Vect_T,C1, Orth,c_Orth);