OPTIMISE THE STALAGMITE FUNCTION AND FIND THE GLOBAL MAXIMA OF THE FUNCTION USING GENETIC ALGORITHM IN MATLAB

AIM

To write a code in MATLAB to optimise the stalagmite function and find the global maxima of the function using the Genetic Algorithm.

GOVERNING EQUATION:

Formula used for run the stalagmite function,

`f(x,y)=f1,x f2,x f1,y f2,y`

`f1.x=[sin(5.1πx+0.5) ]^6`

`f1.y=[sin(5.1πy+0.5) ]^6`

`f2.x=exp[-4ln⁡(2)〖(x-0.0667)〗^2/0.64]`

`f2.y=exp[-4ln⁡(2)〖(y-0.0667)〗^2/0.64]`

OBJECTIVE

1. To define the stalagmite function
2. To find the optimum value of the stalagmite function using various methods
3. Explain the concept of the Genetic Algorithm
4. To find the Global maxima value

METHODS TO FIND THE OPTIMUM VALUE OF THE STALAGMITE FUNCTION

1. Statistical behavior method
2. Statistical behavior with bounds method
3. Statistical behavior with bounds and population size method

THEORY

ABOUT OPTIMIZATION

Optimization is the process of making something better. In any process, we have a set of inputs and a set of outputs as shown in the following figure.

Optimization refers to finding the values of inputs in such a way that we get the “best” output values. The definition of “best” varies from problem to problem, but in mathematical terms, it refers to maximizing or minimizing one or more objective functions, by varying the input parameters.

The set of all possible solutions or values which the inputs can take make up the search space. In this search space, lies a point or a set of points which gives the optimal solution. The aim of optimization is to find that point or set of points in the search space.

Optimisation is a decision making process, it helps to find the best results

under a set of prioritized criteria or constraints. It is a procedure which is executed iteratively by comparing various solutions till an optimum or a satisfactory solution is found. These include maximizing factors such as productivity, strength, reliability, efficiency, utilization, etc

Example:

Let’s take, you are designing a simply supported beam to withstand the load of 1-tonne, in this case you need to consider the following parameters: Length, breadth, width, factor of safety, weight and stress induced. Here optimization helps to find the best value for length, breadth and width of the beam in such a way that the factor of safety, material consumed and stress induced are within the limit.

METHODS OF OPTIMIZATION

Deterministic Method

They use specific rules for moving one solution to another. In this technique, the output of the model is fully determined by the parameter values and the initial conditions.

Example:

Stress analysis of the component, this is a deterministic approach since the output results(stress values) fully depend on the initial conditions(loading condition) and the parameter values(force acting on the component and material properties).

Stochastic Method

These are gaining popularity due to certain properties which deterministic algorithms do not possess. Stochastic models possess some inherent randomness, it is a probabilistic approach.

Example:

Varying the water paths from the mountain can vary the regional weather conditions, but the weather condition is not fully dependent on the water paths, hence there will be many statistical predictions to get the predictive model, since the weather condition depends on various factors such as solar distance, latitudinal location, air pressure, etc. Hence it is a predictive model.

The following is the procedure followed by any optimisation process.

ABOUT GENETIC ALGORITHM

Nature has always been a great source of inspiration to all mankind. Genetic Algorithms (GAs) are search based algorithms based on the concepts of natural selection and genetics. GAs are a subset of a much larger branch of computation known as Evolutionary Computation.

GAs were developed by John Holland and his students and colleagues at the University of Michigan, most notably David E. Goldberg and has since been tried on various optimization problems with a high degree of success.

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.

Advantages of GAs

GAs have various advantages which have made them immensely popular. These include −

• Does not require any derivative information (which may not be available for many real-world problems).
• Is faster and more efficient as compared to the traditional methods.
• Has very good parallel capabilities.
• Optimizes both continuous and discrete functions and also multi-objective problems.
• Provides a list of “good” solutions and not just a single solution.
• Always gets an answer to the problem, which gets better over the time.
• Useful when the search space is very large and there are a large number of parameters involved.

