Diophantine equation¶

This notebook contains an example of solving a discrete optimization problem. The task is to search for integer solution of a simple Diophantine equation: $ x^2 + y^2 = z^2 $

References¶

Wikipedia:
- Discrete optimization
- Diophantine equations
Articles:
- Methods for optimization of nonlinear problems with discrete variables (1994)

[1]:

import alogos as al

Search space¶

Define a grammar¶

The grammar describes a language of strings that represent a triple of integer numbers.

[2]:

ebnf_text = """
triple = "(" integer ", " integer ", " integer ")"
integer = nonzero_digit digit+
nonzero_digit = "1" | "2" | "3" | "4" | "5" | "6" | "7" | "8" | "9"
digit = "0" | "1" | "2" | "3" | "4" | "5" | "6" | "7" | "8" | "9"
"""

grammar = al.Grammar(ebnf_text=ebnf_text)

Generate random strings¶

[3]:

for _ in range(10):
    print(grammar.generate_string())

(84, 4571, 106)
(44, 224, 94)
(43, 25, 5103)
(227, 48, 52)
(7014, 393, 84)
(14, 85, 245251)
(518, 556, 57)
(94998, 73, 262)
(769, 3608, 26)
(56, 2618, 47)

Search goal¶

The problem of finding $ x, y, z $ such that $ x^2 + y^2 = z^2 $ can be expressed as a minimization problem in following way:

$ f(x, y, z) = \left| ; z^2 - (x^2 + y^2) ; \right| $
$ :nbsphinx-math:`min `f(x, y, z) $

Define an objective function¶

[4]:

def objective_function(string):
    x, y, z = eval(string)
    deviation = abs(z**2 - (x**2 + y**2))
    return deviation

Test it on a random string¶

[5]:

string = grammar.generate_string()
fitness = objective_function(string)

print(string)
print(fitness)

(73, 96, 880)
759855

Search for an optimal string¶

Parameterize the search¶

[6]:

ea = al.EvolutionaryAlgorithm(grammar, objective_function, 'min', verbose=True, max_or_min_fitness=0)

Search for an optimal string¶

[7]:

best_individual = ea.run()

print(best_individual.phenotype)

Progress         Generations      Evaluations      Runtime (sec)    Best fitness
..... .....      10               990              0.5              28.0
..... .....      20               1908             0.6              7.0
..... .....      30               2706             0.7              6.0
..... .....      40               3485             0.8              4.0
..... .....      50               4227             0.9              4.0
..... .....      60               4888             1.0              1.0
..... .....      70               5621             1.1              1.0
..... ..

Finished         77               6180             1.2              0.0
(30, 16, 34)

Analyze the result¶

[8]:

x, y, z = eval(best_individual.phenotype)
rest = int(best_individual.fitness)

print('Diophantine equation: {}^2 + {}^2 = {}^2'.format(x, y, z))
print()
print('x^2 = {}^2 = {}'.format(x, x**2))
print('y^2 = {}^2 = {}'.format(y, y**2))
print('z^2 = {}^2 = {}'.format(z, z**2))

Diophantine equation: 30^2 + 16^2 = 34^2

x^2 = 30^2 = 900
y^2 = 16^2 = 256
z^2 = 34^2 = 1156