필터 지우기
필터 지우기

i want some help in my code to do the optimization to the main structure and how to derive the max. response under harmonic excitation with different freq. and optimize it

조회 수: 2 (최근 30일)
clear;clc; warning off;
global Kx Ky Kz m Cx Cy Cz Kxt Kyt Kzt mt Cxt Cyt Czt xg L1x L2x L1y L2y L1z L2z j Beta
wx=2*pi/1.11;
wy=wx;
wz=wx*100;
m=100;
mt=5/100;
Kx=m*wx^2;
Ky=m*wy^2;
Kz=m*wz^2;
Cx=2*0.05*m*wx;
Cy=2*0.05*m*wy;
Cz=2*0.05*m*wz;
sigmax=0.01;
sigmay=0.001;
wxt=4.1/2;
wxt=wx;
wyt=wxt;
wzt=2.9*wxt+sigmay;
%wzt=32^0.5*wxt;
%wxt=3;
%wyt=wxt;
%wyt=4.5;
%mt=0.6;
whar=wxt;
%whar=wyt;
Kxt=mt*wxt^2;
%Kxt=21.6;
%Kyt=21.6/4*2.25;
Kyt=mt*wyt^2;
Kzt=mt*wzt^2;
L1x=0.17;
L2x=0.375-L1x;
L1y=0.17;
L2y=0.375-L1x;
L1z=0.215;
L2z=mt*9.81/Kzt;
%L2z=0.001;
L2=0.01;
L1=1.0;
%wxt=6;
%Lx=1.;
%Lz=1. -mt*9.81/Kyt;
Beta=0.0;
%mt is for TMD
Kisay=0.05;
%Kisay=0.020;
%mt is for TMD
Cxt=2*0.023*(Kxt/mt)^0.5*mt; Cyt=2*0.023*(Kyt/mt)^0.5*mt;Czt=2*0.018*(Kzt/mt)^0.5*mt;
w=1.39;
%--------------------------------------------------------------------------
A=2.1*0.2;
dt=.02;tf=30.;t=0:dt:tf;n=tf/dt;%tsp=time step; tf=final time;
%load('CHICHI0968.mat');%CHICHI0968.mat-max of elcentro=0.3487
%ug=9.81*CHICHI0968(2,1:n+1);%%ELCENTRO_NS0348;A*(wx*w)^2*sin(wx*w*t);
ug=A*sin(wxt*w*t);
xg=ug;
%--------------------------------------------------------------------------
x0N=[0.014 0 0 0 0 0 0 0 0 0 0 0];
x0L=x0N;
xjN(1,:)=x0N;
xjL(1,:)=x0L;
for j=1:n
%for j=1:3500
tint=dt*[j-1 j];%tint=time interval
[tN,xN] = ode45(@nonlinearmodel_structure,tint,x0N);
xjN(j+1,:)=xN(length(tN),:);
x0N=xN(length(tN),:);
%[tL,xL] = ode45(@linearmodel,tint,x0L);
%xjL(j+1,:)=xL(length(tL),:);
%x0L=xL(length(tL),:);
j;
end
%--------------------------------------------------------------------------
ux1=[xjN(:,1)];ux2=[xjN(:,2)];ux3=[xjN(:,3)];
save('Nolin.txt','xjN','-ASCII');
%ux1=[xjL(:,1)];ux2=[xjL(:,2)];
%uy1=[xjN(:,3) xjL(:,3)];uy2=[xjN(:,4) xjL(:,4)];
% clear;clc;
global Kx Ky Kz m Cx Cy Cz Kxt Kyt Kzt mt Cxt Cyt Czt xg L1x L2x L1y L2y L1z L2z j Beta
% Define optimization parameters
nvars = 6; % Number of variables to be optimized
lb = [0.01 0.01 0.01 0.01 , 0.01]; % Lower bounds for the optimization variables
ub = [10 10 10 10 10 10 ]; % Upper bounds for the optimization variables
% Define the GA optimization options
options = optimoptions('ga', 'Display', 'iter', 'PlotFcn' ...
, @gaplotbestf, 'MaxGenerations', 6, 'PopulationSize', 10);
% Define the objective function to be minimized
objFunc = @(x) simulate_structure(x);
% Run the GA optimization
[x_opt, fval_opt] = ga(objFunc, nvars, [], [], [], [], lb, ub, [], options);
% Display the optimized parameters and the corresponding objective function value
disp('Optimized parameters:');
disp(x_opt);
disp('Objective function value:');
disp(fval_opt);
filename = 'Nolin.txt';
formatSpec = '%16f%16f%16f%16f%16f%16f%16f%16f%16f%16f%16f%f%[^\n\r]';
fileID = fopen(filename,'r');
dataArray = textscan(fileID, formatSpec, 'Delimiter', '', 'WhiteSpace', '', 'TextType', 'string', 'ReturnOnError', false);
fclose(fileID);
Nolin = table(dataArray{1:end-1}, 'VariableNames', {'e02','e00','e1','e2','e3','e4','e5','e6','e7','e8','e9','e10'});
clearvars filename formatSpec fileID dataArray ans;
figure
dt = 0.02;
tf = 30.;
t = 0:dt:tf;
subplot(2, 1, 1)
plot(t, Nolin.e2);
subplot(2, 1, 2)
plot(t, x_opt(1,4)*sin(4*w*t/0.6905)/100, t, Nolin.e2);
legend("optimized Model", "Non-linear Model")
function f = simulate_structure(x)
global Kx Ky Kz m Cx Cy Cz Kxt Kyt Kzt mt Cxt Cyt Czt xg L1x L2x L1y L2y L1z L2z j Beta
% Define the system parameters
wx = 2*pi/1.11;
wy = wx;
wz = wx*100;
m = 100;
mt = 5/100;
Kx = m*wx^2;
Ky = m*wy^2;
Kz = m*wz^2;
Cx = 2*0.05*m*wx;
Cy = 2*0.05*m*wy;
Cz = 2*0.05*m*wz;
sigmax = 0.01;
sigmay = 0.001;
wxt = 4.1/2;
wxt = wx;
wyt = wxt;
wzt = 2.9*wxt+sigmay;
Kyt = mt*wyt^2;
Kzt = mt*wzt^2;
L1x = 0.17;
L2x = 0.375-L1x;
L1y = 0.17;
L2y = 0.375-L1x;
L1z = 0.215;
L2z = mt*9.81/Kzt;
L2 = 0.01;
L1 = 1.0;
Beta = 0;
Kisay = 0.05;
Cxt = 2*0.023*(Kxt/mt)^0.5*mt;
Cyt = 2*0.023*(Kyt/mt)^0.5*mt;
Czt = 2*0.018*(Kzt/mt)^0.5*mt;
w = 1.39;
A = 2.1*0.2;
dt = 0.1;
tf = 30.;
t = 0:dt:tf;
n = tf/dt;
ug = A*sin(wxt*w*t);
xg = ug;
% Update the system parameters with the optimized values
Kxt = x(1);
Kyt = x(2);
Kzt = x(3);
L1x = x(4);
L1y = x(5);
% Initialize the simulation
x0N = [0.014 0 0 0 0 0 0 0 0 0 0 0];
x0N(1,:) = x0N;
for j = 1:n
tint = dt*[j-1 j];
[tN,xN] = ode45(@nonlinearmodel_structure,tint,x0N);
xjN(j+1,:) = xN(length(tN),:);
x0N = xN(length(tN),:);
end
% Calculate the objective function value as the maximum displacement of the structure
% from its initial position over the simulation time
f = max(xjN(:,1));
end
function dx=nonlinearmodel_structure(t,x)
global Kx Ky Kz m Cx Cy Cz Kxt Kyt Kzt mt Cxt Cyt Czt xg L1x L2x L1y L2y L1z L2z j Beta
dx = zeros(12,1);
dx(1)=x(7);
dx(2)=x(8);
dx(3)=x(9);
dx(4)=x(10);
dx(5)=x(11);
dx(6)=x(12);
%dx(1)=x(4);
%dx(2)=x(5);
%dx(3)=x(6);
%Vatar_x=((Lx+xx)^2+yy^2)^0.5;
% Vatar_y=((yy+mt*9.81/Kyt+Ly)^2+xx^2)^0.5;
dx(7)=1/m*(-Cx*x(7)-Kx*x(1)+...
-(-Cxt*(x(10)-x(7))...
-Kxt/2*( (x(4)-x(1))+L2x+(x(6)-x(3))^2*0.5*L1x/(L1x+L2x)^2+(x(5)-x(2))^2*0.5*L1x/(L1x+L2x)^2-(x(4)-x(1))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3-(x(4)-x(1))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3)...
+Kxt/2*(-(x(4)-x(1))+L2x+(x(6)-x(3))^2*0.5*L1x/(L1x+L2x)^2+(x(5)-x(2))^2*0.5*L1x/(L1x+L2x)^2+(x(4)-x(1))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3+(x(4)-x(1))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3)...
-Kyt/2*((x(4)-x(1))* (x(5)-x(2)) *L1y/(L1y+L2y)^2+(x(4)-x(1))*L2y/(L1y+L2y)+0.5*(x(4)-x(1))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3- (x(4)-x(1)) *(x(5)-x(2))^2*L1y/(L1y+L2y)^3+0.5*(x(4)-x(1))^3*L1y/(L1y+L2y)^3)...
-Kyt/2*((x(4)-x(1))*(-(x(5)-x(2)))*L1y/(L1y+L2y)^2+(x(4)-x(1))*L2y/(L1y+L2y)+0.5*(x(4)-x(1))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3-(x(4)-x(1))*(-(x(5)-x(2)))^2*L1y/(L1y+L2y)^3+0.5*(x(4)-x(1))^3*L1y/(L1y+L2y)^3)...
-Kzt*((x(4)-x(1))*(x(6)-x(3))*L1z/(L1z+L2z)^2+(x(4)-x(1))*L2z/(L1z+L2z)+0.5*(x(4)-x(1))*(x(5)-x(2))^2*L1z/(L1z+L2z)^3-(x(4)-x(1))*(x(6)-x(3))^2*L1z/(L1z+L2z)^3+0.5*(x(4)-x(1))^3*L1z/(L1z+L2z)^3)))...
-xg(j)*cos(Beta);
dx(8)=1/m*(-Cy*x(8)-Ky*x(2)+...
-(-Cyt*(x(11)-x(8))...
-Kyt/2*( (x(5)-x(2))+L2y+(x(6)-x(3))^2*0.5*L1y/(L1y+L2y)^2+(x(4)-x(1))^2*0.5*L1y/(L1y+L2y)^2-(x(5)-x(2))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3-(x(5)-x(2))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3)...
+Kyt/2*(-(x(5)-x(2))+L2y+(x(6)-x(3))^2*0.5*L1y/(L1y+L2y)^2+(x(4)-x(1))^2*0.5*L1y/(L1y+L2y)^2+(x(5)-x(2))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3+(x(5)-x(2))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3)...
-Kxt/2*((x(5)-x(2))* (x(4)-x(1))*L1x/(L1x+L2x)^2+(x(5)-x(2))*L2x/(L1x+L2x)+0.5*(x(5)-x(2))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3-(x(5)-x(2))* (x(4)-x(1))^2*L1x/(L1x+L2x)^3+0.5*(x(5)-x(2))^3*L1x/(L1x+L2x)^3)...
-Kxt/2*((x(5)-x(2))*(-(x(4)-x(1)))*L1x/(L1x+L2x)^2+(x(5)-x(2))*L2x/(L1x+L2x)+0.5*(x(5)-x(2))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3-(x(5)-x(2))*(-(x(4)-x(1)))^2*L1x/(L1x+L2x)^3+0.5*(x(5)-x(2))^3*L1x/(L1x+L2x)^3)...
-Kzt*((x(5)-x(2))*(x(6)-x(3))*L1z/(L1z+L2z)^2+(x(5)-x(2))*L2z/(L1z+L2z)+0.5*(x(5)-x(2))*(x(4)-x(1))^2*L1z/(L1z+L2z)^3-(x(5)-x(2))*(x(6)-x(3))^2*L1z/(L1z+L2z)^3+0.5*(x(5)-x(2))^3*L1z/(L1z+L2z)^3)))...
-xg(j)*sin(Beta);
dx(9)=1/m*(-Cz*x(9)-Kz*x(3)+...
-(-Czt*(x(12)-x(9))...
-Kzt*((x(6)-x(3))+L2z+(x(4)-x(1))^2*0.5*L1z/(L1z+L2z)^2+(x(5)-x(2))^2*0.5*L1z/(L1z+L2z)^2-(x(6)-x(3))*(x(4)-x(1))^2*L1z/(L1z+L2z)^3-(x(6)-x(3))*(x(5)-x(2))^2*L1z/(L1z+L2z)^3)...
-Kyt/2*((x(6)-x(3))* (x(5)-x(2))*L1y/(L1y+L2y)^2+(x(6)-x(3))*L2y/(L1y+L2y)+0.5*(x(6)-x(3))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3-(x(6)-x(3))* (x(5)-x(2))^2*L1y/(L1y+L2y)^3+0.5*(x(6)-x(3))^3*L1y/(L1y+L2y)^3)...
-Kyt/2*((x(6)-x(3))*(-(x(5)-x(2)))*L1y/(L1y+L2y)^2+(x(6)-x(3))*L2y/(L1y+L2y)+0.5*(x(6)-x(3))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3-(x(6)-x(3))*(-(x(5)-x(2)))^2*L1y/(L1y+L2y)^3+0.5*(x(6)-x(3))^3*L1y/(L1y+L2y)^3)...
-Kxt/2*((x(6)-x(3))* (x(4)-x(1))*L1x/(L1x+L2x)^2+(x(6)-x(3))*L2x/(L1x+L2x)+0.5*(x(6)-x(3))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3-(x(6)-x(3)) *(x(4)-x(1))^2*L1x/(L1x+L2x)^3+0.5*(x(6)-x(3))^3*L1x/(L1x+L2x)^3)...
-Kxt/2*((x(6)-x(3))*(-(x(4)-x(1)))*L1x/(L1x+L2x)^2+(x(6)-x(3))*L2x/(L1x+L2x)+0.5*(x(6)-x(3))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3-(x(6)-x(3))*(-(x(4)-x(1)))^2*L1x/(L1x+L2x)^3+0.5*(x(6)-x(3))^3*L1x/(L1x+L2x)^3)))...
-xg(j)*sin(Beta)*0.+mt/mt*9.81;
dx(10)=1/mt*(-Cxt*(x(10)-x(7))...
-Kxt/2*( (x(4)-x(1))+L2x+(x(6)-x(3))^2*0.5*L1x/(L1x+L2x)^2+(x(5)-x(2))^2*0.5*L1x/(L1x+L2x)^2-(x(4)-x(1))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3-(x(4)-x(1))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3)...
+Kxt/2*(-(x(4)-x(1))+L2x+(x(6)-x(3))^2*0.5*L1x/(L1x+L2x)^2+(x(5)-x(2))^2*0.5*L1x/(L1x+L2x)^2+(x(4)-x(1))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3+(x(4)-x(1))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3)...
-Kyt/2*((x(4)-x(1))* (x(5)-x(2)) *L1y/(L1y+L2y)^2+(x(4)-x(1))*L2y/(L1y+L2y)+0.5*(x(4)-x(1))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3- (x(4)-x(1)) *(x(5)-x(2))^2*L1y/(L1y+L2y)^3+0.5*(x(4)-x(1))^3*L1y/(L1y+L2y)^3)...
-Kyt/2*((x(4)-x(1))*(-(x(5)-x(2)))*L1y/(L1y+L2y)^2+(x(4)-x(1))*L2y/(L1y+L2y)+0.5*(x(4)-x(1))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3-(x(4)-x(1))*(-(x(5)-x(2)))^2*L1y/(L1y+L2y)^3+0.5*(x(4)-x(1))^3*L1y/(L1y+L2y)^3)...
-Kzt*((x(4)-x(1))*(x(6)-x(3))*L1z/(L1z+L2z)^2+(x(4)-x(1))*L2z/(L1z+L2z)+0.5*(x(4)-x(1))*(x(5)-x(2))^2*L1z/(L1z+L2z)^3-(x(4)-x(1))*(x(6)-x(3))^2*L1z/(L1z+L2z)^3+0.5*(x(4)-x(1))^3*L1z/(L1z+L2z)^3))...
-xg(j)*cos(Beta);
dx(11)=1/mt*(-Cyt*(x(11)-x(8))...
-Kyt/2*( (x(5)-x(2))+L2y+(x(6)-x(3))^2*0.5*L1y/(L1y+L2y)^2+(x(4)-x(1))^2*0.5*L1y/(L1y+L2y)^2-(x(5)-x(2))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3-(x(5)-x(2))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3)...
+Kyt/2*(-(x(5)-x(2))+L2y+(x(6)-x(3))^2*0.5*L1y/(L1y+L2y)^2+(x(4)-x(1))^2*0.5*L1y/(L1y+L2y)^2+(x(5)-x(2))*(x(6)-x(3))^2*L1y/(L1y+L2y)^3+(x(5)-x(2))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3)...
-Kxt/2*((x(5)-x(2))* (x(4)-x(1))*L1x/(L1x+L2x)^2+(x(5)-x(2))*L2x/(L1x+L2x)+0.5*(x(5)-x(2))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3-(x(5)-x(2))* (x(4)-x(1))^2*L1x/(L1x+L2x)^3+0.5*(x(5)-x(2))^3*L1x/(L1x+L2x)^3)...
-Kxt/2*((x(5)-x(2))*(-(x(4)-x(1)))*L1x/(L1x+L2x)^2+(x(5)-x(2))*L2x/(L1x+L2x)+0.5*(x(5)-x(2))*(x(6)-x(3))^2*L1x/(L1x+L2x)^3-(x(5)-x(2))*(-(x(4)-x(1)))^2*L1x/(L1x+L2x)^3+0.5*(x(5)-x(2))^3*L1x/(L1x+L2x)^3)...
-Kzt*((x(5)-x(2))*(x(6)-x(3))*L1z/(L1z+L2z)^2+(x(5)-x(2))*L2z/(L1z+L2z)+0.5*(x(5)-x(2))*(x(4)-x(1))^2*L1z/(L1z+L2z)^3-(x(5)-x(2))*(x(6)-x(3))^2*L1z/(L1z+L2z)^3+0.5*(x(5)-x(2))^3*L1z/(L1z+L2z)^3))...
-xg(j)*sin(Beta);
dx(12)=1/mt*(-Czt*(x(12)-x(9))...
-Kzt*((x(6)-x(3))+L2z+(x(4)-x(1))^2*0.5*L1z/(L1z+L2z)^2+(x(5)-x(2))^2*0.5*L1z/(L1z+L2z)^2-(x(6)-x(3))*(x(4)-x(1))^2*L1z/(L1z+L2z)^3-(x(6)-x(3))*(x(5)-x(2))^2*L1z/(L1z+L2z)^3)...
-Kyt/2*((x(6)-x(3))* (x(5)-x(2))*L1y/(L1y+L2y)^2+(x(6)-x(3))*L2y/(L1y+L2y)+0.5*(x(6)-x(3))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3-(x(6)-x(3))* (x(5)-x(2))^2*L1y/(L1y+L2y)^3+0.5*(x(6)-x(3))^3*L1y/(L1y+L2y)^3)...
-Kyt/2*((x(6)-x(3))*(-(x(5)-x(2)))*L1y/(L1y+L2y)^2+(x(6)-x(3))*L2y/(L1y+L2y)+0.5*(x(6)-x(3))*(x(4)-x(1))^2*L1y/(L1y+L2y)^3-(x(6)-x(3))*(-(x(5)-x(2)))^2*L1y/(L1y+L2y)^3+0.5*(x(6)-x(3))^3*L1y/(L1y+L2y)^3)...
-Kxt/2*((x(6)-x(3))* (x(4)-x(1))*L1x/(L1x+L2x)^2+(x(6)-x(3))*L2x/(L1x+L2x)+0.5*(x(6)-x(3))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3-(x(6)-x(3)) *(x(4)-x(1))^2*L1x/(L1x+L2x)^3+0.5*(x(6)-x(3))^3*L1x/(L1x+L2x)^3)...
-Kxt/2*((x(6)-x(3))*(-(x(4)-x(1)))*L1x/(L1x+L2x)^2+(x(6)-x(3))*L2x/(L1x+L2x)+0.5*(x(6)-x(3))*(x(5)-x(2))^2*L1x/(L1x+L2x)^3-(x(6)-x(3))*(-(x(4)-x(1)))^2*L1x/(L1x+L2x)^3+0.5*(x(6)-x(3))^3*L1x/(L1x+L2x)^3))...
-xg(j)*sin(Beta)*0.+mt/mt*9.81;
end
  댓글 수: 2
Star Strider
Star Strider 2023년 5월 26일
Please go into some detail as to what you are doing.
If you are optimising something, what is it?
If you are fitting data, which one of your differential equations are you fitting it to?
LD = load('CHICHI0968.mat');
Data = LD.CHICHI0968;
t = Data(1,:);
s = Data(2,:);
figure
plot(t, s)
grid
xlabel('Time')
ylabel('Signal')
It is best not to use global variables. Instead use the approach described in Passing Extra Parameters.
.
Doha Ali
Doha Ali 2023년 5월 26일
I just want to put the structure under harmonic excitation then optimize the response of the main structure .

댓글을 달려면 로그인하십시오.

답변 (0개)

카테고리

Help CenterFile Exchange에서 Assembly에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by