CONCEPT OF THE 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. It is frequently used to solve optimization problems, in research, and in machine learning. We need to choose an optimization algorithm to perform the optimisation process. Genetic algorithm is one of the optimisation processes.

The following flow chart describes how the genetic algorithm works to find the optimal solution.

Five phases are considered in a genetic algorithm.

1. Initial population
2. Fitness function
3. Selection
4. Crossover
5. Mutation

Initial Population

The process begins with a set of individuals which is called a Population. Each individual is a solution to the problem you want to solve.

Fitness Function

The fitness function determines how fit an individual is (the ability of an individual to compete with other individuals). It gives a fitness score to each individual. The probability that an individual will be selected for reproduction is based on its fitness score.

Selection

The idea of selection phase is to select the fittest individuals and let them pass their genes to the next generation.

Crossover

Crossover is the most significant phase in a genetic algorithm. For each pair of parents to be mated, a crossover point is chosen at random from within the genes.

Mutation

In certain new offspring formed, some of their genes can be subjected to a mutation with a low random probability. This implies that some of the bits in the bit string can be flipped.

Termination

The algorithm terminates if the population has converged (does not produce offspring which are significantly different from the previous generation). Then it is said that the genetic algorithm has provided a set of solutions to our problem.

SYNTAX FOR GA

1) x = ga(fun,nvars)

Finds a local unconstrained minimum, x, to the objective function, fun. nvars is the dimension (number of design variables) of fun.

2) x = ga(fun,nvars,A,b):

Finds a local minimum x to fun, subject to the linear inequalities A*x ≤ b. ga evaluates the matrix product A*x as if x is transposed (A*x').

3) x = ga(fun,nvars,A,b,Aeq,beq):

Finds a local minimum x to fun, subject to the linear equalities Aeq*x = beq and A*x ≤ b. (Set A=[] and b=[] if no linear inequalities exist.) ga evaluates the matrix product Aeq*x as if x is transposed (Aeq*x').

4) x = ga(fun,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. (Set Aeq=[] and beq=[] if no linear equalities exist.)

5) x = ga(fun,nvars,A,b,Aeq,beq,lb,ub,nonlcon):

Subjects the minimization to the constraints defined in nonlcon. The function nonlcon accepts x and returns vectors C and Ceq, representing the nonlinear inequalities and equalities respectively. ga minimizes the fun such that C(x) ≤ 0 and Ceq(x) = 0. (Set lb=[] and ub=[] if no bounds exist.)

6) x = ga(fun,nvars,A,b,Aeq,beq,lb,ub,nonlcon,options):

Minimizes with the default optimization parameters replaced by values in options. (Set nonlcon=[] if no nonlinear constraints exist.) Create options using optimoptions.

7) x = ga(fun,nvars,A,b,[],[],lb,ub,nonlcon,IntCon):

Requires that the variables listed in IntCon take integer values.

8) x = ga(problem):

Finds the minimum for problem, where problem is a structure.

9) [x,fval] = ga(___):

For any previous input arguments, also returns fval, the value of the fitness function at x.

10) [x,fval,exitflag,output] = ga(___):

Also returns exitflag, an integer identifying the reason the algorithm terminated, and output, a structure that contains output from each generation and other information about the performance of the algorithm.

11) [x,fval,exitflag,output,population,scores] = ga(___):

Also returns a matrix population, whose rows are the final population, and a vector scores, the scores of the final population.

CODE FOR STALAGMITE FUNCTION(SUB PROGRAM)

%Program to define the stalagmite function
function[f] = stalagmite(input_vector)
% Inputs for stalagmite function
x = input_vector(1);
y = input_vector(2);
% Formula of finding optimum values
f1x =  (sin((5.1*pi*x)+0.5))^6;
f1y =  (sin((5.1*pi*y)+0.5))^6;
f2x = exp((-4*log(2))*(((x-0.0667)^2)/0.64));
f2y = exp((-4*log(2))*(((y-0.0667)^2)/0.64));
f = -(f1x*f2x*f1y*f2y);

