Aller au contenu

Mathc matrices/h05b

Ajouter des liens
Un livre de Wikilivres.


Bibliothèque



Installer ce fichier dans votre répertoire de travail.

vdotu.h
/* ------------------------------------ */
/*  Save as :   vdotu.h                */
/* ------------------------------------ */

/* ------------------------------------ */
/* ------- u.v = v^t u ---------------- */
/* ------------------------------------ */
double dot_R(
double **u,
double **v
)
{
double **vt  = i_mR(R1, rsize_R(v));
double **vtu = i_mR(R1, C1);

double udotv;
  
  transpose_mR(v,vt);
        mul_mR(vt,u,vtu);

  udotv = vtu[R1][C1];
  
  f_mR(vt);
  f_mR(vtu);
  
  return (udotv);
}
/* ------------------------------------ */
/* ------  ||u|| =  (u.u)^(1/2) ------- */
/* ------------------------------------ */
double norm_R(
double **u
)
{
    return (pow(dot_R(u,u),(1./2.)));
}
/* ------------------------------------  */
/* -- d(u,v) = ||u-v||                -- */
/* ------------------------------------  */
double dist_R(
double **u,
double **v
)
{
double **umnsv  =  i_mR(rsize_R(v),csize_R(v));

double d;
  
  sub_mR(u,v,umnsv);

  d = norm_R(umnsv);
  
  f_mR(umnsv);
  
  return (d);
}
/* ------------------------------------ */
/*        (<u,v>/||v||^2)  * v          */
/* ------------------------------------ */
double **proj_mR(
double **u,
double **v,
double **projuv
)
{
           /* <u,v>      /  ||v||^2  */
double s = dot_R(u,v) / dot_R(v,v);


    smul_mR(s,v,projuv);
    
  return(projuv); 
}
/* ------------------------------------ */
/* ------------------------------------ */
double **normalize_mR(
double **U,
double **Q
)
{
double **u  = i_mR(rsize_R(U),C1);
double **q  = i_mR(rsize_R(U),C1);	

int c_U     = csize_R(U);
int c;

		for( c = C1; c <= c_U; c++)
          {		
		    c_c_mR(U,c, u,C1); 
		    		           
		    smul_mR(1./norm_R(u),u,q);

            c_c_mR(q,C1, Q,c); 
	      }

  f_mR(u);  
  f_mR(q);  
  
 return(Q);
}
/* ------------------------------------ */
/* ------------------------------------ */
double **Normalize_mR(
double **Q
)
{
double **u  = i_mR(rsize_R(Q),C1);
double **q  = i_mR(rsize_R(Q),C1);	

int c_Q   = csize_R(Q);
int c;

		for( c = C1; c <= c_Q; c++)
          {		
		    c_c_mR(Q,c, u,C1); 
		    		           
		    smul_mR(1./norm_R(u),u, q);

            c_c_mR(q,C1, Q,c); 
	      }

  f_mR(u);  
  f_mR(q);  
  
 return(Q);
}
/* ------------------------------------ */
/* ------------------------------------ */
double **orth_mR(
double **U,
double **Q,
int      Normal  /* YES = 1 NO = 0 */
)
{
double **u  = i_mR(rsize_R(U),C1);
double **v  = i_mR(rsize_R(U),C1);

double **q  = i_mR(rsize_R(U),C1);

double **projuv = i_mR(rsize_R(U),C1);
double **tq     = i_mR(rsize_R(U),C1);


int c_U;	
int c_Q;

  c_c_mR(U,C1, Q,C1);    
 
  for( c_U = C2; c_U < U[C_SIZE][C0]; c_U++)
      {
		c_c_mR(U,c_U, u,C1);
		c_c_mR(U,c_U, q,C1);

		for( c_Q = C1; c_Q < c_U; c_Q++)
          {		
		    c_c_mR(Q,c_Q, v,C1); 

		    proj_mR(u,v, projuv);      /* [<u,v> / ||v||^2] * v */		                		        
		        sub_mR(q,   projuv, tq);		    
		          c_mR(             tq, q);    
	      }
	      c_c_mR(q,C1, Q,c_Q);		
	  }
	  
  if(Normal) Normalize_mR(Q); 
	   
  f_mR(u);
  f_mR(v);
  f_mR(q); 
  
  f_mR(projuv);  
  f_mR(tq);
           
 return(Q);
}
/* ------------------------------------ */
/* ------------------------------------ */
/* ------------------------------------ */
/* ------------------------------------ 
double **xorth_mR(
double **U,
double **Q
)
{
double **u  = i_mR(rsize_R(U),C1);
double **v  = i_mR(rsize_R(U),C1);

double **q  = i_mR(rsize_R(U),C1);

double **projuv = i_mR(rsize_R(U),C1);

int c_U;	
int c_Q;

  c_c_mR(U,C1, Q,C1);    
 
  for( c_U = C2; c_U < U[C_SIZE][C0]; c_U++)
      {
		c_c_mR(U,c_U, u,C1);
		c_c_mR(U,c_U, q,C1);

		for( c_Q = C1; c_Q < c_U; c_Q++)
          {		
		    c_c_mR(Q,c_Q, v,C1); 

		    proj_mR(u,v, projuv); 		                		        
		        sub_mR(q,projuv, v);
		    
		    c_mR(v,q);    
	      }
	      c_c_mR(q,C1, Q,c_Q);		
	  }
	  
  Normalize_mR(Q); 
	    
  f_mR(u);
  f_mR(v);
  f_mR(q); 
  
  f_mR(projuv);  
         
 return(Q);
}*/
/* ------------------------------------ */
/* ------------------------------------ */
void QR_mR(
double **U,
double **Q,
double **R)
{
double **Q_T = i_mR(csize_R(Q),rsize_R(Q));

  orth_mR(U,Q,YES); 
  transpose_mR(Q,Q_T);
  mul_mR(Q_T,U, R);
  
  f_mR(Q_T);    
}
/* ------------------------------------ */
double **eigs_mR(
double **A,
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(A,T);
 
 for(i=0;i<LOOP_EIGSUV;i++)
     { 
	   QR_mR(T,Q,R);  
        mul_mR(R,Q,T);
     }
     
  for(rc=R1;rc<=r;rc++)
     { 
	   EigsValue[rc][C1] = T[rc][rc];
     }
     
  f_mR(T);
  f_mR(Q);
  f_mR(R); 
  
  return(EigsValue);    
}
/* ------------------------------------ */
double **sqrt_svd_mR(
double **A,
double **sqrtA
)
{
int r;
int c;

  dif_sizes_mR(A,sqrtA,"svd_sqrt_mR();","(A or sqrtA)");

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

            sqrtA[r][c] = sqrt(A[r][c]);
            
 return(sqrtA);
}
/* ------------------------------------ */
double **inv_svd_mR(
double **S,
double **invS
)
{
int r;
int c;

  dif_sizes_mR(S,invS,    "inv_svd_mR();","(S,invS)");

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

         if(r==c)
            {
			 if(fabs(S[r][c]) < ERROR_E) 
			             invS[r][c] = 0;
             else        invS[r][c] = (1./S[r][c]);
		    }
         else
            invS[r][c] = 0;
         
 return(invS);
}
/* ------------------------------------ */
/* ------------------------------------ */
/* ------------------------------------ */
double **svds_mR(
double **A,
double **SvdValue)
{
int rA = rsize_R(A);
int cA = csize_R(A);

double **A_T    =      i_mR(cA,rA);
double **S      =      i_mR(cA,cA);

double **SvdValue_2 =  i_mR(cA,C1);

  transpose_mR(A,A_T);
  mul_mR(A_T,A, S);   
  eigs_mR(S,SvdValue_2);
  sqrt_svd_mR(SvdValue_2,SvdValue);
  
  f_mR(A_T);
  f_mR(S);
  f_mR(SvdValue_2);
  
  return(SvdValue);    
}
/* ------------------------------------ */
/* ------------------------------------ */



Déclaration des fichiers h.
Récupérée de « https://fr.wikibooks.org/w/index.php?title=Mathc_matrices/h05b&oldid=706150 »