i am not getting the exact data for finite element model updating

10 visualizaciones (últimos 30 días)
Susanta  Barik
Susanta Barik el 29 de Mzo. de 2022
Comentada: Susanta Barik el 29 de Mzo. de 2022
clear;
close all;
clc;
% __________ M3=40tons
% | | K3=50*10^6
% ---------- M2=45tons
% | | K2=50*10^6
% ---------- M1=50tons
% | | K1=50*10^6
% --- ---
set-up for stiffness matrix and other
m1=50000;
m2=50000;
m3=50000;
mass=diag([m1 m2 m3]);
fprintf('The mass matrix \n');
disp(mass)
k1=50*10^6;
k2=50*10^6;
k3=50*10^6;
stiffness=[k1+k2 -k2 0;-k2 k2+k3 -k3;0 -k3 k3];
fprintf('The stiffness matrix \n');
disp(stiffness)
Natural frequency(eigen value) and Mode shape(eigen vector) calculation
[v,d]=eig(stiffness,mass);
fprintf('lambda value \n');
disp(d);
w=[sqrt(d)];
w1=w(1,1);
w2=w(2,2);
w3=w(3,3);
W=[w1,w2,w3];
fprintf('The natural frequencies of this structure are as follows(HZ) \n');
disp(W)
fprintf('The mode shape of the structure \n');
disp(v)
fprintf('lambda value \n');
disp(d);
w=[sqrt(d)];
w1=w(1,1);
w2=w(2,2);
w3=w(3,3);
W=[w1,w2,w3];
fprintf('The natural frequencies of this structure are as follows(HZ) \n');
disp(W)
fprintf('The mode shape of the structure \n');
disp(v)
Individusal mode shapes vector
fprintf('individusal mode shapes \n');
v1=v(:,1)
v2=v(:,2)
v3=v(:,3)
Plotting of mode shape
% height=[0;3;6;9];
% for i=1:5
% subplot(1,3,i)
% plot([0; v(:,i)],height);
% xlabel('mode shape','FontSize',12);
% ylabel('Height of the structure','FontSize',12);
% title(['mode shape', num2str(i)],'FontSize',18);
% end
Sensitivity matrix
l1 = [1 0 0;0 0 0;0 0 0]; % DIFFERENTIATION OF STIFFNESS WRT K1
l2 = [1 -1 0;-1 1 0;0 0 0]; % DIFFERENTIATION OF STIFFNESS WRT K2
l3 = [0 0 0;0 1 -1;0 -1 1]; % DIFFERENTIATION OF STIFFNESS WRT K3
Mode shape decomposition
a111=0;
a112=(v2'*l1*v1)/(w1-w2);
a113=(v3'*l1*v1)/(w1-w3);
a121=0;
a122=(v2'*l2*v1)/(w1-w2);
a123=(v3'*l2*v1)/(w1-w3);
a131=0;
a132=(v2'*l3*v1)/(w1-w2);
a133=(v3'*l3*v1)/(w1-w3);
a211=(v1'*l1*v2)/(w2-w1);
a212=0;
a213=(v3'*l1*v2)/(w2-w3);
a221=(v1'*l2*v2)/(w2-w1);
a222=0;
a223=(v3'*l2*v2)/(w2-w3);
a231=(v1'*l3*v2)/(w2-w1);
a232=0;
a233=(v3'*l3*v2)/(w2-w3);
a311=(v1'*l1*v3)/(w3-w1);
a312=(v2'*l1*v3)/(w3-w2);
a313=0;
a321=(v1'*l2*v3)/(w3-w1);
a322=(v2'*l2*v3)/(w3-w2);
a323=0;
a331=(v1'*l3*v3)/(w3-w1);
a332=(v2'*l3*v3)/(w3-w2);
a333=0;
SENSITIVITY MATRIX FORMATION (CONSIDERING BOTH EIGEN VALUE AND MODE-SHAPE DECOMPOSITION)
s11 = v1'*l1*v1;
s12 = v1'*l2*v1;
s13 = v1'*l3*v1;
s21=a111*v1+a112*v2+a113*v3;
s22=a121*v1+a122*v2+a123*v3;
s23=a131*v1+a132*v2+a133*v3;
s31 = v2'*l1*v2;
s32 = v2'*l2*v2;
s33 = v2'*l3*v2;
s41=a211*v1+a212*v2+a213*v3;
s42=a221*v1+a222*v2+a223*v3;
s43=a231*v1+a232*v2+a233*v3;
s51 = v3'*l1*v3;
s52 = v3'*l2*v3;
s53 = v3'*l3*v3;
s61=a311*v1+a312*v2+a313*v3;
s62=a321*v1+a322*v2+a323*v3;
s63=a331*v1+a332*v2+a333*v3;
S = [s11 s12 s13;s21 s22 s23;s31 s32 s33;s41 s42 s43;s51 s52 s53;s61 s62 s63]
SENSITIVITY MATRIX FORMATION 2 (CONSIDERING ONLY EIGEN VALUE DECOMPOSITION)
% s11 = v1'*l1*v1;
% s12 = v1'*l2*v1;
% s13 = v1'*l3*v1;
% s21 = v2'*l1*v2;
% s22 = v2'*l2*v2;
% s23 = v2'*l3*v2;
% s31 = v3'*l1*v3;
% s32 = v3'*l2*v3;
% s33 = v3'*l3*v3;
% S=[s11 s12 s13;s21 s22 s23;s31 s32 s33];
UPDATING PARAMETER (CONSIDERING BOTH EIGEN VALUE AND MODE-SHAPE DECOMPOSITION)
K = [k1 k2 k3];
WA = [w1 v1' w2 v2' w3 v3'];
WM = [14.7257 -0.0015 -0.0027 -0.0038 35.4870 0.0029 0.0021 -0.0031 55.2413 -0.0031 0.0033 -0.0011];
K1 = K'+inv(S'*S)*S'*(WM'-WA')
UPDATING PARAMETER 2 (CONSIDERING ONLY EIGEN VALUE DECOMPOSITION)
% K = [k1 k2 k3];
% WA = [w1 w2 w3];
% WM = [13.3221 36.7237 51.7097];
% K1 = K'+inv(S*S')*S'*(WM'-WA')
  2 comentarios
Steven Lord
Steven Lord el 29 de Mzo. de 2022
It's not clear what problem you're experiencing or what question you'd like help with.
Susanta  Barik
Susanta Barik el 29 de Mzo. de 2022
I AM USING SENSITIVITY METHOD TO UPDATE FINITE ELEMENT DATA BUT WHAT THE RESULT I AM GETTING THAT IS NOT EXACT DATA . ALSO WHEN I AM COMPAIRING MY EXPERIMENTAL DATA WITH UPDATED DATA I AM GETTING NOT EXACT DATA

Iniciar sesión para comentar.

Respuestas (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by