[PYTHON] Optimisation globale à usage général avec Z3

Aperçu

Une méthode d'optimisation globale "numérique" utilisant Z3 Prover est présentée.

Qu'est-ce que Z3 Prover?

En gros, il s'agit d'un certificateur de théorème sous licence MIT créé par Microsoft Research.

https://github.com/Z3Prover/z3

Z3 is a theorem prover from Microsoft Research. It is licensed under the MIT license.

If you are not familiar with Z3, you can start here.

Z3 can be built using Visual Studio, a Makefile or using CMake. It provides bindings for several programming languages.

See the release notes for notes on various stable releases of Z3.

Le charme du Z3

Même si j'écris jusqu'à présent, je ne pense pas pouvoir comprendre la signification de Z3, alors je voudrais donner un exemple.

code

from z3 import *

x = Real("x")
y = Real("y")
c = Real("c")

s = Solver()

s.add(
    x > 0,
    y > 0,
    c > 0,
    c * c == 2.0,
    x*x + y*y == 1,
    y == -x + c
)

print s.check()
print s.model()

Résultat de sortie

$ python z3_test.py 
> sat
> [x = 0.7071067811?, c = 1.4142135623?, y = 0.7071067811?]

Le mot «sat» est sorti. De là, on peut voir que cette expression logique de contrainte peut être satisfaite. Ensuite, la combinaison de variables qui satisfait cette contrainte est sortie.

L'exemple le plus simple est c. c est

c > 0
c * c == 2.0

Résolvez la contrainte de. est ce qu'il lit. De là, vous pouvez voir que c = √2. La valeur de s.model () est également 1,4142 ..., ce qui montre qu'elle est correcte. Il est plus rapide de voir la figure de la relation entre x et y.

image

