|
|
Line 1: |
Line 1: |
| Dear All,
| | Please ask questions at Octave's mailing-list, see [[FAQ#Where_do_I_get_additional_help.3F]]. |
| I am trying to write an transient analysis solver for a multi degree of freedom system. it consists of mass, stiffness, damping, gyroscopic damping matrices and so on. my initial condition is declared by zeros(2*nr,1). after I run the code, the output y which is a matrix of 6x41 has all it's rows and columns just zero. It would be really apprectiated if some one has an advice / opinion about what am I doing wrong?
| |
| for the main .m file it is as below:
| |
| | |
| %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
| |
| alpha = [0.0020*2*pi 8*pi 0];
| |
| tspan = [0 4] ;
| |
| frunb=zeros(size(K,1),1); %%% unbalance vector
| |
| frunb(5,1) = 1; %%%
| |
| nr = 3;
| |
| options = odeset;
| |
| [t,y] = ode45(@deriv,tspan,zeros(2*nr,1),options,M,K,C,G,alpha,frunb);
| |
| %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
| |
| | |
| | |
| and below is the function deriv (saved in a different file named as deriv.m)
| |
| | |
| | |
| %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
| |
| function dz = deriv(t,z,M,K,C,G,alpha,ubforce)
| |
|
| |
| phi = alpha(1)*t*t + alpha(2)*t + alpha(3);
| |
| d_phi = 2*alpha(1)*t + alpha(2);
| |
| dd_phi = 2*alpha(1);
| |
| A = [(-inv(M)*(C-j*phi*G)) (-inv(M)*K); eye(size(M)) zeros(size(M))];
| |
| [U,D] = eig(A);
| |
| uncoupled_A = (inv(U))*A*(U);
| |
| B = [inv(M)*ubforce;zeros(size(M))*ubforce];
| |
| uncoupled_B = inv(U) * B;
| |
| double_uncoupled_A_1 = [(uncoupled_A(1,1) + uncoupled_A(2,2)) (-uncoupled_A(1,1)*uncoupled_A(2,2)); 1 0];
| |
| double_uncoupled_A_2 = [(uncoupled_A(3,3) + uncoupled_A(4,4)) (-uncoupled_A(3,3)*uncoupled_A(4,4)); 1 0];
| |
| double_uncoupled_A_3 = [(uncoupled_A(5,5) + uncoupled_A(6,6)) (-uncoupled_A(5,5)*uncoupled_A(6,6)); 1 0];
| |
| U11 = [uncoupled_A(1,1) uncoupled_A(2,2); 1 1];
| |
| U22 = [uncoupled_A(3,3) uncoupled_A(4,4); 1 1];
| |
| U33 = [uncoupled_A(5,5) uncoupled_A(6,6); 1 1];
| |
| U = U(1:6,1:6);
| |
| U1 = blkdiag(U11,U22,U33);
| |
| double_uncoupled_A = blkdiag(double_uncoupled_A_1,double_uncoupled_A_2,double_uncoupled_A_3);
| |
| double_uncoupled_B = U1 * inv(U) * B(1:6,:);
| |
| gain = [1 0 0 0 0 0] * (phi^2)*exp(j*phi*t)*z;
| |
| dz = (double_uncoupled_A * z) + double_uncoupled_B * gain';
| |
| endfunction
| |