필터 지우기
필터 지우기

I need to know the value of parameter of this code and do the optimization for x,y,z ? can you the parameter in z direction and the properties of MODEL TMD?

조회 수: 1 (최근 30일)

clear
clc
%A = eye(4);
A=eye(5);
%b = [80; 320; 50;50];
%UB=[80; 320; 50;50];

b = [80*2; 290*2; 200;200;200];
UB=[80*2; 290*2; 200;200;200];

b = [65*2/2/2; 570/2/2; 200;200;220];
UB=[65*2/2/2; 570/2/2; 200;200;220];

%b = [107; 453; 94;163;200];
%UB=[107; 453; 94;163;200];
%lb = [50;200;30;30];
%LB=[50;200;30;30];

lb = [50;220;50;80;50];
LB=[50;220;50;80;50];

lb = [48*2/2/2;215*2/2/2;30;50;20];
LB=[48*2/2/2;215*2/2/2;30;50;20];

%lb = [107;453;94;163;25];
%LB=[107;453;94;163;25];
%pause(5)
%[x,fval,exitflag] = ga(@TMD_nonlineardamploop,4,A,b,[],[],lb)
[x,fval,exitflag] = ga(@TMD_nonlineardamploop,5,A,b,[],[],LB,UB,[],[1,2,3,4,5])
function y = TMD_nonlineardamploop(x)
if ~isvector(x)
error('Input must be a vector')
end

Kxt=x(1);
Kyt=x(2);
Kisayx=x(3)/1000;
Kisayy=x(4)/1000;
Kxt
Kyt
Kisayx
Kisayy
Ly=x(5)/100
m=30.;
Kx=1320; Ky=1320*90.0;
Kx=1320; Ky=1320*1.00;

wx=(Kx/m)^0.5;
%Kyt=13.5;
mt=1.5*2/2/2;

Lx=1.;
%Ly=1.;
Cx=2*0.025*4*(Kx/m)^0.5*m;
Cy=2*0.025*4*(Ky/m)^0.5*m;
%Kisayy=0.065;
%mt is for TMD
Cxt=2*Kisayx*(Kxt/mt)^0.5*mt; Cyt=2*Kisayy*(Kyt/mt)^0.5*mt;
%%stiffness & mass & damping
Beta=0.0;
%6063 7310 6160 6456
%309+15
%1270642
M=[m 0 0 0;0 m 0 0;0 0 mt 0;0 0 0 mt];
K=[Kx+Kxt 0 -Kxt 0;0 Ky+Kyt 0 -Kyt; -Kxt 0 Kxt 0; 0 -Kyt 0 Kyt];
C=[Cx+Cxt 0 -Cxt 0;0 Cy+Cyt 0 -Cyt; -Cxt 0 Cxt 0; 0 -Cyt 0 Cyt];

%eig(K*M^(-1))
wy=(Ky/m)^0.5;
%AAa=0.99;
www=1;
for F_wy=0.65*wx:0.1/2:1.25*wx
for jjjj=1:2000.
time=jjjj*0.02;

AAa=0.65;

%U_dotdotg(jjjj)=AAa*sin(1.*wy*time);
U_dotdotg(jjjj)=AAa*sin(F_wy*time);
end
g=1;
dt=0.02;

%%---------------------------------------------
%these parameters are for inital displacement, acceleration, velocity
x=zeros(length(U_dotdotg),4);
v=zeros(length(U_dotdotg),4);
a=zeros(length(U_dotdotg),4);
%time
time=0:dt:dt*(length(U_dotdotg)-1);
%these parameters are used for Newmark integration
%aa = K*dt^2/2+C*dt;
%bb = K*dt;
%Keff=M+0.5*dt*C+dt^2/4*K;

length(time);
length(U_dotdotg);

for i=1:length(U_dotdotg)-1

i ;
%initial step calcultion
disp(1,:)=x(i,:)+dt*v(i,:)+0.5*dt^2*a(i,:);
velocity(1,:)=v(i,:)+dt*a(i,:);
accel(1,:)=a(i,:);
%initial guess of displacement, velocity, acceleration
DD=disp(i,:)+dt*velocity(i,:)+dt^2/2*accel(i,:);
VV=velocity(i,:)+dt*accel(i,:);
AA=accel(i,:);
criterion=10;
j=1;
%%ccc(1,:)=[0,0];
Keff=M+0.5*dt*C+dt^2/4*K;
eig(K^(-1)*M);
while criterion>0.00005
for ii=1:4
xx=DD(3)-DD(1) ;
yy=DD(4)-DD(2);
VVx=VV(3)-VV(1);
VVy=VV(4)-VV(2);
fx1=Kxt*( Lx+xx -(Lx^2+Lx*xx)/ ((Lx+xx )^2+yy ^2)^0.5)-Kxt*xx ; fx2=Kyt*(xx- (Ly*xx) /((Ly+yy) ^2+xx^2 )^0.5) ;

fxdamp1=Cxt*(VVx*(-yy^2/Lx^2+2*xx*yy^2/Lx^3)+VVy*(yy/Lx-xx*yy/Lx^2-yy^3/Lx^3+yy*xx^2/Lx^3));
fxdamp2=Cyt*(VVx*(+xx^2/Ly^2-2*xx^2*yy/Ly^3)+VVy*(xx/Ly-xx*yy/Ly^2-xx^3/Lx^3+yy^2*xx/Ly^3));

%%fXdamp1=Cxt*(VVx*(Lx+xx)^2+VVy*yy*(Lx+xx))/((Lx+xx)^2+yy^2)-Cxt*VVx;
%%fXdamp2=Cyt*(VVy*(Ly+yy)*xx+VVx*xx^2)/((Ly+yy)^2+xx^2);
%fx1=Kxt*( xx +yy^2/2/Lx -xx*yy^2/Lx^2)-Kxt*xx ; fx2=Kyt*(xx*yy/Ly +xx^3/2/Ly^2 -xx*yy^2/Ly^2) ;

