Weight Optimization of 10-bar truss structure with MATLAB

To demonstrate the application of OpenSees.NET, we optimize weight the 10-bar truss structure that study by several researchers before [?]. MATLAB software capability to interact with .NET assemblies. In this problem, the fitness function is the weight and constraints of the optimization problem are depend on stress and the minimum cross-section of elements. Figure () shows the properties of the 10-bar truss structure. Figure () shows the input codes of the model in MATLAB software.

MATLAB Optimization Toolbox 

optimtool opens the Optimization app. Use the Optimization app to select a solver, optimization options, and run problems. See Optimization App for a complete description of the Optimization app. The Optimization app can be used to run any Optimization Toolbox™ solver except intlinprog, and any Global Optimization Toolbox solver except GlobalSearch and MultiStart. Results can be exported to a file or to the MATLAB® workspace as a structure.

optimtool(optstruct) starts the Optimization app and loads optstructoptstruct can be either optimization options or an optimization problem structure. Create optimization options with the optimoptions or optimset function, or by using the export option from the Optimization app. Create a problem structure by exporting the problem from the Optimization app to the MATLAB workspace. If you have Global Optimization Toolbox, you can create a problem structure for fminconfminunclsqnonlin, or lsqcurvefit using the createOptimProblem function.

optimtool('solver') starts the Optimization app with the specified solver, identified as a character vector, and the corresponding default options and problem fields. All Optimization Toolbox and Global Optimization Toolbox solvers are valid inputs to the optimtool function, except for intlinprogGlobalSearch, and MultiStart.    more info

MATLAB Optimization toolbox

Fig (1) - optimtool-MATLAB Optimization Toolbox

 

Building Problem (Fitness & Constraint Functions)

In this example, the truss model is constructed in OpenSees.NET. Structural specifications are shown in Fig. (2). The minimum cross-sectional area is 0.1 in2. The maximum stress of the members are 25 ksi, but the member 9 has a maximum stress of 75ksi. The length of the horizontal and vertical members are 360 inches. The force on nodes 2 and 4 of the structure is 100Kips.

10-bar truss structures

Fig (2) - 10-bar truss structure

Fitness Function 

Fitness function (fitness.m) calculates the weight of the structure.


function weight = fitness(A)
l=360;
r= 0.1;
lp = l*sqrt(2);
L = [l,l,l,l,l,l,lp,lp,lp,lp];
weight =0;
for i=1:10
weight = weight + L(i)*A(i)*r;
end
end
    

Constraint Function

Constraints define in constraints function (constraints.m) and contain 20 non-equal constraints for limit minimum section area and maximum stress in each element.


function [ c,ceq ] = constraints(A)
l=360;
e=3000.0;
p=100;
ceq = [];
NET.addAssembly('System');
NET.addAssembly('System.IO');
NET.addAssembly(strcat(pwd,'\OpenSees.NET.X64.dll'));
%domain
theDomain = OpenSees.Components.DomainWrapper();

% nodes
node6 = OpenSees.Components.NodeWrapper(6,2,0,0);
node4 = OpenSees.Components.NodeWrapper(4,2,l*1,0);
node2 = OpenSees.Components.NodeWrapper(2,2,l*2,0);
node5 = OpenSees.Components.NodeWrapper(5,2,0,l);
node3 = OpenSees.Components.NodeWrapper(3,2,l*1,l);
node1 = OpenSees.Components.NodeWrapper(1,2,l*2,l);
theDomain.AddNode(node1);
theDomain.AddNode(node2);
theDomain.AddNode(node3);
theDomain.AddNode(node4);
theDomain.AddNode(node5);
theDomain.AddNode(node6);

%material
theMaterial = OpenSees.Materials.Uniaxials.ElasticMaterialWrapper(1, e,0);


%elements
ndm = 2;
truss1 = OpenSees.Elements.TrussWrapper(1, ndm, 5, 3, theMaterial, A(1), 0, 0, 0);
truss2 = OpenSees.Elements.TrussWrapper(2, ndm, 3, 1, theMaterial, A(2), 0, 0, 0);
truss3 = OpenSees.Elements.TrussWrapper(3, ndm, 6, 4, theMaterial, A(3), 0, 0, 0);
truss4 = OpenSees.Elements.TrussWrapper(4, ndm, 4, 2, theMaterial, A(4), 0, 0, 0);
truss5 = OpenSees.Elements.TrussWrapper(5, ndm, 3, 4, theMaterial, A(5), 0, 0, 0);
truss6 = OpenSees.Elements.TrussWrapper(6, ndm, 1, 2, theMaterial, A(6), 0, 0, 0);
truss7 = OpenSees.Elements.TrussWrapper(7, ndm, 5, 4, theMaterial, A(7), 0, 0, 0);
truss8 = OpenSees.Elements.TrussWrapper(8, ndm, 6, 3, theMaterial, A(8), 0, 0, 0);
truss9 = OpenSees.Elements.TrussWrapper(9, ndm, 3, 2, theMaterial, A(9), 0, 0, 0);
truss10 = OpenSees.Elements.TrussWrapper(10, ndm, 4, 1, theMaterial, A(10), 0, 0, 0);
theDomain.AddElement(truss1);
theDomain.AddElement(truss2);
theDomain.AddElement(truss3);
theDomain.AddElement(truss4);
theDomain.AddElement(truss5);
theDomain.AddElement(truss6);
theDomain.AddElement(truss7);
theDomain.AddElement(truss8);
theDomain.AddElement(truss9);
theDomain.AddElement(truss10);

