Azzera filtri
Azzera filtri

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 visualizzazione (ultimi 30 giorni)

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

Risposte (0)

Categorie

Scopri di più su Get Started with Optimization Toolbox 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!

Translated by