fy1=Kxt*(yy -(Lx*yy) /((Lx+xx)^2+yy^2 )^0.5) ; fy2= Kyt*( Ly+yy -(Ly^2+Ly*yy) /((Ly+yy) ^2+ xx^2)^0.5)-Kyt*yy ;

fydamp1=Cxt*(VVy*(+yy^2/Lx^2-2*yy^2*xx/Lx^3)+VVx*(yy/Lx-yy*xx/Lx^2-yy^3/Ly^3+xx^2*yy/Lx^3));
fydamp2=Cyt*(VVy*(-xx^2/Ly^2+2*yy*xx^2/Ly^3)+VVx*(xx/Ly-yy*xx/Ly^2-xx^3/Ly^3+xx*yy^2/Ly^3));

%%fYdamp1=Cxt*(VVx*(Lx+xx)*yy+VVy*yy^2)/((Lx+xx)^2+yy^2);
fYdamp2=Cyt*(VVy*(Ly+yy)^2+VVx*xx*(Ly+yy))/((Ly+yy)^2+xx^2)-Cyt*VVy;
%fy1= Kxt*(xx*yy/Lx +yy^3/2/Lx^2 -yy*xx^2/Lx^2) ;fy2=Kyt*( yy +xx^2/2/Ly -yy*xx^2/Ly^2)-Kyt*yy;
if ii==1
F_N_S(ii)=-fx1-fx2+(-fxdamp1-fxdamp2);
Fe(i,ii)=-U_dotdotg(i)*g*cos(Beta);
Fe(i+1,ii)=-U_dotdotg(i+1)*g*cos(Beta);
elseif ii==2
F_N_S(ii)=-fy1-fy2+(-fydamp1-fydamp2);
Fe(i,ii)=-U_dotdotg(i)*g*sin(Beta);
Fe(i+1,ii)=-U_dotdotg(i+1)*g*sin(Beta);
elseif ii==3
F_N_S(ii)=fx1+fx2+(fxdamp1+fxdamp2);
Fe(i,ii)=-U_dotdotg(i)*g*cos(Beta);
Fe(i+1,ii)=-U_dotdotg(i+1)*g*cos(Beta);
elseif ii==4
F_N_S(ii)=fy1+fy2+(fydamp1+fydamp2);
Fe(i,ii)=-U_dotdotg(i)*g*sin(Beta);
Fe(i+1,ii)=-U_dotdotg(i+1)*g*sin(Beta);
end
end
%ccc(i+1,:)=[F_N_S(1)+Kxt*xx,F_N_S(2)+Kyt*yy];
%ccc(i+1,:)=[Kxt*( Lx+xx -(Lx^2+Lx*xx)/ ((Lx+xx )^2+yy ^2)^0.5),Kyt*(xx- (Ly*xx) /((Ly+yy) ^2+xx^2 )^0.5)];
%ccc(i+1,:)=[Kxt*(yy -(Lx*yy) /((Lx+xx)^2+yy^2 )^0.5),Kyt*( Ly+yy -(Ly^2+Ly*yy) /((Ly+yy) ^2+ xx^2)^0.5)];
%%ccc(i+1,:)=[fydamp2,fYdamp2];
Fe(i+1,:)=Fe(i+1,:)*M;

deltaPhat=-AA*M +Fe(i+1,:) -VV*C -DD*K -F_N_S;

deltaA=deltaPhat*Keff^(-1);
DD=DD+dt^2*0.25*deltaA;
VV=VV+dt*0.5*deltaA;
AA=AA+deltaA;
criterion=max(abs(deltaA));
j=j+1;
if j==95
criterion=0.0;
end
NORM1(i)=max(abs(deltaPhat));
end

j ;
iterationnum(i)=j;
accel(i+1,:)=AA;
velocity(i+1,:)=VV;
disp(i+1,:)=DD;

end

MMM(www,:)=[max(abs(disp(1500:1999,1))),max(abs(disp(1500:1999,2)))];
C_MMM(www)=www;
C_MMM1(www)=F_wy;

www=www+1;
if mod(www,10)==0
www;
end
%Axt=max(abs(disp(1500:1999,3)-disp(1500:1999,1)));
%Ayt=max(abs(disp(1500:1999,4)-disp(1500:1999,2)));
% AMP(www,:)=[max(abs(disp(1500:1999,1))),max(abs(disp(1500:1999,2))) ,Axt ,Ayt];
end
%plot(C_MMM1,[MMM(:,1),MMM(:,2)]),grid
%load 'disp_linear.txt'
non_linansw=[C_MMM1',MMM(:,1),MMM(:,2)];
%save('non_linansw_1.txt','non_linansw','-ASCII');

%load 'linansw.txt'
%plot(linansw(:,1)/6.63,[linansw(:,2),non_linansw(:,2),linansw(:,3),non_linansw(:,3)]),grid
%plot(AMP(:,1))
%plot(disp(2000:2500,3)-disp(2000:2500,1))
%plot(linansw(:,1)/6.63,[non_linansw(:,2),non_linansw(:,3)]),grid
%plot(time(:),[disp(:,2)-disp(:,4)]),grid
%plot(time(:),[ccc(:,1),ccc(:,2)]),grid
%plot(time(:),[ccc(:,1)]),grid
y =max( max([non_linansw(:,2),non_linansw(:,3)]));
y
end

답변 (0개)

카테고리

Help CenterFile Exchange에서 Get Started with Optimization Toolbox에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by