www.gusucode.com > globaloptim 案例源码程序 matlab代码 > globaloptim/steppedCantileverExample.m
%% Solving a Mixed Integer Engineering Design Problem Using the Genetic Algorithm % % This example shows how to solve a mixed integer engineering design % problem using the Genetic Algorithm (|ga|) solver in Global % Optimization Toolbox. % % The problem illustrated in this example involves the design of a stepped % cantilever beam. In particular, the beam must be able to carry a % prescribed end load. We will solve a problem to minimize the beam volume % subject to various engineering design constraints. % % In this example we will solve two bounded versions of the problem % published in [1]. % Copyright 2012-2015 The MathWorks, Inc. %% Stepped Cantilever Beam Design Problem % % A stepped cantilever beam is supported at one end and a load is applied % at the free end, as shown in the figure below. The beam must be able to % support the given load, $P$, at a fixed distance $L$ from the support. % Designers of the beam can vary the width ($b_i$) and height ($h_i$) % of each section. We will assume that each section of the cantilever has % the same length, $l$. % %% % <<../steppedCantileverFigure.png>> % %% % % _Volume of the beam_ % % The volume of the beam, $V$, is the sum of the volume of the individual % sections % % $$V = l(b_1h_1 + b_2h_2 + b_3h_3 + b_4h_4 + b_5h_5)$$ % %% % % _Constraints on the Design : 1 - Bending Stress_ % % Consider a single cantilever beam, with the centre of coordinates at the % centre of its cross section at the free end of the beam. The bending % stress at a point $(x, y, z)$ in the beam is given by the following % equation % % $$\sigma_b = M(x)y/I$$ % % where $M(x)$ is the bending moment at $x$, $x$ is the distance % from the end load and $I$ is the area moment of inertia of the beam. % % Now, in the stepped cantilever beam shown in the figure, the maximum % moment of each section of the beam is $PD_i$, where $D_i$ is the maximum % distance from the end load, $P$, for each section of the beam. Therefore, % the maximum stress for the $i$-th section of the beam, $\sigma_i$, is % given by % % $$\sigma_i = PD_i(h_i/2)/I_i$$ % % where the maximum stress occurs at the edge of the beam, $y = h_i/2$. % The area moment of inertia of the $i$-th section of the beam is given by % % $$I_i = b_ih_i^3/12$$ % % Substituting this into the equation for $\sigma_i$ gives % % $$\sigma_i = 6PD_i/b_ih_i^2$$ % % The bending stress in each part of the cantilever should not exceed the % maximum allowable stress, $\sigma_{max}$. Consequently, we can finally % state the five bending stress constraints (one for each step of the % cantilever) % % $$\frac{6Pl}{b_5h_5^2} \leq \sigma_{max}$$ % % $$\frac{6P(2l)}{b_4h_4^2} \leq \sigma_{max}$$ % % $$\frac{6P(3l)}{b_3h_3^2} \leq \sigma_{max}$$ % % $$\frac{6P(4l)}{b_2h_2^2} \leq \sigma_{max}$$ % % $$\frac{6P(5l)}{b_1h_1^2} \leq \sigma_{max}$$ % %% % % _Constraints on the Design : 2 - End deflection_ % % The end deflection of the cantilever can be calculated using % Castigliano's second theorem, which states that % % $$\delta = \frac{\partial U}{\partial P}$$ % % where $\delta$ is the deflection of the beam, $U$ is the energy stored % in the beam due to the applied force, $P$. % % The energy stored in a cantilever beam is given by % % $$U = \int_0^L \! M^2/2EI \, \mathrm{d} x$$ % % where $M$ is the moment of the applied force at $x$. % % Given that $M = Px$ for a cantilever beam, we can write the above % equation as % % $$U = % P^2/2E \int_0^l \! [(x+4l)^2/I_1 \, % + (x+3l)^2/I_2 \, % + (x+2l)^2/I_3 \, % + (x+l)^2/I_4 \, % + x^2/I_5 ]\, \mathrm{d} x$$ % % where $I_n$ is the area moment of inertia of the $n$-th part of the % cantilever. Evaluating the integral gives the following expression for % $U$. % % $$U = (P^2/2)(l^3/3E)(61/I_1 + 37/I_2 + 19/I_3 + 7/I_4 + 1/I_5)$$ % % Applying Castigliano's theorem, the end deflection of the beam is given % by % % $$\delta = Pl^3/3E(61/I_1 + 37/I_2 + 19/I_3 + 7/I_4 + 1/I_5)$$ % % Now, the end deflection of the cantilever, $\delta$, should be less than % the maximum allowable deflection, $\delta_{max}$, which gives us the % following constraint. % % $$Pl^3/3E(61/I_1 + 37/I_2 + 19/I_3 + 7/I_4 + 1/I_5) \leq % \delta_{max}$$ %% % % _Constraints on the Design : 3 - Aspect ratio_ % % For each step of the cantilever, the aspect ratio must not exceed a % maximum allowable aspect ratio, $a_{max}$. That is, % % $h_i/b_i \leq a_{max}$ for $i = 1, ..., 5$ %% State the Optimization Problem % % We are now able to state the problem to find the optimal parameters for % the stepped cantilever beam given the stated constraints. % % Let $x_1 = b_1$, $x_2 = h_1$, $x_3 = b_2$, $x_4 = % h_2$, $x_5 = b_3$, $x_6 = h_3$, $x_7 = b_4$, $x_8 = % h_4$, $x_9 = b_5$ and $x_{10} = h_5$ % % Minimize: % % $$V = l(x_1x_2 + x_3x_4 + x_5x_6 + x_7x_8 + x_9x_{10})$$ % % Subject to: % % $$\frac{6Pl}{x_9x_{10}^2} \leq \sigma_{max}$$ % % $$\frac{6P(2l)}{x_7x_8^2} \leq \sigma_{max}$$ % % $$\frac{6P(3l)}{x_5x_6^2} \leq \sigma_{max}$$ % % $$\frac{6P(4l)}{x_3x_4^2} \leq \sigma_{max}$$ % % $$\frac{6P(5l)}{x_1x_2^2} \leq \sigma_{max}$$ % % $$\frac{Pl^3}{E}(\frac{244}{x_1x_2^3} + \frac{148}{x_3x_4^3} + % \frac{76}{x_5x_6^3} + \frac{28}{x_7x_8^3} + % \frac{4}{x_9x_{10}^3}) \leq \delta_{max}$$ % % $$x_2/x_1 \leq 20$$, $$x_4/x_3 \leq 20$$, $$x_6/x_5 \leq 20$$, % $$x_8/x_7 \leq 20$$ and $$x_{10}/x_9 \leq 20$$ % % The first step of the beam can only be machined to the nearest % centimetre. That is, $x_1$ and $x_2$ must be integer. The remaining % variables are continuous. The bounds on the variables are given below:- % % $$1 \leq x_1 \leq 5$$ % % $$30 \leq x_2 \leq 65$$ % % $$2.4 \leq x_3, x_5 \leq 3.1$$ % % $$45 \leq x_4, x_6 \leq 60$$ % % $$1 \leq x_7, x_9 \leq 5$$ % % $30 \leq x_8, x_{10} \leq 65$ %% % % _Design Parameters for this Problem_ % % For the problem we will solve in this example, the end load that the beam % must support is $P = 50000 N$. % % The beam lengths and maximum end deflection are: % % * Total beam length, $L = 500 cm$ % * Individual section of beam, $l = 100 cm$ % * Maximum beam end deflection, $\delta_{max} = 2.7 cm$ % % The maximum allowed stress in each step of the beam, $\sigma_{max} = % 14000 N/cm^2$ % % Young's modulus of each step of the beam, $E = 2\times10^{7} N/cm^2$ %% Solve the Mixed Integer Optimization Problem % % We now solve the problem described in _State the Optimization % Problem_. %% % % _Define the Fitness and Constraint Functions_ % % Examine the MATLAB files |cantileverVolume.m| and % |cantileverConstraints.m| to see how the fitness and constraint % functions are implemented. % % _A note on the linear constraints_: When linear constraints are specified % to |ga|, you normally specify them via the |A|, |b|, |Aeq| and |beq| % inputs. In this case we have specified them via the nonlinear constraint % function. This is because later in this example, some of the variables % will become discrete. When there are discrete variables in the problem it % is far easier to specify linear constraints in the nonlinear constraint % function. The alternative is to modify the linear constraint matrices to % work in the transformed variable space, which is not trivial and maybe % not possible. Also, in the mixed integer |ga| solver, the linear % constraints are not treated any differently to the nonlinear constraints % regardless of how they are specified. %% % % _Set the Bounds_ % % Create vectors containing the lower bound (|lb|) and upper bound % constraints (|ub|). lb = [1 30 2.4 45 2.4 45 1 30 1 30]; ub = [5 65 3.1 60 3.1 60 5 65 5 65]; %% % % _Set the Options_ % % To obtain a more accurate solution, we increase the |PopulationSize|, and % |MaxGenerations| options from their default values, and decrease the % |EliteCount| and |FunctionTolerance| options. These settings cause |ga| % to use a larger population (increased PopulationSize), to increase the % search of the design space (reduced EliteCount), and to keep going until % its best member changes by very little (small FunctionTolerance). We % also specify a plot function to monitor the penalty function value as % |ga| progresses. % % Note that there are a restricted set of |ga| options available when % solving mixed integer problems - see Global Optimization Toolbox User's % Guide for more details. opts = optimoptions(@ga, ... 'PopulationSize', 150, ... 'MaxGenerations', 200, ... 'EliteCount', 10, ... 'FunctionTolerance', 1e-8, ... 'PlotFcn', @gaplotbestf); %% % % _Call |ga| to Solve the Problem_ % % We can now call |ga| to solve the problem. In the problem statement $x_1$ % and $x_2$ are integer variables. We specify this by passing the index % vector |[1 2]| to |ga| after the nonlinear constraint input and before % the options input. We also seed and set the random number generator here % for reproducibility. rng(0, 'twister'); [xbest, fbest, exitflag] = ga(@cantileverVolume, 10, [], [], [], [], ... lb, ub, @cantileverConstraints, [1 2], opts); %% % % _Analyze the Results_ % % If a problem has integer constraints, |ga| reformulates it internally. In % particular, the fitness function in the problem is replaced by a penalty % function which handles the constraints. For feasible population members, % the penalty function is the same as the fitness function. % % The solution returned from |ga| is displayed below. Note that the section % nearest the support is constrained to have a width ($x_1$) and height % ($x_2$) which is an integer value and this constraint has been honored by % GA. display(xbest); %% % We can also ask |ga| to return the optimal volume of the beam. fprintf('\nCost function returned by ga = %g\n', fbest); %% Add Discrete Non-Integer Variable Constraints % % The engineers are now informed that the second and third steps of the % cantilever can only have widths and heights that are chosen from a % standard set. In this section, we show how to add this constraint to the % optimization problem. Note that with the addition of this constraint, % this problem is identical to that solved in [1]. % % First, we state the extra constraints that will be added to the above % optimization % % * The width of the second and third steps of the beam must be chosen from % the following set:- [2.4, 2.6, 2.8, 3.1] cm % * The height of the second and third steps of the beam must be chosen % from the following set:- [45, 50, 55, 60] cm % % To solve this problem, we need to be able to specify the variables $x_3$, % $x_4$, $x_5$ and $x_6$ as discrete variables. To specify a component % $x_j$ as taking discrete values from the set $S = {v_1,\ldots,v_k}$, % optimize with $x_j$ an integer variable taking values from 1 to $k$, and % use $S(x_j)$ as the discrete value. To specify the range (1 to $k$), set % 1 as the lower bound and $k$ as the upper bound. % % So, first we transform the bounds on the discrete variables. Each set has % 4 members and we will map the discrete variables to an integer in the % range [1, 4]. So, to map these variables to be integer, we set the lower % bound to 1 and the upper bound to 4 for each of the variables. lb = [1 30 1 1 1 1 1 30 1 30]; ub = [5 65 4 4 4 4 5 65 5 65]; %% % % Transformed (integer) versions of $x_3$, $x_4$, $x_5$ and $x_6$ will now % be passed to the fitness and constraint functions when the |ga| solver is % called. To evaluate these functions correctly, $x_3$, $x_4$, $x_5$ and % $x_6$ need to be transformed to a member of the given discrete set in % these functions. To see how this is done, examine the MATLAB files % |cantileverVolumeWithDisc.m|, |cantileverConstraintsWithDisc.m| and % |cantileverMapVariables.m|. %% % % Now we can call |ga| to solve the problem with discrete variables. In % this case $x_1, ..., x_6$ are integers. This means that we pass the index % vector |1:6| to |ga| to define the integer variables. rng(0, 'twister'); [xbestDisc, fbestDisc, exitflagDisc] = ga(@cantileverVolumeWithDisc, ... 10, [], [], [], [], lb, ub, @cantileverConstraintsWithDisc, 1:6, opts); %% % _Analyze the Results_ % % |xbestDisc(3:6)| are returned from |ga| as integers (i.e. in their % transformed state). We need to reverse the transform to retrieve the % value in their engineering units. xbestDisc = cantileverMapVariables(xbestDisc); display(xbestDisc); %% % % As before, the solution returned from |ga| honors the constraint that % $x_1$ and $x_2$ are integers. We can also see that $x_3$, $x_5$ are % chosen from the set [2.4, 2.6, 2.8, 3.1] cm and $x_4$, $x_6$ are chosen % from the set [45, 50, 55, 60] cm. %% % Recall that we have added additional constraints on the variables |x(3)|, % |x(4)|, |x(5)| and |x(6)|. As expected, when there are additional % discrete constraints on these variables, the optimal solution has a % higher minimum volume. Note further that the solution reported in [1] has % a minimum volume of $64558 cm^3$ and that we find a solution which is % approximately the same as that reported in [1]. fprintf('\nCost function returned by ga = %g\n', fbestDisc); %% Summary % % This example illustrates how to use the genetic algorithm solver, |ga|, % to solve a constrained nonlinear optimization problem which has integer % constraints. The example also shows how to handle problems that have % discrete variables in the problem formulation. %% References % % [1] Survey of discrete variable optimization for structural design, P.B. % Thanedar, G.N. Vanderplaats, J. Struct. Eng., 121 (3), 301-306 (1995)