using fmincon with extra parameters
조회 수: 19 (최근 30일)
이전 댓글 표시
I am in the process of writing an optimisation code for the design of beams. I am using the fmincon function to achieve the minimal cross sectional area. I have completed the code for 1 iteration, which is shown below:
The objective function:
objective = @(x) (2*x(1)*x(2)) + ((x(3)-(2*x(2)))*x(4)); %Minimise cross sectional area
x0 = [210,17.4,536.7,10.8]; %Initial guess
disp(['Initial Objective: ' num2str(objective(x0))]) %Initial guess C.A
%fmincon inputs
A =[];
b = [];
Aeq = [];
beq = [];
lb = [50; 5; 100; 4];
ub = [450; 70; 1100; 36];
%non linear inequalities and equalities constraints
nonlincon = @nlcon;
%options for the analysis
options = optimoptions('fmincon', 'Algorithm', 'sqp', 'Display', 'iter-detailed',...
'MaxFunctionEvaluations', 1000000, 'MaxIterations', 20000,...
'FunctionTolerance', 1e-1);
%optimisation
x = fmincon(objective, x0, A,b,Aeq,beq,lb,ub,nonlincon,options);
%Display final results with constraints shown
disp(['Final Objective: ' num2str(objective(x))])
[c,ceq] = nlcon(x)
the non-linear constraints:
function [c,ceq] = nlcon(x)
% Initial information
w = 51.1; % N/mm load
L = 4300; % mm Span of beam
E = 210000; % N/mm2 Youngs Modulus
v = 0.3; % Poisson's Ratio
ymo = 1; % Material Partial Factor
...
% Section Classification Flange width to thickness ratio in compression
% Class 2 taken to be the upper limit as same procedure as Class 1
c1 = (((x(1)-x(4))/2)/(x(2)))-(10*eps);
% Section Classification: Web width to thickness ratio in bending
c2 = (x(3)-(2*(x(2))))/(x(4))-(83*eps);
% Bending Moment Constraint
Emx = (w*1.5)*L^2/2; % Nmm Bending Moment of a UDL Cantilever Beam wl2/2
Rmx = (zp*fy)/ymo; % Nmm Bending Resistance in accordance to Eurocode 3 for Class 1 or 2 beams
c3 = Emx-Rmx; % Inequality constraint for Bending Moment
% Shear Force Constraint
Evx = (w*1.5)*L; % N Shear force of a UDL Cantilever beam
Rvx = (Av*(fy*(3^(1/2))))/ymo; % N Shear Resistance in accordance to Eurocode 3 for class 1 or 2
c4 = Evx-Rvx;
% Deflection Constraint
Edf = ((w*L^4)/(8*E*Ix));
Rdf = L/250;
c5 = Edf-Rdf;
...
...
...
...
% Imperfection factor alphaLT
if x(3)/x(1)<=2
alphaLT = 0.49;
else
alphaLT = 0.79;
end
%Non-dimensional Slenderness lambdaLT calculation
lambdaLT = ((zp*fy)/Mcr)^(1/2);
%phiLT
lambdaLT0 = 0.2; %From NA to BS EN 1993-1-1 cl. NA.2.17b for welded sections
BetaLT = 1; % From NA to BS EN 1993-1-1 cl. NA.2.17b for welded sections
phiLT = 0.5*(1+(alphaLT*(lambdaLT-lambdaLT0))+(BetaLT*lambdaLT^2));
%Reduction Factor
XLT0 = 1/(phiLT+(((phiLT^2)-(BetaLT*(phiLT^2)))^(1/2)));
XLT = min([XLT0,1]);
%Buckling Resistance
Mbrd = XLT*zp*fy/ymo;
c6 = Emx-Mbrd;
c7 = (x(3)/x(1))-3.1
%% contraint inequality
c = [c1, c2, c3, c4, c5, c6, c7];
%% constraint equality
ceq=[];
My question is how would I go about doing multiple iteration of this with different Length of beam (L)?
I want to find the optimal design variables for each length of beam between 100:100:15000. I was thinking about using a for loop but not sure how and where to implement this within the code. Any help or advice would be greatly appreciated. Thank you for your time!
댓글 수: 0
채택된 답변
Cris LaPierre
2020년 3월 31일
You'd probably need a loop to iterate through all the lengths. Update your call to fmincon like so
...
L = 4300;
%optimisation
x = fmincon(objective, x0, A,b,Aeq,beq,lb,ub,@(x)nlcon(x,L),options);
...
function [c,ceq] = nlcon(x,L)
...
end
댓글 수: 5
추가 답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Systems of Nonlinear Equations에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!