end

MAIN CODE FOR FIND LOCAL MINIMA OF STALAGMITE FUNCTION

%Program for finding Local minima of a Stalagmite function Using Genetic Algorithm

% Basic start commands
clear all
close all
clc

% Inputs of one dimensional array
% Defining our search place
x = linspace(0,0.6,150);% x axis values of array
y = linspace(0,0.6,150);% y axis values of array

% Creating 2 Dimensional mesh
[xx yy] = meshgrid(x,y);% Creating 2d array for stalagmite function

% Evaluating the stalagmite function
for i = 1:length(xx)% value starts from frist value to end of xx value
    for j = 1:length(yy)% value starts from frist value to end of yy value
        input_vector(1)=xx(i,j);% random vvalue selected from xx for as input of stalagmite function
        input_vector(2)=yy(i,j);% random vvalue selected from yy for as input of stalagmite function
        f(i,j) = stalagmite(input_vector);% function of finding function maximum
    end
end
surf(xx,yy,f)% command for creating 3D plot structure
shading interp % command for shading and smoothening the plot to avoid grid lines
xlabel('xvalue')% name of x axis
ylabel('yvalue')% name of y axis

% Syntax to find the Local unconstrained minima
[inputs,Fopt] = ga(@stalagmite,2);% output of Fopt is optimum Function of GA minimum value, inputs are x and y values

WRITING PROGRAM WITH DIFFERENT OPTIMIZATION STUDIES

STUDY-1 STATISTICAL BEHAVIOR - UNBOUNDED INPUTS

• Initially, we use the one Dimensional arrays x and y with linspace command and as we know the x and y are input array it does not define the workspace for the genetic algorithm. Here we are using the component used in the program as follows.
• num_cases = 50 indicates that running the genetic algorithm 50 times.
• And next , we are just calling the (ga) function and passing a number of variables as an input which is equal to 2 and storing the x and y value also the (fopt) optimum value which function is produced.
• after that opening new figure window [figure(1)] and creating the subplot (which creating the multiple plots in the same figure window) for the two plot windows in the same figure.
• the first plot is surfc(x,y,F) and using the plot3 command for x, y and the optimum value of (ga) predict's over the 50-studies. which help us to know how the algorithm is performing.
• Now, the second subplot in which directly optimum value is get plotted as a function of iterations.
• To measure how long it takes to measure the study we are using tic and study1_time=toc.

CODE FOR STUDY 1 STATISTICAL BEHAVIOR OF UNBOUNDED INPUTS

%Program for finding Global maxima of a Stalagmite function Using Genetic Algorithm
clear all
close all
clc

% Defining our search place
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);
% Total num of Iteration process
n=50;

% Creating 2 Dimensional mesh
[xx,yy] = meshgrid(x,y);

% 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

%Study1 Statical behaviour with unbounded inputs
tic % strating time to record the process
for i = 1:n
    %Syntax to find all minima and maxima .
    [inputs,Fopt(i)] = ga(@stalagmite,2);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
Study1_time = toc% strating time to record the process
subplot(2,1,1)% create multi plots in same figure
hold on 
surf(xx,yy,-f)
shading interp
xlabel('x values')
ylabel('y values')
title('Unbounded Inputs') % title of the plot
plot3(xopt,yopt,-Fopt,'marker','o','markersize',5,'markerfacecolor','r')% 3D plot for function maximum and minimum values 
subplot(2,1,2)
plot(-Fopt)% plot for optimum function maximum and iteration 
xlabel('Iteration')
ylabel('Function maximum')

OUTPUT SCREENSHOT FOR STUDY 1

OUTPUT PLOT FOR STUDY 1

Now in the first study, we get multiple results, due to random sampling and maximum value changes from between 0.2 to 0.8, also time is taken to the first study is 7.6273 seconds.

study1_time = 5.4318

Hence, Further optimization is required.

STUDY 2 STATISTICAL BEHAVIOUR OF UPPER AND LOWER BOUND

