Plotting functions with matrices
4 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
function [Om2, Pp] = GeneralizedEigenProblem(Kk,Mm);
% Solves the generalized eigen problem (Kk - W^2 Mm) u = 0
% Inputs
% Kk: stiffness matrix
% Mm: mass matrix
%Torque Force
%Excitation Frequency
%Force array
Mm= [3.2 0 0 0 0; 0 5.4 0 0 0; 0 0 6.5 0 0; 0 0 0 12.5 0; 0 0 0 0 3.7];
Kk= [23 -12.2 -4.4 0 0; -12.2 27.8 -15.6 0 0; -4.4 -15.6 30.5 -7.9 -2.6; 0 0 -7.9 16 0; 0 0 -2.6 0 26.9];
%----- Construct dynamic matrix
%
K1 = inv(Kk); Dd = K1*Mm;
%
% Compute Eigenvalues and Eigenvectors
%
[Pp, Lam] = eig(Dd);
%
%----- Adjust eigen values
%
nn = size(Kk,1);
for i=1:nn
Om2(i) = 1/Lam(i,i);
end
%
%----- Order eigen values
%
swap = 1;
while (swap == 1)
swap = 0;
for i=1:nn-1
if (Om2(i) > Om2(i+1))
swap = 1;
tr0 = Om2(i); Om2(i) = Om2(i+1); Om2(i+1) = tr0;
tr1 = Pp(:,i); Pp(:,i) = Pp(:,i+1); Pp(:,i+1) = tr1;
end
end
end
%
%----- Normalize eigen modes
%
for i=1:nn
ui = Pp(:,i);
mu = sqrt(transpose(ui)*Mm*ui);
Pp(:,i) = Pp(:,i)/mu;
end
disp('Eigenvalues 1 through 3')
u1=Pp(:,1)
freq1=sqrt(Om2(:,1))
u2=Pp(:,2)
freq2=sqrt(Om2(:,2))
u3=Pp(:,3)
freq3=sqrt(Om2(:,3))
u4=Pp(:,4)
freq4=sqrt(Om2(:,4))
u5=Pp(:,5)
freq5=sqrt(Om2(:,5))
disp('P matrix')
Pp;
%Finding the Rayleigh's Quotient
%Randomly picked numbers for alpha
a1=5;
a3=9;
a4=6;
a5=-4;
e=(a1.*u1)+(a3.*u3)+(a4.*u4)+(a5.*u5);
ep=-linspace(-.1,.1);
rq= ((u2'.*Kk.*u2)+(2*ep.*u2'.*Kk.*e)+((ep^2).*e'.*Kk.*e))/((u2'.*Mm.*u2)+(2*ep.*u2'.*Mm.*e)+((ep^2).*e'.*Mm.*e));
figure(1)
plot(ep,rq)
grid;
5 Commenti
madhan ravi
il 18 Ott 2018
size(ep) is 1 by 100 And the rest is 5 by 1 , ep size should be 5 by 5 or 5 by 1
Risposte (1)
madhan ravi
il 18 Ott 2018
Modificato: madhan ravi
il 18 Ott 2018
function [Om2, Pp] = GeneralizedEigenProblem(Kk,Mm);
% Solves the generalized eigen problem (Kk - W^2 Mm) u = 0
% Inputs
% Kk: stiffness matrix
% Mm: mass matrix
%Torque Force
%Excitation Frequency
%Force array
Mm= [3.2 0 0 0 0; 0 5.4 0 0 0; 0 0 6.5 0 0; 0 0 0 12.5 0; 0 0 0 0 3.7];
Kk= [23 -12.2 -4.4 0 0; -12.2 27.8 -15.6 0 0; -4.4 -15.6 30.5 -7.9 -2.6; 0 0 -7.9 16 0; 0 0 -2.6 0 26.9];
%----- Construct dynamic matrix
%
K1 = inv(Kk); Dd = K1*Mm;
%
% Compute Eigenvalues and Eigenvectors
%
[Pp, Lam] = eig(Dd);
%
%----- Adjust eigen values
%
nn = size(Kk,1);
for i=1:nn
Om2(i) = 1/Lam(i,i);
end
%
%----- Order eigen values
%
swap = 1;
while (swap == 1)
swap = 0;
for i=1:nn-1
if (Om2(i) > Om2(i+1))
swap = 1;
tr0 = Om2(i); Om2(i) = Om2(i+1); Om2(i+1) = tr0;
tr1 = Pp(:,i); Pp(:,i) = Pp(:,i+1); Pp(:,i+1) = tr1;
end
end
end
%
%----- Normalize eigen modes
%
for i=1:nn
ui = Pp(:,i);
mu = sqrt(transpose(ui)*Mm*ui);
Pp(:,i) = Pp(:,i)/mu;
end
disp('Eigenvalues 1 through 3')
u1=Pp(:,1)
freq1=sqrt(Om2(:,1))
u2=Pp(:,2)
freq2=sqrt(Om2(:,2))
u3=Pp(:,3)
freq3=sqrt(Om2(:,3))
u4=Pp(:,4)
freq4=sqrt(Om2(:,4))
u5=Pp(:,5)
freq5=sqrt(Om2(:,5))
disp('P matrix')
Pp;
%Finding the Rayleigh's Quotient
%Randomly picked numbers for alpha
a1=5;
a3=9;
a4=6;
a5=-4;
e=(a1.*u1)+(a3.*u3)+(a4.*u4)+(a5.*u5);
ep=-linspace(-.1,.1,5);
rq= ((u2'.*Kk.*u2)+(2*ep.*u2'.*Kk.*e)+((ep.^2).*e'.*Kk.*e))/((u2'.*Mm.*u2)+(2*ep.*u2'.*Mm.*e)+((ep.^2).*e'.*Mm.*e));
figure(1)
plot(ep,rq)
grid;
end
2 Commenti
Vedere anche
Categorie
Scopri di più su Computational Geometry in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!