Mathc matrices/h12c
Bibliothèque
Installer ce fichier dans votre répertoire de travail.
/* Save as : vdsvdcn.h */
/* ------------------------------------ */
void svd_U_Cn_uv_mR(
double **A,
double **U
int i; 
int r = rsize_R(A);
int c = csize_R(A);
double **A_T = i_mR(c,r); 
double **AA_T = i_mR(r,r);
double **EigsValue = i_mR(r,C1);
double **eignV = i_mR(r,r);
double **Ab;
double **b = m0_mR(i_mR(r,C1));
double **Ide = eye_mR(i_mR(r,r));
double **sIde = i_mR(r,r);
double **AA_TmnssIde = i_mR(r,r);
 eigs_value_AA_T_uv_mR(A, eignValue); if(c > r)
 eigs_value_AA_T_uv_mR(A_T,eignValue); if(r > c)
 eigs_value_AA_T_uv_mR(A,EigsValue); 
 transpose_mR(A,A_T) ; 
 mul_mR(A,A_T,AA_T);
 for(i = R1; i <= rsize_R(EigsValue); i++)
 smul_mR(EigsValue[i][C1],Ide,sIde);
 MmnsD_mR(AA_T,sIde,AA_TmnssIde);
 Ab = c_A_b_Ab_mR(AA_TmnssIde,b,i_Abr_Ac_bc_mR(r,r,C1));
 gj_freevariable_eignvector_uv_mR(Ab,eignV);
 c_c_mR(eignV,C2,U,i); 
 } 
 Normalize_uv_mR(U); 
 f_mR(EigsValue);
 f_mR(eignV);
 f_mR(A_T); 
 f_mR(AA_T);
 f_mR(Ab);
 f_mR(b);
 f_mR(Ide);
 f_mR(sIde);
 f_mR(AA_TmnssIde);
