solutions for every iteration

3 visualizzazioni (ultimi 30 giorni)
sura naji
sura naji il 23 Gen 2022
Risposto: Image Analyst il 23 Gen 2022
format long
clear;
clc;
fun=@testcost
lb=[0.1,0.0003];
ub=[1.0,0.01];
%lb=[ 1,1,1,1,1,1]; %lower bound of variables
%ub=[10,10,10,10,10,10]; % upper bound of variables
nvars=2 % number rof variables
options = optimoptions('particleswarm')
options.SwarmSize=2000
options.Display='iter'
options = optimoptions(options,'PlotFcn',@pswplotbestf);
[x,fval,exitflag,output,MaxIterations] = particleswarm(fun,nvars,lb,ub,options);
%As=pi*(x(1)/2)^2-pi*((x(1)-2*x(2))/2)^2;
%Ac=pi*((x(1)-2*x(2))/2)^2;
%h=3.0;
%csteel=0.75*7850;
%cconc=180;
%cost=As*h*csteel+Ac*h*cconc
%****************************************************************
Ninp=6000;
Minp=100;
%****************************************************************
clc
[penalty,Mcalc,cost,x,Nmax,lambda, slenderness] = testcostbest(x,Ninp,Minp)
%diameter=x(1)
%thickness=x(2)
h=3;
fc=35;
fy=360;
Es=200;
Ec=((fc+8)/10)^0.3*22;
%-----------------calculations for global buckling ---------------
Asb=(pi*(x(1)/2)^2-pi*((x(1)-2*x(2))/2)^2);
Acb=(pi*((x(1)-2*x(2))/2)^2);
Npl=(fy*Asb+fc*Acb)*1000;
moicb=(1/4*pi*((x(1)-2*x(2))/2)^4);
moisb=((1/4*pi*(x(1)/2)^4)-(1/4*pi*((x(1)-2*x(2))/2)^4));
k=0.7;
EIeff=Es*moisb+0.6*Ec*moicb;
Ncr=pi^2*EIeff/(k*h)^2*1000000;
lambda=sqrt(Npl/Ncr);
%--------------slenderness-----------------------------------------
slenderness=x(1)/x(2);
% objective function
function f = testcost(x,Ninp,Minp)
% input axşal load (kN) and bending moment(kN.m)
%***********************************************************
Ninp=6000;
Minp=100;
%**********************************************************
% costs of steel and concrete
csteel=0.75*7850;
cconc=180;
% cross-sectional areas (mm2), moment of inertias (mm4), section modulus (mm3)
As=pi*(x(1)/2)^2-pi*((x(1)-2*x(2))/2)^2;
Ac=pi*((x(1)-2*x(2))/2)^2;
%moic=1/4*pi*((x(1)-2*x(2))/2)^4;
%mois=(1/4*pi*(x(1)/2)^4)-(1/4*pi*((x(1)-2*x(2))/2)^4);
sm=pi/32*(x(1)-2*x(2))^3;
% material properties: yield streses (MPa), leastic modulus (GPA) and height(m)
fc=35;
fy=360;
Es=200;
Ec=((fc+8)/10)^0.3*22;
h=3.0;
% other variables and constants
a1=817.2;
a2=5528;
a3=128.3;
a4=1282;
b1=595.2;
b2=5924;
b3=275.5;
b4=2103;
c1=22.5;
c2=-6358;
c3=6.827;
c4=-2257;
d1=385.7;
d2=1961;
d3=108.4;
d4=213.9;
teta=As*fy/(Ac*fc);
Nmax=Ac*fc*(a1+a2*x(2)/x(1)+a3*fc/fy+a4*teta);
Mmax=sm*fc*(b1+b2*x(2)/x(1)+b3*fc/fy+b4*teta);
M0=sm*fc*(c1+c2*x(2)/x(1)+c3*fc/fy-c4*teta);
NMmax=Ac*fc*(d1+d2*x(2)/x(1)+d3*fc/fy+d4*teta);
% coefficients of the polynomial
k1=M0;
k2=((2*Nmax^4-4*Nmax^2*NMmax^2)*(Mmax-M0)-2*M0*NMmax^2)/(Nmax*NMmax*(2*NMmax^3-3*Nmax*NMmax^2+Nmax^3));
k3=((4*Nmax*NMmax^3-Nmax^4)*(Mmax-M0)+3*M0*NMmax^4)/(Nmax*NMmax^2*(2*NMmax^3-3*Nmax*NMmax^2+Nmax^3));
k4=((Nmax^2-2*Nmax*NMmax)*(Mmax-M0)-M0*NMmax^2)/(Nmax*NMmax^2*(2*NMmax^3-3*Nmax*NMmax^2+Nmax^3));
Mcalc=k1+k2*Ninp+k3*Ninp^2+k4*Ninp^4;
%
% ---------penalties-------------
penalty=0;
%--------penalty #1 ------------(d-2*t>0)--------
cc1=x(1)-2*x(2);
if cc1<0
penalty=abs(cc1);
end
%-------penalty #2 -----------
cc2=x(1)/x(2);
if cc2>90*235/fy
penalty=penalty+100*abs(cc2);
end
slenderness=x(1)/x(2);
%--------penalty #3 ----------
if Minp>Mcalc
cc3=abs(Minp-Mcalc);
penalty=penalty+cc3*10;
end
%--------penalty #4 ----------
Asb=(pi*(x(1)/2)^2-pi*((x(1)-2*x(2))/2)^2);
Acb=(pi*((x(1)-2*x(2))/2)^2);
Npl=(fy*Asb+fc*Acb)*1000;
moicb=(1/4*pi*((x(1)-2*x(2))/2)^4);
moisb=((1/4*pi*(x(1)/2)^4)-(1/4*pi*((x(1)-2*x(2))/2)^4));
k=0.7;
EIeff=Es*moisb+0.6*Ec*moicb;
Ncr=pi^2*EIeff/(k*h)^2*1000000;
lambda=sqrt(Npl/Ncr);
cc4=lambda-2;
if cc4>0
penalty=penalty+abs(cc4)*1000;
end
%--------------------------------
% objective function
cost=As*h*csteel+Ac*h*cconc;
%cost1=abs(Minp-Mcalc)
f=cost*(1+penalty*100000);
if As<0
f=1000000000;
end
if Nmax<Ninp
f=1000000000;
penalty=f;
end
penalty;
end
% objective function
function [penalty,Mcalc,cost,x,Nmax, lambda, slenderness] = testcostbest(x, Ninp, Minp)
% costs of steel and concrete
csteel=0.75*7850;
cconc=180;
% cross-sectional areas (mm2), moment of inertias (mm4), section modulus (mm3)
As=pi*(x(1)/2)^2-pi*((x(1)-2*x(2))/2)^2;
Ac=pi*((x(1)-2*x(2))/2)^2;
%moic=1/4*pi*((x(1)-2*x(2))/2)^4;
%mois=(1/4*pi*(x(1)/2)^4)-(1/4*pi*((x(1)-2*x(2))/2)^4);
sm=pi/32*(x(1)-2*x(2))^3;
% material properties: yield streses (MPa), leastic modulus (GPA) and height(m)
fc=35;
fy=360;
Es=200;
Ec=((fc+8)/10)^0.3*22;
h=3.0;
% other variables and constants
a1=817.2;
a2=5528;
a3=128.3;
a4=1282;
b1=595.2;
b2=5924;
b3=275.5;
b4=2103;
c1=22.5;
c2=-6358;
c3=6.827;
c4=-2257;
d1=385.7;
d2=1961;
d3=108.4;
d4=213.9;
teta=As*fy/(Ac*fc);
Nmax=Ac*fc*(a1+a2*x(2)/x(1)+a3*fc/fy+a4*teta);
Mmax=sm*fc*(b1+b2*x(2)/x(1)+b3*fc/fy+b4*teta);
M0=sm*fc*(c1+c2*x(2)/x(1)+c3*fc/fy-c4*teta);
NMmax=Ac*fc*(d1+d2*x(2)/x(1)+d3*fc/fy+d4*teta);
% coefficients of the polynomial
k1=M0;
k2=((2*Nmax^4-4*Nmax^2*NMmax^2)*(Mmax-M0)-2*M0*NMmax^2)/(Nmax*NMmax*(2*NMmax^3-3*Nmax*NMmax^2+Nmax^3));
k3=((4*Nmax*NMmax^3-Nmax^4)*(Mmax-M0)+3*M0*NMmax^4)/(Nmax*NMmax^2*(2*NMmax^3-3*Nmax*NMmax^2+Nmax^3));
k4=((Nmax^2-2*Nmax*NMmax)*(Mmax-M0)-M0*NMmax^2)/(Nmax*NMmax^2*(2*NMmax^3-3*Nmax*NMmax^2+Nmax^3));
Mcalc=k1+k2*Ninp+k3*Ninp^2+k4*Ninp^4;
%
% ---------penalties-------------
penalty=0;
%--------penalty #1 ------------(d-2*t>0)--------
cc1=x(1)-2*x(2);
if cc1<0
penalty=abs(cc1);
end
%-------penalty #2 -----------
cc2=x(1)/x(2);
if cc2>90*235/fy
penalty=penalty+100*abs(cc2);
end
slenderness=x(1)/x(2);
%--------penalty #3 ----------
if Minp>Mcalc
cc3=abs(Minp-Mcalc);
penalty=penalty+cc3*10;
end
%--------penalty #4 ----------
Asb=(pi*(x(1)/2)^2-pi*((x(1)-2*x(2))/2)^2);
Acb=(pi*((x(1)-2*x(2))/2)^2);
Npl=(fy*Asb+fc*Acb)*1000;
moicb=(1/4*pi*((x(1)-2*x(2))/2)^4);
moisb=((1/4*pi*(x(1)/2)^4)-(1/4*pi*((x(1)-2*x(2))/2)^4));
k=0.7;
EIeff=Es*moisb+0.6*Ec*moicb;
Ncr=pi^2*EIeff/(k*h)^2*1000000;
lambda=sqrt(Npl/Ncr);
cc4=lambda-2;
if cc4>0
penalty=penalty+abs(cc4)*1000;
end
%--------------------------------
% objective function
cost=As*h*csteel+Ac*h*cconc;
%cost1=abs(Minp-Mcalc)
f=cost*(1+penalty*100000);
if As<0
f=1000000000;
end
if Nmax<Ninp
f=1000000000;
penalty=f;
end
penalty;
end
Can anyone help me please, when l do run l get just the optimum solutions, l need the cost and valued x1 and x2 for evey iteration ?

Accedi per rispondere a questa domanda.

Risposte (1)

Image Analyst
Image Analyst il 23 Gen 2022
You need to index your variables to make them vectors (lists of values) rather than just overwriting them every iterations, like instead of
myVariable = whatever;
do
myVariable(index) = whatever;
where index is your "for" or "while" loop iterator variable that gets incremented on each iteration.

Categorie

Scopri di più su Whos in Help Center e File Exchange

Tag

Community Treasure Hunt

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

Start Hunting!

Translated by