In study 2 we are going to use the upper and lower limit which is known as upper and lower bounds. the syntax use is based on our referred syntax section.

Where:-

[0;0] = upper and lower limit of the x.

[1;1] = upper and lower limit of the y.

CODE FOR STUDY 2 STUDY 2 STATISTICAL BAHAVIOR OF UPPER AND LOWER BOUND

%Program for finding Global maxima of a Stalagmite function Using Genetic Algorithm

clear all
close all
clc

% Defining our search place
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

n=50;

% Creating 2 Dimensional mesh
[xx,yy] = meshgrid(x,y);

% 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

%Study2 Statical behaviour with upper and lower bounded inputs
tic
for i = 1:n
    % Syntax to find upper and lower bound maxima, minima values
    [inputs,Fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1]);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
Study2_time = toc
subplot(2,1,1)
hold on
surf(xx,yy,-f)
shading interp
xlabel('x values')
ylabel('y values')
title('Bounded Inputs')
plot3(xopt,yopt,-Fopt,'marker','o','markersize',5,'markerfacecolor','r')
subplot(2,1,2)
plot(-Fopt)
xlabel('Iteration')
ylabel('Function maximum')

OUTPUT SCREENSHOT FOR STUDY 2

OUTPUT PLOT FOR STUDY 2

Now, we get results as more clear from the previous one.

Here, Maximum value is between 0.7 to 0.85 which is much less range. also the time is taken to the second study is 6.1794 seconds.

study2_time = 6.1794

But still, optimization is needed.

STUDY 3 - INCREASING GA ITERATIONS

To optimize the results and getting more accuracy we are increasing the iteration in genetic algorithm. Because by default the (ga) takes 50 value. Here we are going to increase the population size because if initial population is too small then it not going to find out global maxima (or) minima, it only finds the local maxima (or) minima.

The (ga) value is 50 which is changeable by initially (i.e we can use the increased value by default) in the individual program, but in the combined program we need to use the command as follows,

options=optimoptions('ga')
options=optimoptions(options,'PopulationSize',250)

here we are increasing the size to 250 but this is going to slow down the (ga).

CODE FOR INCREASING GA ITERATIONS

%Program for finding Global maxima of a Stalagmite function Using Genetic Algorithm
clear all
close all
clc

% Defining our search place
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

n=50;

% Creating 2 Dimensional mesh
[xx,yy] = meshgrid(x,y);

% 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

%Study-3 Increasing GA Iteration process
options = optimoptions('ga') % option command to changes in ga
options = optimoptions(options,'PopulationSize',250)% increase the population size of iteration process

tic
for i = 1:n
    % syntax of finding maximum value and increase the iteration process
    [inputs,Fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1],[],[],options);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
Study3_time = toc
subplot(2,1,1)
hold on
surf(xx,yy,-f)
shading interp
xlabel('x values')
ylabel('y values')
title('Increased Iteration')
plot3(xopt,yopt,-Fopt,'marker','o','markersize',5,'markerfacecolor','r')
subplot(2,1,2)
plot(-Fopt)
xlabel('Iteration')
ylabel('Function maximum')

OUTPUT SCREENSHOT FOR STUDY 3

OUTPUT PLOT FOR STUDY 3

Now, we get the most accurate result which is between 0.9999997 to 0.9999998.

time is taken to the third study is 23.5356 seconds.

study3_time = 23.5356

CODE FOR RESULT AND COMBAINED STUDY PROGRAM

To get Final combined result we have combined all the Three studies in the same program as follows,

%Program for finding Global maxima of a Stalagmite function Using Genetic Algorithm
clear all
close all
clc

% Defining our search place
x = linspace(0,0.6,150);
y = linspace(0,0.6,150);

n=50;

% Creating 2 Dimensional mesh
[xx,yy] = meshgrid(x,y);

% 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

%Study1 Statical behaviour with bounded inputs
tic
for i = 1:n
    [inputs,Fopt(i)] = ga(@stalagmite,2);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
