Optimization of Stalagmite function using Genetic Algorithm in MATLAB

Genetic Algorithm(GA)

AIM:- To optimise the stalagmite function and find the global maxima using “Genetic Algorithm”.

And explain the concept of genetic algorithm and also explain the syntax for ga in MATLAB.

Governing Equations:-

Stalagmite Function: `f(x,y)= f_(1,x)*f_(2,x)*f_(1,y)*f_(2,y)`

Where,

`f_(1,x)=[sin⁡(5.1*π*x+0.5)]^6`

`f_(1,y)=[sin⁡(5.1*π*y+0.5)]^6`

`f_(2,x)= e^([-4*ln⁡(2)(x-0.0667)^2/0.64])`

`f_(2,y)= e^([-4*ln⁡(2)(y-0.0667)^2/0.64])`

Objectives of the Project :-

To optimize stalagmite function to find its global maximum using GA.

Theory about ‘Genetic Algorithm’ :-

Genetic Algorithm (GA) is a search-based optimization technique based on the principles of Genetics and Natural Selection. It is frequently used to find optimal or near-optimal solutions to difficult problems which otherwise would take a lifetime to solve.

In GAs, we have a pool or a population of possible solutions to the given problem. These solutions then undergo recombination and mutation (like in natural genetics), producing new children, and the process is repeated over various generations. Each individual (or candidate solution) is assigned a fitness value (based on its objective function value) and the fitter individuals are given a higher chance to mate and yield more “fitter” individuals. This is in line with the Darwinian Theory of “Survival of the Fittest”.

In this way we keep “evolving” better individuals or solutions over generations, till we reach a stopping criterion.

Genetic Algorithms are sufficiently randomized in nature, but they perform much better than random local search (in which we just try various random solutions, keeping track of the best so far), as they exploit historical information as well.

Basic Terminology in GA,

• Population− It is a subset of all the possible (encoded) solutions to the given problem. The population for a GA is analogous to the population for human beings except that instead of human beings, we have Candidate Solutions representing human beings.
• Chromosomes− A chromosome is one such solution to the given problem.
• Gene− A gene is one element position of a chromosome.
• Allele− It is the value a gene takes for a particular chromosome.
• Genotype− Genotype is the population in the computation space. In the computation space, the solutions are represented in a way which can be easily understood and manipulated using a computing system.
• Phenotype− Phenotype is the population in the actual real world solution space in which solutions are represented in a way they are represented in real world situations.
• Decoding and Encoding− For simple problems, the phenotype and genotype spaces are the same. However, in most of the cases, the phenotype and genotype spaces are different. Decoding is a process of transforming a solution from the genotype to the phenotype space, while encoding is a process of transforming from the phenotype to genotype space. Decoding should be fast as it is carried out repeatedly in a GA during the fitness value calculation.

For example, consider the 0/1 Knapsack Problem. The Phenotype space consists of solutions which just contain the item numbers of the items to be picked.

However, in the genotype space it can be represented as a binary string of length n (where n is the number of items). A 0 at position x represents that x^th item is picked while a 1 represents the reverse. This is a case where genotype and phenotype spaces are different.

• Fitness Function− A fitness function simply defined is a function which takes the solution as input and produces the suitability of the solution as the output. In some cases, the fitness function and the objective function may be the same, while in others it might be different based on the problem.
• Genetic Operators− These alter the genetic composition of the offspring. These include crossover, mutation, selection, etc.

The basic structure of a GA is as follows −

We start with an initial population (which may be generated at random or seeded by other heuristics), select parents from this population for mating. Apply crossover and mutation operators on the parents to generate new off-springs. And finally these off-springs replace the existing individuals in the population and the process repeats. In this way genetic algorithms actually try to mimic the human evolution to some extent.

Syntax for GA in MATLAB® :-

X = ga(FITNESSFCN,NVARS) finds a local unconstrained minimum X to the FITNESSFCN using ga. NVARS is the dimension (number of design variables) of the FITNESSFCN.

FITNESSFCN accepts a vector X of size 1-by-NVARS, and returns a scalar evaluated at X.

X = ga(FITNESSFCN,NVARS,A,b) finds a local minimum X to the function FITNESSFCN, subject to the linear inequalities A*X <= B. Linear constraints are not satisfied when the PopulationType option is set to 'bitString' or 'custom'. See the documentation for details.

X = ga(FITNESSFCN,NVARS,A,b,Aeq,beq) finds a local minimum X to the function FITNESSFCN, subject to the linear equalities Aeq*X = beq as well as A*X <= B. (Set A=[] and B=[] if no inequalities exist.) Linear constraints are not satisfied when the PopulationType option is set to 'bitString' or 'custom'. See the documentation for details.

