www.pudn.com > simulator.zip > General_eigen.cpp

#include "General_ext.h" 
 
void eigen_finder(float[][MAX_CTRL_CHAN], float[], float[][MAX_CTRL_CHAN], int); 
 
/**************************************************************************/ 
/* Function: mode_solution                                                */ 
/* Purpose:	 solve natural periods and mode shapes of an [M]-[K] system   */ 
/**************************************************************************/ 
void mode_solution(float eig_val[MAX_CTRL_CHAN]) 
{ 
const float pi=3.14159f; 
float Stiff_i [MAX_CTRL_CHAN][MAX_CTRL_CHAN], 
	  temp    [MAX_CTRL_CHAN][MAX_CTRL_CHAN], 
      period  [MAX_CTRL_CHAN], 
      eig_vec [MAX_CTRL_CHAN][MAX_CTRL_CHAN]; 
 
//  Mass[][]*v[]=lambda*Stif[][]*v[] to find eigs(lambda,v[])  
//     with lambda=(1/omega)^2, equals to A[][]*v[]=lambda*v[]  
//     with A[][]=Flex[][]*Mass[][], where Flex[][]=Stif[][]^(-1) 
 
    invert_matrix(Stiff, Stiff_i); 
    mpy_matrix_matrix(Stiff_i, Mass, temp); 
    eigen_finder(temp, eig_val, eig_vec, 10); 
 
    printf("\n  No.      Periods                        Mode Shapes\n"); 
    for(int i=0; i4) && (fabs(d[ip])+g==fabs(d[ip])) && (fabs(d[iq])+g==fabs(d[iq])) ) 
			{ 
               a[ip][iq]=0.f; 
			} 
            else if (fabs(a[ip][iq])>tresh) 
			{ 
               h=d[iq]-d[ip]; 
 
               if (fabs(h)+g==fabs(h)) 
			   { 
                  t=a[ip][iq]/h; 
			   } 
               else 
			   { 
                  theta=0.5f*h/a[ip][iq]; 
                  t=1.f/((float)fabs(theta)+(float)sqrt(1.f+theta*theta)); 
			      { 
                     if (theta<0.f) t=-t; 
			   	  } 
			   }; 
 
               c=1.f/(float)sqrt(1.f+t*t); 
               s=t*c; 
               tau=s/(1.f+c); 
               h=t*a[ip][iq]; 
               z[ip]-=h; 
               z[iq]+=h; 
               d[ip]-=h; 
               d[iq]+=h; 
               a[ip][iq]=0.f; 
 
               for( j=0 ; j