C'est le problème de trouver l'intersection du cercle unité et y = -x + c. La réponse est très simple, et vous pouvez voir la réponse comme x = y = √2 / 2. (Vous pouvez voir à partir du résultat de sortie qu'il est également numériquement correct.) En général, il est nécessaire que les humains réfléchissent à la solution en remplaçant y = -x + c par x * x + y * y = 1 et ainsi de suite. Il n'y a pas besoin.

De ce qui précède, nous pouvons voir que ce Z3 a trois points merveilleux.

Avec juste cela, vous pouvez faire la plupart des choses.

Procédure d'optimisation globale avec Z3

Idée basique

Z3 lui-même n'a pas la capacité d'optimiser les valeurs. Si l'expression logique donnée est satisfaite ou non. Quelle valeur prend-il quand il est satisfait? Je le sais seulement. Par conséquent, faites attention à l'erreur entre le modèle et les données. Dans quelle mesure l'erreur peut-elle être minimisée? Je vais résoudre la proposition. En conséquence, l'erreur maximale que la proposition ne peut pas être résolue. Trouvez la valeur de la plus petite erreur qui peut résoudre la proposition par dichotomie. Cela permet d'obtenir une optimisation globale avec Z3. Dans l'exemple suivant, nous allons expliquer l'utilisation d'un modèle qui renvoie 2x + 3 avec du bruit par ax + b.

Script de génération de données de test (data.py)

import sys
from random import gauss,seed

seed(0)
n = int(sys.argv[1])

for i in range(n):
    print i,gauss(2.0 * i + 3.0, 0.1)

Procédure d'exécution

$python data.py 50 > data.csv

Optimisation globale par résiduel

  1. Préparez les données
  2. Définir la formule du modèle et les résidus de données
  3. Définir la somme des carrés des résidus (delta)
  4. Vérifiez le modèle
  5. Soit max_delta le delta obtenu en 4.
  6. Définissez min_delta = 0.0
  7. tmp_delta = (max_delta + min_delta) / 2
  8. Définissez delta <max_delta
  9. Si le modèle peut être vérifié, définissez max_delta = tmp_delta.
  10. Si le modèle ne peut pas être vérifié, annulez la définition en 8. et définissez min_delta = tmp_delta.
  11. max_delta --min_delta <Terminer lorsque la précision souhaitée est atteinte. Sinon, passez à 8.
from z3 import *

if __name__ == "__main__":
    data = map(lambda x:map(float,x.split()),open("data.csv")) # 1.

    s = Solver()
    a = Real("a")
    b = Real("b")
    delta = Real("delta")

    deltas = []
    for i,(x,y) in enumerate(data):
        d = Real("d%d"%i)
        deltas.append(d)
        s.add(
            d == y - ( a * x + b )  # 2.
        )
        s.check()

    s.add(delta == sum(d * d for d in deltas)) # 3.
    
    print s.check() # 4.
    print s.model()

    max_delta = float(s.model()[delta].as_decimal(15)[:-1]) # 5.
    min_delta = 0                                           # 6.

    while 1:
        tmp_delta = (max_delta + min_delta) / 2.0           #7.
        s.push()

        s.add( delta < tmp_delta )   # 8.

        if s.check() == sat: # 9.
            print "sat:",tmp_delta,min_delta,max_delta
            s.push()
            max_delta = tmp_delta
        else: # 10.
            print "unsat:",tmp_delta,min_delta,max_delta
            s.pop()
            min_delta = tmp_delta

        if max_delta - min_delta < 1.0e-6: # 11.
            break

    print s.check()
    model = s.model()
    
    print delta,model[delta].as_decimal(15)
    print a,model[a].as_decimal(15)
    print b,model[b].as_decimal(15)

résultat

delta 0.604337347396090?
a 2.001667022705078?
b 2.963314133483435?

Vous pouvez voir que des nombres proches de a = 2,0 et b = 3,0 sont requis.

Optimisation globale pour le périmètre des cas individuels

  1. Préparez les données
  2. Définissez epsilon> 0
  3. Définissez y --epsilon <= f (x) <= y + epsilon, où la variable de contrôle (x), la variable objectif (y) et le modèle sont y = f (x), selon le cas de chaque donnée. ..
  4. Vérifiez le modèle
  5. Soit max_epsilon l'epsilon obtenu en 4.
  6. Définissez min_epsilon = 0.0
  7. tmp_epsilon = (max_epsilon + min_epsilon) / 2
  8. Définissez epsilon <max_epsilon
  9. Si le modèle peut être vérifié, définissez max_epsilon = tmp_epsilon.
  10. Si le modèle ne peut pas être vérifié, annulez la définition en 8. et définissez min_epsilon = tmp_epsilon.
  11. max_epsilon --min_epsilon <Fin lorsque la précision souhaitée est atteinte. Sinon, passez à 8.
from z3 import *

if __name__ == "__main__":
    data = map(lambda x:map(float,x.split()),open("data.csv")) # 1.

    s = Solver()
    a = Real("a")
    b = Real("b")
    epsilon = Real("epsilon")

    s.add(epsilon >= 0.0) # 2.

    for i,(x,y) in enumerate(data):
        s.add(
            y - epsilon <=  ( a * x + b ) , (a*x+b) <= y + epsilon # 3.
        )
        s.check() 

    print s.check() # 4.
    print s.model()

    max_epsilon = float(s.model()[epsilon].as_decimal(15)[:-1]) # 5.
    min_epsilon = 0 # 6.

    while 1:
        tmp_epsilon = (max_epsilon + min_epsilon) / 2.0 # 7.
        s.push()

        s.add( epsilon < tmp_epsilon ) # 8.

        if s.check() == sat: # 9.
            print "sat:",tmp_epsilon,min_epsilon,max_epsilon
            s.push()
            max_epsilon = tmp_epsilon
        else: # 10.
            print "unsat:",tmp_epsilon,min_epsilon,max_epsilon
            s.pop()
            min_epsilon = tmp_epsilon

        if max_epsilon - min_epsilon < 1.0e-6: # 11.
            break

    print s.check()
    model = s.model()
    
    print epsilon,model[epsilon].as_decimal(15)
    print a,model[a].as_decimal(15)
    print b,model[b].as_decimal(15)

résultat

epsilon 0.223683885406060?
a 2.000628752115151?
b 3.006668013951515?

Vous pouvez voir que des nombres proches de a = 2,0 et b = 3,0 sont également requis ici.

Quelle méthode doit être utilisée

Ce dernier est "l'optimisation globale pour la région limite des cas individuels". L'un est la vitesse. Le premier prend beaucoup de temps à vérifier une fois. C'est probablement l'histoire de l'implémentation interne de Z3, donc je ne connais pas les détails, mais il semble que cela prend du temps car la relation entre delta et x, y a une forte non-linéarité. D'autre part, epsilon et x, y ont une faible non-linéarité, donc cela ne prendra pas beaucoup de temps. Le second est la facilité de juger de la rationalité du modèle. En effet, en ajoutant d'abord une contrainte telle que «epsilon <précision spécifiée», lorsque le modèle ne répond pas à la précision spécifiée, il devient insatiable, de sorte que la détection d'un mauvais modèle devient plus rapide. De plus, il y a d'autres avantages à mettre les restrictions en premier, et même si de nouvelles données sont ajoutées, la précision spécifiée ici est garantie, donc c'est très facile à utiliser. La méthode existante peut être utilisée pour "l'optimisation globale par résidu". C'est une explication du degré. Personnellement, la première méthode était seulement facile pour les humains d'organiser des formules mathématiques, et je pense que la dernière méthode est meilleure pour le traitement avec Z3.

Discussion sur l'optimisation globale

Cette méthode n'est pas une optimisation globale. Vous pourriez penser cela. Moitié Oui, moitié Non. Cela dépend de la position, de la position de l'ordinateur ou de la position des mathématiques. Je pense qu'il y a probablement deux revendications. L'un est le point d'erreur. Par exemple, il y a le nombre √2. Ceci est exprimé mathématiquement, mais ne peut pas être exprimé sur un ordinateur général. Les ordinateurs ne peuvent gérer que des nombres rationnels. Par conséquent, il n'y a généralement aucun moyen de gérer parfaitement le nombre √2. Par conséquent, dans le cas d'une virgule flottante de 64 bits, √2 est approximé et traité dans la plage d'erreur d'environ 1e-15. Par conséquent, il n'y a pas de √2 mathématiquement exact sur l'ordinateur. Cependant, je ne me soucie pas personnellement de cette histoire. Peu importe à quel point vous organisez manuellement les formules à la limite et les rendez faciles à lire, vous ne pouvez pas les utiliser par programme à moins de les réduire à des valeurs numériques. Par conséquent, l'erreur entre les nombres rationnels et les nombres irrationnels, ou les valeurs numériques exprimées mathématiquement et les valeurs numériques exprimées sur un ordinateur doit être autorisée par programme (doit être configurée pour être autorisée). Jusqu'à présent, je n'ai pas d'autre choix que de le tolérer, ou jusqu'à ce qu'on me le dise, je ne savais pas que c'était toléré, donc je ne pense pas qu'il y ait de problème pratique. Le second n'est pas officiel. à propos de ça. Si vous regardez le titre, vous pensez probablement à quelque chose comme une formule pour la solution d'une équation quadratique. Il est possible de trouver une solution limitée dans une plage limitée de valeurs numériques sur une formule mathématique limitée. Je pense que c'est la formule générale de la solution. C'est plus un algorithme, c'est plus une procédure qu'une formule. Bien sûr, comme il a été créé sur l'hypothèse que la mise en œuvre de Z3 Prover est correcte, je ne sais pas si cette procédure conduira nécessairement à la valeur correcte. Cependant, je pense qu'il est logique de pouvoir présenter la solution optimale avec quelques exemples. Peut-être que cette méthode a de nombreux points non naturels mathématiquement, mais il semble qu'il n'y ait pas de problème particulier d'utilisation pratique en termes d'ordinateur. Par conséquent, il est exprimé ici comme une optimisation globale "numérique" et non comme une optimisation globale "mathématique".

Attractivité de l'optimisation globale par Z3

Comme vous pouvez le voir dans le code, il peut être optimisé simplement en définissant la formule. Je pense que les points difficiles de l'apprentissage automatique sont "la formulation de formules mathématiques" et "l'analyse de formules mathématiques". Le premier peut être fait assez facilement en cas de problème. Cependant, «l'analyse des formules mathématiques» est également requise en tant qu'ensemble. Comment calculer le gradient de la fonction objectif ... Parce que le gradient ne peut pas être calculé, il est nécessaire de changer la forme de la fonction objectif ... etc., Et il est nécessaire de transformer la formule mathématique boueuse elle-même, ce qui est un problème douloureux depuis de nombreuses années. Cependant, depuis que j'ai mis au point la méthode d'optimisation pour Z3 Prover, je me demande si c'est le moment. Je me sens comme ça.

Différence avec le système de traitement des formules

Les systèmes de traitement mathématiques célèbres sont Mathematica, Maple, OSS est Maxima et les mineurs sont PySim et Sage. J'ai moi-même utilisé Maxima, mais je n'en suis pas satisfait personnellement. En particulier, Maxima ne renvoie pas de réponse avec une erreur après une attente de 3 jours si la gestion des inégalités est faible ou si des expressions trop compliquées sont saisies. Il y avait une telle chose. Aussi, même si une réponse était retournée, Maxima essayait parfois de la manipuler avec une formule mathématique et retournait une réponse incompréhensible. En tant que programmeur, je ne suis pas intéressé par la vérité des mathématiques et je suis souvent satisfait si seule la solution numérique d'une formule donnée est donnée. En ce sens, même si vous entrez une formule compliquée, vous pouvez renvoyer la réponse, et Z3 avec liaison Python est très reconnaissant. (Je ne veux plus écrire le processus d'analyse des formules de Maxima et de les convertir en programmes)