X = ga(FITNESSFCN,NVARS,A,b,Aeq,beq,lb,ub) defines a set of lower and upper bounds on the design variables, X, so that a solution is found in the range lb <= X <= ub. Use empty matrices for lb and ub if no bounds exist. Set lb(i) = -Inf if X(i) is unbounded below; set ub(i) = Inf if X(i) is unbounded above. Linear constraints are not satisfied when the PopulationType option is set to 'bitString' or 'custom'. See the documentation for details.

X = ga(FITNESSFCN,NVARS,A,b,Aeq,beq,lb,ub,NONLCON) subjects the minimization to the constraints defined in NONLCON. The function NONLCON accepts X and returns the vectors C and Ceq, representing the nonlinear inequalities and equalities respectively. ga minimizes FITNESSFCN such that C(X)<=0 and Ceq(X)=0. (Set lb=[] and/or ub=[] if no bounds exist.) Nonlinear constraints are not satisfied when the PopulationType option is set to 'bitString' or 'custom'.

X = ga(FITNESSFCN,NVARS,A,b,Aeq,beq,lb,ub,NONLCON,options) minimizes with the default optimization parameters replaced by values in OPTIONS. OPTIONS can be created with the OPTIMOPTIONS function.

X = ga(FITNESSFCN,NVARS,A,b,[],[],lb,ub,NONLCON,INTCON) requires that the variables listed in INTCON take integer values. Note that ga does not solve problems with integer and equality constraints. Pass empty matrices for the Aeq and beq inputs if INTCON is not empty.

X = ga(FITNESSFCN,NVARS,A,b,[],[],lb,ub,NONLCON,INTCON,options)

minimizes with integer constraints and the default optimization parameters replaced by values in OPTIONS. OPTIONS can be created with the OPTIMOPTIONS function. See OPTIMOPTIONS for details.

X = ga(PROBLEM) finds the minimum for PROBLEM. PROBLEM is a structure that has the following fields:

fitnessfcn:

nvars:

Aineq:

bineq:

Aeq:

beq:

lb:

ub:

nonlcon:

intcon:

options:

rngstate:

[X,FVAL] = ga(FITNESSFCN, ...) returns FVAL, the value of the fitness function FITNESSFCN at the solution X.

[X,FVAL,EXITFLAG] = ga(FITNESSFCN, ...) returns EXITFLAG which describes the exit condition of ga. Possible values of EXITFLAG and the corresponding exit conditions are,

1 Average change in value of the fitness function over

options.MaxStallGenerations generations less than

options.FunctionTolerance and constraint violation less than

options.ConstraintTolerance.

3 The value of the fitness function did not change in

options.MaxStallGenerations generations and constraint violation

less than options.ConstraintTolerance.

4 Magnitude of step smaller than machine precision and constraint

violation less than options.ConstraintTolerance. This exit

condition applies only to nonlinear constraints.

5 Fitness limit reached and constraint violation less than

options.ConstraintTolerance.

0 Maximum number of generations exceeded.

-1 Optimization terminated by the output or plot function.

-2 No feasible point found.

-4 Stall time limit exceeded.

-5 Time limit exceeded.

[X,FVAL,EXITFLAG,OUTPUT] = ga(FITNESSFCN, ...) returns a structure OUTPUT with the following information:

rngstate:

generations:

funccount:

maxconstraint: , if any

message:

[X,FVAL,EXITFLAG,OUTPUT,POPULATION] = ga(FITNESSFCN, ...) returns the final POPULATION at termination.

[X,FVAL,EXITFLAG,OUTPUT,POPULATION,SCORES] = ga(FITNESSFCN, ...) returns the SCORES of the final POPULATION.

Function Code :-

function [f] = stalagmite(input_vector)

x = input_vector(1);
y = input_vector(2);

f_1_x = (sin(5.1*pi*x + 0.5))^6;
f_1_y = (sin(5.1*pi*y + 0.5))^6;

f_2_x = exp(-4*log(2)*(x-0.0667)^2/0.64);
f_2_y = exp(-4*log(2)*(y-0.0667)^2/0.64);

f = -(f_1_x)*(f_1_y)*(f_2_x)*(f_2_y);
end

Main Code :-

close all 
clear all
clc

%Stalagmite Model:

%Defining search space.
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

% Creating a 2D-Mesh 
[xx yy] = meshgrid(x,y);
input_vector = zeros(1,1);

f = zeros(150,150);
% Evaluating the stalagmite function 
for i = 1:length(xx)
    for j= 1 :length(yy)
        input_vector(1) = xx(i,j);
        input_vector(2) = yy(i,j);
        f(i,j)= stalagmite(input_vector);
    end