/* ------------------------------------ */
void svd_V_Cn_uv_mR(
double **A,
double **V
int i;
int j; 
int r = rsize_R(A);
int c = csize_R(A);
double **A_T = i_mR(c,r); 
double **A_TA = i_mR(c,c);
double **EigsValue = i_mR(c,C1);
double **eignV = i_mR(c,c);
double **Ab;
double **Ab0;
double **b = m0_mR(i_mR(c,C1));
double **Ide = eye_mR(i_mR(c,c));
double **sIde = i_mR(c,c);
double **A_TAmnssIde = i_mR(c,c);
 eigs_value_AA_T_uv_mR(A, eignValue); if(c > r)
 eigs_value_AA_T_uv_mR(A_T,eignValue); if(r > c)
 eigs_value_AA_T_uv_mR(A,EigsValue);
 transpose_mR(A,A_T);
 mul_mR(A_T,A,A_TA);
 for(i = R1; i <= rsize_R(EigsValue); i++)
 if(EigsValue[i][C1])
 smul_mR(EigsValue[i][C1],Ide,sIde);
 MmnsD_mR(A_TA,sIde,A_TAmnssIde);
 Ab = i_Abr_Ac_bc_mR(rsize_R(A_TAmnssIde),
 csize_R(A_TAmnssIde),
 csize_R(b));
 c_A_b_Ab_mR(A_TAmnssIde,b, Ab); 
 gj_freevariable_eignvector_uv_mR(Ab,eignV);
 c_c_mR(eignV,C2,V,i); 
 else
 smul_mR(EigsValue[i][C1],Ide,sIde);
 MmnsD_mR(A_TA,sIde,A_TAmnssIde);
 Ab0 = i_Abr_Ac_bc_mR(rsize_R(A_TAmnssIde), 
 csize_R(A_TAmnssIde), 
 (C1+c-r)

c_A_b_Ab_mR(A_TAmnssIde,b, Ab0);
 gj_free_V_Cn_value0_uv_mR(Ab0,eignV); 
 for(j=C2; i <= rsize_R(EigsValue); i++,j++)
 c_c_mR(eignV,j,V,i);
 } 
 Normalize_uv_mR(V); 
 f_mR(EigsValue);
 f_mR(eignV);
 f_mR(A_TA);
 f_mR(Ab);
 f_mR(Ab0);
 f_mR(b);
 f_mR(Ide);
 f_mR(sIde);
 f_mR(A_TAmnssIde);
/* ------------------------------------ */
double **pseudo_Cn_uv_mR(
double **A,
double **Pinv,
double factor_E
int r = rsize_R(A);
int c = csize_R(A);
double **a = smul_mR(factor_E,A,i_mR(r,c));
double **U = i_mR(r,r);
double **U_T = i_mR(r,r);
double **V = i_mR(c,c);
double **U_TA = i_mR(r,c);
double **S = i_mR(r,c);
double **S_T = i_mR(c,r);
double **invS_T = i_mR(c,r); 
double **VS_T = i_mR(c,r);
 svd_U_Cn_uv_mR(a,U); 
 svd_V_Cn_uv_mR(a,V);
/* S = U_T * B * V */
 transpose_mR(U,U_T); 
 mul_mR(U_T,A,U_TA);
 mul_mR(U_TA,V,S); 
/* invS_T = */ 
 transpose_mR(S,S_T);
 svd_inv_D_uv_mR(S_T,invS_T); 
/* PseudoInverse = V * invS_T * U_T */ 
 mul_mR(V,invS_T,VS_T);
 mul_mR(VS_T,U_T,Pinv); 
 f_mR(a); 
 f_mR(U);
 f_mR(V);
 f_mR(U_T); 
 f_mR(U_TA);
 f_mR(S);
 f_mR(S_T);
 f_mR(invS_T);
 f_mR(VS_T); 
 return(Pinv); 
/* ------------------------------------ */
 Déclaration des fichiers h.

Mathc matrices/h12d
Bibliothèque
Installer ce fichier dans votre répertoire de travail.
/* Save as : v_d.h */
/* ------------------------------------ */

Mathc matrices/c12an4
Application
Installer et compiler ces fichiers dans votre répertoire de travail.
/* ------------------------------------ */
/* Save as : cgn0f.c */
/* ------------------------------------ */
void fun(int r,int c)
double **B = r_mR(i_mR(c,r),9.);
double **B_T = transpose_mR(B,i_mR(r,c));
double **U_TB = i_mR(r,r); /* rc cr*/
double **U = i_mR(c,r);
double **V = i_mR(r,r);
double **S = i_mR(r,r);
double **S_T = i_mR(r,r);
double **invS_T = i_mR(r,r);
double **U_T = i_mR(r,c);
double **VS_T = i_mR(r,r); /* rr rr */
double **Pinv = i_mR(r,c); /* rr rr rc */
double **PinvB = i_mR(r,r); /* rr rr rc cr */
 clrscrn();
 printf(" Copy/Past into the octave windows \n\n\n");
 p_Octave_mR(B,"B",P2);
 printf("pinv(B)\n\n\n");
 stop();
 clrscrn(); 
 svd_U_Rn_uv_mR(B_T,U); 
 svd_V_Rn_uv_mR(B_T,V); 
 /* S = U_T * B * V */ 
 transpose_mR( U,U_T);
 mul_mR(U_T,B,U_TB);
 mul_mR(U_TB,V,S);
 transpose_mR(S,S_T);
 svd_inv_D_uv_mR(S_T,invS_T);
 /* Pinv = V * invS_T * U_T */ 
 mul_mR(V,invS_T,VS_T);
 mul_mR(VS_T,U_T,Pinv);
 printf(" PseudoInverse = V * invS_T * U_T "); 
 p_mR(Pinv,S10,P6,C10); 
 printf(" Ide = PseudoInverse * B "); 
 mul_mR(Pinv,B,PinvB); 
 p_mR(PinvB,S10,P6,C10); 
 f_mR(B);
 f_mR(B_T);
 f_mR(U_TB); 
 f_mR(U);
 f_mR(V);
 f_mR(S); 
 f_mR(S_T);
 f_mR(invS_T);
 f_mR(U_T); 
 f_mR(VS_T); 
 f_mR(Pinv);
 f_mR(PinvB); 
/* ------------------------------------ */
int main(void)
time_t t;
 srand(time(&t));
int i;
do
 i = rp_I(R3)+R1;
 fun(i,i+R0);
} while(stop_w());
 return 0;
/* ------------------------------------ */

Exemple de sortie écran :
 Copy/Past into the octave windows 
 B=[
-5.00,+8.00,-1.00,+4.00;
-5.00,-5.00,-9.00,-3.00;
+6.00,-3.00,+6.00,-3.00;
-1.00,-5.00,-7.00,-7.00]
pinv(B)
 Press return to continue. 
 PseudoInverse = V * invS_T * U_T 
 -0.418605 -0.569767 -0.612403 +0.267442 
 +0.034884 -0.244186 -0.143411 +0.186047 
 +0.302326 +0.383721 +0.534884 -0.220930 
 -0.267442 -0.127907 -0.344961 -0.093023 
 Ide = PseudoInverse * B 
 +1.000000 -0.000000 +0.000000 +0.000000 
 +0.000000 +1.000000 +0.000000 +0.000000 
 +0.000000 +0.000000 +1.000000 -0.000000 
 -0.000000 -0.000000 -0.000000 +1.000000 
 Press return to continue
 Press X to stop

Mathc matrices/c12b8
Application
Installer ce fichier dans votre répertoire de travail.
/* Save as : ccn1f.c */
/* ------------------------------------ */
void fun(int r,int c)
double **A = r_mR(i_mR(r,c),9.);
double **U = i_mR(r,r);
double **V = i_mR(c,c);
double **U_T = i_mR(r,r);
double **U_TA = i_mR(r,c); /* rr rc */
double **S = i_mR(r,c);
double **S_T = i_mR(c,r);
double **invS_T = i_mR(c,r); 
double **VS_T = i_mR(c,r); /* cc cr VS_T */
double **Pinv = i_mR(c,r); /* cc cr rr VS_TU_T */
double **Ide = i_mR(r,r); /* rc cc cr rr AVS_TU_T */
 clrscrn();
 printf(" Copy/Past into the octave windows \n\n\n");
 p_Octave_mR(A,"A",P2);
 printf("pinv(A)\n\n\n");
 stop();
 clrscrn();
 svd_U_Cn_uv_mR(A,U); 
 svd_V_Cn_uv_mR(A,V); 
 transpose_mR( U,U_T);
/* S = U_T * B * V */ 
 mul_mR(U_T,A,U_TA);
 mul_mR(U_TA,V,S);
 transpose_mR(S,S_T);
 svd_inv_D_uv_mR(S_T,invS_T); 
/* PseudoInverse = V * invS_T * U_T */ 
 mul_mR(V,invS_T,VS_T);
 mul_mR(VS_T,U_T,Pinv);
 printf(" PseudoInverse = V * invS_T * U_T "); 
 pE_mR(Pinv,S10,P6,C10); 
 printf(" Ide = A * PseudoInverse ");
 mul_mR(A,Pinv,Ide); 
 p_mR(Ide,S10,P6,C10); 
 f_mR(A); 
 f_mR(U);
 f_mR(V); 
 f_mR(U_T); 
 f_mR(U_TA);
 f_mR(S); 
 f_mR(S_T);
 f_mR(invS_T);
 f_mR(VS_T);
 f_mR(Pinv); 
 f_mR(Ide); 
/* ------------------------------------ */
int main(void)
time_t t;
 srand(time(&t));
int i;
do
 i = rp_I(R3)+R1;
 fun(i,i+C1);
} while(stop_w());
 return 0;
