Optimization

321

Appendix B: MATLAB

®

Code

Chapter 1

Code Name Details

graph_examp12.m Solves Example 1.2 using the graphical method

graph_examp14.m Solves Example 1.4 using the graphical method

convexity.m Plots some convex functions

derivative.m Computes and plots rst and second derivatives of a function

grad.m Plots the gradient vector

positive_denite.m Checks whether the square matrix is positive denite

quadr.m Linear and quadratic approximations of a function

quadr_examp12.m Linear and quadratic approximations of a function given in Example 1.2

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% MATLAB code graph_examp12.m

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

% x1 -> radius of can

% x2 -> height of can

% area -> area of can

% pi -> MATLAB variable

% x,y,z -> array of design points

% vv -> user identified contour values

% cons_x2-> value of x2 when constraint is active

% contour -> MATLAB function to generate contours

% xlabel, ylabel, legend, hold on -> MATLAB functions

%

for x1 = 1:100

for x2 = 1:200

area = 2*pi*x1*x2 + 2*pi*x1*x1;

x(x1,x2) = x1;

y(x1,x2) = x2;

z(x1,x2) = area;

end

end

vv = [15000;26436;50000;70000;200000];

[c, h] = contour(x,y,z,vv); clabel(c, h);

hold on

322 Appendix B

for x1 = 10:100

cons_x2 = 330000/(pi*x1*x1);

plot(x1,cons_x2,'*')

hold on

end

xlabel('x_1, mm')

ylabel('x_2, mm')

legend('Objective Function','Constraint')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% MATLAB code graph_examp14.m

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

% L -> length of rod

% rho -> density of rod material

% d -> diameter of rod

% m -> mass of rod

% I -> moment of inertia

% k -> mass per unit length

% f1 -> frequency

% plot, xlabel, ylabel, legend -> MATLAB function

%

L = 1;

rho = 7800;

E = 2e11;

for d = 0.0:0.001:0.05

mass = (pi/4)*d*d*L*rho;

k = mass/L;

I = (pi/64)*d^4;

f1 = (1/(2*pi))*(3.5156/(L*L))*sqrt((E*I)/k);

plot(d,mass,'+')

hold on

plot(d,f1,'*')

hold on

end

xlabel('d, m')

ylabel('objective function (kg), constraint(Hz)')

legend('Objective Function','Constraint')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% MATLAB code convexity.m

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

% x -> independent variable

% y -> dependent variable

323Appendix B

% plot, xlabel, ylabel, meshgrid -> MATLAB functions

%

x = meshgrid(-2:0.01:2);

y = x.^2;

subplot(2,2,1), plot(x,y)

xlabel('x')

ylabel('y')

hold on

y = exp(x);

subplot(2,2,2), plot(x,y)

xlabel('x')

ylabel('y')

hold on

y = exp(y);

subplot(2,2,3), plot(x,y)

xlabel('x')

ylabel('y')

hold on

y = exp(x.^2);

subplot(2,2,4), plot(x,y)

xlabel('x')

ylabel('y')

hold on

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% MATLAB code derivative.m

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

% delx -> delta-x

% fx -> f(x)

% deriv -> derivative of the function at xd

% sderiv -> second derivative of the function at xdd

% signchange -> change of derivative sign

% locatepoints -> point at which derivative changes sign

% plot, subplot, xlabel, ylabel, hold -> MATLAB functions

%

delx = 0.01;

x=0.1:delx:1.0;

fx = @(x)2*sin(5*x)+3*x.^3-2*x.^2+3*x-5;

subplot(311), plot(x,fx(x),'LineWidth',2)

hold on

ylabel('f(x)')

grid on

for i = 2:length(x)-1

xd(i-1) = x(i);

deriv(i-1) = (fx(x(i+1))-fx(x(i-1)))/(2*delx);

end

subplot(312), plot(xd,deriv,'LineWidth',2)

grid on

324 Appendix B

hold on

ylabel('f''(x)')

signchange = deriv(1:length(deriv)-1).* deriv(2:length(deriv));

locatepoints = xd(find(signchange<0))

subplot(311), plot(locatepoints,fx(locatepoints),'r*')

subplot(312), plot(xd(find(signchange<0)),deriv(find(signchange<0)),'r*')

for ii = 2:length(xd)-1

xdd(ii-1) = xd(ii);

sderiv(ii-1) = (fx(xd(ii+1))+fx(xd(ii-1))-2*fx(xd(ii)))/

(delx*delx);

end

subplot(313), plot(xdd,sderiv,'LineWidth',2)

grid on

hold on

subplot(313), plot(xdd(find(signchange<0)),sderiv(find(signchange<0)),'r*')

xlabel('x')

ylabel('f''''(x)')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% MATLAB code grad.m

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%

% x1 -> radius of can

% x2 -> height of can

% area -> area of can

% pi -> MATLAB variable

% x,y,z -> array of design points

% vv -> user identified contour values

% xit, x2t -> identified point at which gradient required

% contour -> MATLAB function to generate contours

% xlabel, ylabel, legend, plot, hold on -> MATLAB functions

%

clear all

clc

for x1 = 1:100

for x2 = 1:100

area = 2*pi*x1*x2 + 2*pi*x1*x1;

x(x1,x2) = x1;

y(x1,x2) = x2;

z(x1,x2) = area;

end

end

vv = [5000,15000,30000,50000,70000,90000];

[c, h] = contour(x,y,z,vv); clabel(c, h);

hold on

x1t = 25;

x2t = 70.493;

slope = (x2t+2*x1t)/x1t;

i = 1;

325Appendix B

for delx1 = -10:10

delx2 = -slope*delx1;

x11(i) = x1t+delx1;

x22(i) = x2t+delx2;

i = i+1;

end

plot(x11,x22,'r--')

hold on

i = 1;

for delx1 = -10:10

delx2 = (1/slope)*delx1;

x11(i) = x1t+delx1;

x22(i) = x2t+delx2;

i = i+1;

end

plot(x11,x22,'b+')

xlabel('x_1, mm')

ylabel('x_2, mm')

legend('Objective Function','Tangent','Gradient')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% MATLAB code positive_definite.m

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% H -> hessian matrix

% eig, fprintf -> MATLAB function

% eigenvalues -> of the hessian matrix

%

H = [2 1 1;

1 2 1;

0 2 3];

eigenvalues = eig(H);

eigenvalues

if eigenvalues >= 0

fprintf('The matrix is positive definite\n')

else

fprintf('The matrix is not positive definite\n')

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%