end

CASE: 1 Unbounded Inputs
% Number of iterations
num_cases_1 = 10;
tic
fopt_1 = zeros(num_cases_1);
for i = 1:num_cases_1
    [inputs, fopt_1(i),exitFlag, out] = ga(@stalagmite,2);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
study1_time = toc

figure(1)
subplot(2,1,1)
hold on 
surfc(x,y,-f)
xlabel('X-Axis')
ylabel('Y-Axis')
zlabel('Function Output')
shading interp
colorbar
plot3(xopt,yopt,-fopt_1,'marker','o','markersize',5,'markerfacecolor','r')
title('Unbounded Inputs')
subplot(2,1,2)
plot(-fopt_1)
axis([0 num_cases_1 0 1.5])
xlabel('Iterations')
ylabel('Function Maximium')


% CASE: 2 Bounded Inputs

% Number of iterations
num_cases_2 = 10;
tic
fopt_2 = zeros(num_cases_2);
for i = 1:num_cases_2
    [inputs, fopt_2(i),exitFlag, out] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1]);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end

study2_time = toc

figure(2)
subplot(2,1,1)
hold on 
surfc(x,y,-f)
xlabel('X-Axis')
ylabel('Y-Axis')
zlabel('Function Output')
shading interp
colorbar
plot3(xopt,yopt,-fopt_2,'marker','o','markersize',5,'markerfacecolor','r')
title('Bounded Inputs')
subplot(2,1,2)
plot(-fopt_2)
axis([0 num_cases_2 0 1.5])
xlabel('Iterations')
ylabel('Function Maximium')

% CASE: 3 Bounded Outputs & Increase in Population_Size

% Number of iterations
num_cases_3 = 10;
option_1 = optimoptions('ga');
option_1 = optimoptions(option_1,'Populationsize',300);

tic
fopt_3 = zeros(num_cases_3);
for i = 1:num_cases_3
    [inputs, fopt_3(i),exitFlag, out] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1],[],[],option_1);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end


study3_time = toc

figure(3)
subplot(2,1,1)
hold on 
surfc(x,y,-f)
xlabel('X-Axis')
ylabel('Y-Axis')
zlabel('Function Output')
shading interp
colorbar
plot3(xopt,yopt,-fopt_3,'marker','o','markersize',5,'markerfacecolor','r')
title('Bounded Inputs With Increased Population Size')
subplot(2,1,2)
plot(-fopt_3)
axis([0 num_cases_3 0 1.5])
xlabel('Iterations')
ylabel('Function Maximium')
inputs

Explanation of Code :-

1. A program for the Stalagmite function is written first.

Solving the Genetic algorithm for stalagmite function gives us a minimum value of the function. To find the maximum value of the function, A subtraction sign is added to the final answer of stalagmite function so that it will give value for the reverse case of minima, i.e. maxima.

1. Now starting with the main code,The values of x & y in the stalagmite function are defined using ‘linspace()’ command to be between 0 to 0.6 with 150 intervals. Next meshgrid(x,y) is used as the inputs to stalagmite function.
2. Stalagmite function is evaluated using nested for loop, & the function outputs are stored in variable ‘ f ’.
3. The next step is to obtain the results for 3 cases of the genetic algorithm in MATLAB. In Study_1, maximum value of the function is found out with no restrictions, no bounds, and no limits.
4. To get the same value each time the program is executed, Study_2 The program is restricted within specific bounds namely upper bound and lower bound. The whole procedure for programming this case the same as that of the previous one, exception being for GA
5. Even the results from the bounded cases are not the same every time the program is executed. Study_3 Hence we now need to increase GA iterations. We can raise it to any value, we have chosen it as 300. In this case, syntax for GA is different than the previous case. GA has to be given population size via an option. Results for this case will be just one value of function which does not deviate with no of iterations.

Outputs Screens :-

Study_1: No restrictions, no bounds, and no limits.

Time for Study_1 = 0.6142 s

Study_2: Restricted Boundary

Time for Study_2 = 0.7679 s

Study_3 : Increse in Population size with restricted boundary.

Time for Study_3 = 2.5597 s

Conclusion:-

By observing above plots, It can be concluded that with additions of restrictions in upper/lower bound limits & increasing the population size a constant value of (x,y) = (0.0668,0.0668) is obtained for global maxima(fopt_3 =1)

It can also be observed that the time taken for Study_3 is more than that of 1 or 2, this limits the number of iterations which can be performed for Study_3 as it requirees more computational power.

References :-

https://www.tutorialspoint.com/genetic_algorithms/index.htm

https://in.mathworks.com/help/gads/examples/genetic-algorithm-options.html

https://in.mathworks.com/help/gads/ga.html