# Adding a new optimization problem¶

In this Tutorial we will learn how to code simple optimization problems
(continuous, single objective, unconstrained), so that PyGMO can then
apply all of its algorithmic power to solve it. In a nutshell, we will
write a class deriving from **PyGMO.problem.base** and reimplement some
of its ‘virtual’ methods.

## Simple single objective problem¶

Let us start with defining one of the classic textbook examples of an optimization problem.

```
from PyGMO.problem import base
class my_problem(base):
"""
De Jong (sphere) function implemented purely in Python.
USAGE: my_problem(dim=10)
* dim problem dimension
"""
def __init__(self, dim=10):
# First we call the constructor of the base class telling PyGMO
# what kind of problem to expect ('dim' dimensions, 1 objective, 0 contraints etc.)
super(my_problem,self).__init__(dim)
# We set the problem bounds (in this case equal for all components)
self.set_bounds(-5.12, 5.12)
# Reimplement the virtual method that defines the objective function.
def _objfun_impl(self, x):
# Compute the sphere function
f = sum([x[i] ** 2 for i in range(self.dimension)])
# Note that we return a tuple with one element only. In PyGMO the objective functions
# return tuples so that multi-objective optimization is also possible.
return (f, )
# Finally we also reimplement a virtual method that adds some output to the __repr__ method
def human_readable_extra(self):
return "\n\t Problem dimension: " + str(self.__dim)
```

Note that by default PyGMO will assume one wants to minimize the objective function. In the second part of this tutorial we will also see how it is possible to change this default behaviour.

To solve our problem we will use Artificial Bee Colony algorithm with 20 individuals.

```
from PyGMO import algorithm, island
prob = my_problem(dim=10) # Create a 10-dimensional problem
algo = algorithm.bee_colony(gen=500) # 500 generations of bee_colony algorithm
isl = island(algo, prob, 20) # Instantiate population with 20 individuals
isl.evolve(1) # Evolve the island once
isl.join()
print(isl.population.champion.f)
```

```
(5.800103096468531e-25,)
```

And we are done! Objective value in the order of \(10^{-25}\), no big deal for a sphere problem.

## Maximization problem¶

Let’s consider now a maximization problem. To solve such a problem, two possibilities are available to the PaGMO/PyGMO user. The first one is to code the original problem as a minimization problem by premultiplying the objective function by \(-1\) (a technique wich is often used and requires no particular effort). If such a method is used, the final fitness value obtained with PyGMO has to be multiplied by \(-1\) to get back to the correct value.

A second method, more elegant and most of all serving the purpose to show the use of another virtual method which can be reimplemented in python objects deriving from base, is to override the function that compares two fitness vectors. This function is used by all pagmo algorithms to compare performances of individuals. By default, this function compares the fitness \(f_1\) to a fitness \(f_2\) and returns true if \(f_1\) dominates \(f_2\) (which is single objective optimization correspond to minimization). Let us see how...

```
class my_problem_max(base):
"""
Analytical function to maximize.
USAGE: my_problem_max()
"""
def __init__(self):
super(my_problem_max,self).__init__(2)
self.set_bounds(-10, 10)
# We provide a list of the best known solutions to the problem
self.best_x = [[1.0, -1.0], ]
# Reimplement the virtual method that defines the objective function
def _objfun_impl(self, x):
f = -(1.0 - x[0]) ** 2 - 100 * (-x[0] ** 2 - x[1]) ** 2 - 1.0
return (f, )
# Reimplement the virtual method that compares fitnesses
def _compare_fitness_impl(self, f1, f2):
return f1[0] > f2[0]
# Add some output to __repr__
def human_readable_extra(self):
return "\n\tMaximization problem"
```

Additionally in the constructor we provide a list of all known global
minima (we will use those later for testing). The list of corresponding
objective function values will be then computed and accessible through
**best_f** of the problem’s instance.

As before, we use our favorite optimization algorithm:

```
from math import sqrt
prob = my_problem_max()
algo = algorithm.de(gen=20)
isl = island(algo, prob, 20)
isl.evolve(10)
isl.join()
print("Best individual:")
print(isl.population.champion)
print("Comparison of the best found fitness with the best known fitness:")
for best_fitness in prob.best_f:
print(best_fitness[0] - isl.population.champion.f[0])
print("L2 distance to the best decision vector:")
for best_decision in prob.best_x:
l2_norm = 0
for n in range(0, len(best_decision)):
l2_norm += (best_decision[n] - isl.population.champion.x[n]) ** 2
l2_norm = sqrt(l2_norm)
print(l2_norm)
```

```
Best individual:
Decision vector: [1.0000035381312899, -1.000007979785372]
Constraints vector: []
Fitness vector: [-1.0000000000941514]
Comparison of the best found fitness with the best known fitness:
9.41513533803e-11
L2 distance to the best decision vector:
8.72899465051e-06
```

Note here that we used the **best_f** and **best_x** methods which
return the best known fitness and decision vectors. The **best_f**
vector is automatically available as we defined **best_x** in the
problem. With these vectors, we can have an idea of the optimizer
performances. The result of this optimization should be in order of
\(10^{-11}\) for the comparison with the best fitness and
\(10^{-6}\) for the distance to the best decision vector.

## Multi-objective problem¶

As hinted before, users can also define their own multi-objective problem. In that case we need to overload the the base constructor with third argument stating the desired objective function dimension and return a tuple or a list with more than one element in the objective function implementation (both dimensions must agree).

```
class my_mo_problem(base):
"""
A multi-objective problem.
(This is actually a Python implementation of 2-dimensional ZDT-1 problem)
USAGE: my_mo_problem()
"""
def __init__(self, dim=2):
# We call the base constructor as 'dim' dimensional problem, with 0 integer parts and 2 objectives.
super(my_mo_problem,self).__init__(dim, 0, 2)
self.set_bounds(0.0, 1.0)
# Reimplement the virtual method that defines the objective function
def _objfun_impl(self, x):
f0 = x[0]
g = 1.0 + 4.5 * x[1]
f1 = g * (1.0 - sqrt(f0 / g))
return (f0, f1, )
# Add some output to __repr__
def human_readable_extra(self):
return "\n\tMulti-Objective problem"
```

We instantiate our problem as before, but this time we use one of the multi-objective algorithms available in PaGMO:

```
from PyGMO import population
prob = my_mo_problem()
algo = algorithm.sms_emoa(gen=2000) # 2000 generations of SMS-EMOA should solve it
pop = population(prob, 30)
pop = algo.evolve(pop)
```

Since in the Multi-Objective world the idea of a single ‘champion’ solution is not very well defined, we plot the Pareto front of the whole population, i.e., the two objectives \(f_i^{(1)}\) and \(f_i^{(2)}\) of each individual \(i \in 1,\ldots,30\).

```
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
# Fold each objectives into vectors and print the Pareto front
F = np.array([ind.cur_f for ind in pop]).T
plt.scatter(F[0], F[1])
plt.xlabel("$f^{(1)}$")
plt.ylabel("$f^{(2)}$")
plt.show()
```

NOTE1: This problems of tutorial are implemented in PyGMO under the name
**PyGMO.problem.py_example** and **PyGMO.problem.py_example_max**

NOTE2: When evolve is called from an island, the process is forked and
transferred to another python or ipython instance. As a consequence,
when writing your ***obj*fun_impl** you cannot use stuff like
matplotlib to make interactive plots and alike. If you need, during
development, to have this kind of support, use the algorithm evolve
method (see the optimization of the Multi-Objective problemabove)

NOTE3: If performance is your goal, you should implement your problem in C++, and then expose it into Python.