Les défis de l'optimisation globale avec Z3

C'est lent de toute façon. Je l'ai essayé plusieurs fois, mais cela peut être difficile. J'ai utilisé ce Z3 pour optimiser toutes les variables du réseau neuronal, mais parfois il ne revenait pas même après 3 jours. De plus, il semble qu'il dépasse parfois la limite d'expression des fractions, auquel cas il peut mourir avec une exception. Par conséquent, même s'il est possible d'effectuer une optimisation globale non linéaire, il est souvent impossible de même vérifier le modèle s'il est trop fait. Bien qu'il puisse en être au stade de développement comme celui-ci, je pense que c'est une méthode très intéressante en tant que méthode qui met des formules mathématiques dans le programme telles quelles et optimise l'ensemble de la zone sans autorisation.

Recommended Posts

Optimisation globale à usage général avec Z3
Aménagement routier par optimisation
Introduction à l'optimisation
Essayez l'optimisation des fonctions avec Optuna
Enregistrez une adresse IP globale avec python
Jeux de regroupement avec optimisation des combinaisons
Restaurez les photos décousues avec l'optimisation!
Ajuster les hyper paramètres avec l'optimisation bayésienne
Ne changez pas avec pyenv global!
Résolvez le problème des 4 couleurs grâce à l'optimisation des combinaisons
Maximisez les ventes des restaurants grâce à l'optimisation combinée
Allez voir les baleines avec l'optimisation des combinaisons
Paver la route avec l'optimisation des combinaisons
Optimisation bayésienne très simple avec Python
Optimisation apprise avec OR-Tools Part0 [Introduction]
Résolution de la théorie des jeux avec l'optimisation des combinaisons