Study1_time = toc
figure(1)
subplot(2,1,1)
hold on
surf(xx,yy,-f)
shading interp
xlabel('x values')
ylabel('y values')
title('Unbounded Inputs')
plot3(xopt,yopt,-Fopt,'marker','o','markersize',5,'markerfacecolor','r')
subplot(2,1,2)
plot(-Fopt)
xlabel('Iteration')
ylabel('Function maximum')

%Study2 Statical behaviour with upper and lower bounded inputs
tic
for i = 1:n
    [inputs,Fopt(i)] = 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
surf(xx,yy,-f)
shading interp
xlabel('x values')
ylabel('y values')
title('Bounded Inputs')
plot3(xopt,yopt,-Fopt,'marker','o','markersize',5,'markerfacecolor','r')
subplot(2,1,2)
plot(-Fopt)
xlabel('Iteration')
ylabel('Function maximum')

%Study-3 Increasing GA Iteration process
options = optimoptions('ga')
options = optimoptions(options,'PopulationSize',250)

tic
for i = 1:n
    [inputs,Fopt(i)] = ga(@stalagmite,2,[],[],[],[],[0;0],[1;1],[],[],options);
    xopt(i) = inputs(1);
    yopt(i) = inputs(2);
end
Study3_time = toc
figure(3)
subplot(2,1,1)
hold on
surf(xx,yy,-f)
shading interp
xlabel('x values')
ylabel('y values')
title('Increased Iteration')
plot3(xopt,yopt,-Fopt,'marker','o','markersize',5,'markerfacecolor','r')
subplot(2,1,2)
plot(-Fopt)
xlabel('Iteration')
ylabel('Function maximum')

OUTPUT FOR COMBAINED PROGRAM

Now, we get the most accurate result which is between 0.9999997 to 0.9999998.

To get the maximum value of the function we are using the command in the program of study-3 as,

max_value = Fopt

Finally, we get the final result as,

max_value = 1.0000

EXPLAINATION

• Start the program with 'clear all','close all', and 'clc' commands
• First, we have to write a program for the Stalagmite function. Solving the Genetic algorithm for stalagmite function gives us a minimum value of the function. To find the maximum value of the function, we need to change the sign of function, so that it will give value for the reverse case of minima, i.e. maxima.
• Now we need to define the search space for finding the values of x & y in the stalagmite function. The search space is chosen to be between 0 to 0.6 with 150 intervals. For this we have used 'linspace()'. Genetic Algorithm can take no. of iterations each time we run the code giving different values each time. To save the time of calculations we have restricted our script to run for 50 iterations only. It doesn't mean function will give maximum value in 50 iterations. We will consider whatever maximum value that came within those 50 iterations. We have restricted our area of search. These no of iterations is given by 'n=50'.
• As 'linspace()' creates the series of nos, x & y will be the column vector. As input to the stalagmite function, we need to convert the column vector to a row vector. The first step of this conversion is to create the 'meshgrid(x,y)'. Next, we use 'for' loop for conversion.
• The next step is to obtain the results for 3 cases of the genetic algorithm in Matlab. In the first case, maximum value of the function is found out with no restrictions, no bounds, and no limits. Genetic algorithm is free to choose the population then reproduce and mutate.
• To get the same value each time we run the program, we need to restrict the program 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.
• Even the results from the bounded cases are not the same every time the program terminates. 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 give this population size via an option.
• Results for this case will be just one value of function which does not deviate with no of iterations.

REFERENCES

https://skill-lync.com/knowledgebase/genetic-algorithm-explained

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

https://mech.iitm.ac.in/nspch52.pdf

https://towardsdatascience.com/introduction-to-genetic-algorithms-including-example-code-e396e98d8bf3

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

https://in.mathworks.com/help/gads/how-the-genetic-algorithm-works.html

https://in.mathworks.com/help/gads/ga.html#mw_f69991d0-5b23-4023-8c8e-e2890f2474d5

CONCULUTION

Hence, we are able to optimise the stalagmite function and find the global maxima of the function.

The (ga) is only used for the minimum function therefore to Maximize the function we inversed the minimum function.