%restraints
sp1 =  OpenSees.Components.Constraints.SP_ConstraintWrapper(6, 0, 0.0, true);
sp2 =  OpenSees.Components.Constraints.SP_ConstraintWrapper(6, 1, 0.0, true);
sp3 =  OpenSees.Components.Constraints.SP_ConstraintWrapper(5, 0, 0.0, true);
sp4 =  OpenSees.Components.Constraints.SP_ConstraintWrapper(5, 1, 0.0, true);
theDomain.AddSP_Constraint(sp1);
theDomain.AddSP_Constraint(sp2);
theDomain.AddSP_Constraint(sp3);
theDomain.AddSP_Constraint(sp4);

%loading
theSeries = OpenSees.Components.Timeseries.LinearSeriesWrapper();
theLoadPattern = OpenSees.Components.LoadPatterns.LoadPatternWrapper(1);
theLoadPattern.SetTimeSeries(theSeries);
theDomain.AddLoadPattern(theLoadPattern);

theLoadValues = OpenSees.VectorWrapper(2);
theLoadValues.Set(0, 0);
theLoadValues.Set(1, -p);
nl2 = OpenSees.Components.Loads.NodalLoadWrapper(2, 2, theLoadValues, false);
nl4 = OpenSees.Components.Loads.NodalLoadWrapper(1, 4, theLoadValues, false);
theDomain.AddNodalLoad(nl2, 1);
theDomain.AddNodalLoad(nl4, 1);

%setup analysis
theModel =  OpenSees.AnalysisModelWrapper();
theSolnAlgo =  OpenSees.Algorithms.LinearWrapper();
theIntegrator =  OpenSees.Integrators.Static.LoadControlWrapper(1.0, 1, 1.0, 1.0);
theHandler =  OpenSees.Handlers.PlainHandlerWrapper();
theRCM =  OpenSees.GraphNumberers.RCMWrapper(false);
theNumberer =  OpenSees.Numberers.DOF_NumbererWrapper(theRCM);
theSolver =  OpenSees.Systems.Linears.BandSPDLinLapackSolverWrapper();
theSOE =  OpenSees.Systems.Linears.BandSPDLinSOEWrapper(theSolver);
theTest = OpenSees.ConvergenceTests.CTestNormDispIncrWrapper(1e-8, 6, 2, 2, 1.0e10);
theAnalysis =  OpenSees.Analysis.StaticAnalysisWrapper(...
theDomain,...
theHandler,...
theNumberer,...
theModel,...
theSolnAlgo,...
theSOE,...
theIntegrator,...
theTest...
);
ret = theAnalysis.Analyze(1);
nodedisp = node2.GetDisp();
disp(ret);
stressId = 40;
x = truss1.GetDoubleResponse(stressId);
s1 = abs(double(truss1.GetDoubleResponse(stressId)))-25;
s2 = abs(double(truss2.GetDoubleResponse(stressId)))-25;
s3 = abs(double(truss3.GetDoubleResponse(stressId)))-25;
s4 = abs(double(truss4.GetDoubleResponse(stressId)))-25;
s5 = abs(double(truss5.GetDoubleResponse(stressId)))-25;
s6 = abs(double(truss6.GetDoubleResponse(stressId)))-25;
s7 = abs(double(truss7.GetDoubleResponse(stressId)))-25;
s8 = abs(double(truss8.GetDoubleResponse(stressId)))-25;
s9 = abs(double(truss9.GetDoubleResponse(stressId)))-75;
s10 = abs(double(truss10.GetDoubleResponse(stressId)))-25;
constresses = [s1,s2,s3,s4,s5,s6,s7,s8,s9,s10];
conareas = zeros(10,1);
for i=1:10
conareas(i) = 0.1 - A(i);
end
c=zeros(20,1);
c(1:10) = constresses(1:10);
c(11:20) = conareas(1:10);

end


    

Results

In MATLAB software, the optimization of a nonlinear constrained function is performed using the fmincon function. This function finds minimum of a nonlinear multivariate bounded function. For more information on the theory and description, go to the following address. https://www.mathworks.com/help/optim/ug/fmincon.html