/* ------------------------------------ */
 Dans cet exemple on vérifie que PseudoInverse = V * invS_T * U_T
 Exemple de sortie écran :
 Copy/Past into the octave windows 
 A=[
+4.00,-3.00,+4.00,+2.00,-5.00;
+4.00,+8.00,+4.00,+6.00,-3.00;
-7.00,+8.00,+6.00,-7.00,-3.00;
+6.00,-3.00,-1.00,-5.00,+8.00]
pinv(A)
 Press return to continue. 
 PseudoInverse = V * invS_T * U_T 
+4.820956e-02 +4.597405e-02 -1.430147e-02 +6.523338e-02 
-6.775831e-02 +7.282082e-02 +2.292894e-02 +8.884062e-03 
+7.812878e-02 +1.295387e-02 +4.675018e-02 +3.143548e-02 
-2.877289e-02 +4.001486e-02 -5.762949e-02 -5.036124e-02 
-6.978350e-02 +1.945579e-02 -1.085020e-02 +5.186015e-02 
 Ide = A * PseudoInverse 
 +1.000000 +0.000000 +0.000000 -0.000000 
 +0.000000 +1.000000 -0.000000 -0.000000 
 +0.000000 +0.000000 +1.000000 -0.000000 
 -0.000000 -0.000000 +0.000000 +1.000000 
 Press return to continue
 Press X to stop