[PYTHON] Résolvez le problème de minimisation parabolique avec OpenMDAO

Le problème de la minimisation de la surface parabolique

Résolvez le problème de minimisation suivant à l'aide d'OpenMDAO.

\begin{align}
    {\rm min} \: \: \:& f(x,y) = (x-3)^2 + xy + (y+4)^2 - 3 \\
\\
    {\rm subject \: to} \: \: \:&  -50.0\leq x \leq 50.0 \\
                                &  -50.0\leq y \leq 50.0 \\
\\
    {\rm answer} \: \: \: &  f(x,y)=-27.333  \: \:  {\rm at}\:   x=6.667, \: y=-7.333 \\
\\
\end{align}

Préparation des composants

Définissez une classe paraboloïde qui hérite de la classe de composant comme indiqué ci-dessous.

paraboloid.py

from openmdao.api import Component

class Paraboloid(Component):
    """ Evaluates the equation f(x,y) = (x-3)^2 + xy + (y+4)^2 - 3 """
    def __init__(self):
        super(Paraboloid, self).__init__()
        self.add_param('x', val=0.0)
        self.add_param('y', val=0.0)
        self.add_output('f_xy', shape=1)

    def solve_nonlinear(self, params, unknowns, resids):
        """f(x,y) = (x-3)^2 + xy + (y+4)^2 - 3
        """
        x = params['x']; y = params['y']
        unknowns['f_xy'] = (x-3.0)**2 + x*y + (y+4.0)**2 - 3.0

    def linearize(self, params, unknowns, resids):
        """ Jacobian for our paraboloid."""
        x = params['x']; y = params['y']
        J = {}
        J['f_xy', 'x'] = 2.0*x - 6.0 + y
        J['f_xy', 'y'] = 2.0*y + 8.0 + x
        return J

Dans la méthode __init__, ajoutez les variables d'entrée de x et y avec la valeur initiale de 0,0. Définissez une variable de sortie inconnue f_xy avec shape = 1 (type numérique).

La méthode solve_nonlinear est utilisée pour calculer $ f (x, y) $. Cependant, le x, y utilisé dans le calcul utilise la valeur passée dans le dictionnaire appelée params de l'argument, ainsi que l'argument ʻinconnus. ʻMise à jour du contenu du dictionnaire

Elle peut être calculée sans la méthode linéariser. Elle renvoie un dictionnaire équivalent à la matrice de Jacobi. J ['f_xy', 'x'] est le calcul réel de $ \ frac {\ partial f (x, y)} {\ partial x} $.

Paramètres de problème

Paraboloid est un composant (classe) qui représente la fonction $ f (x, y) $. Afin d'effectuer l'optimisation, il est nécessaire de déterminer la méthode d'optimisation pour les variables de conception et les fonctions objectives. Exécutez le code suivant avec un script ou un interpréteur dans le répertoire où paraboloid.py est enregistré

opt_paraboloid.py

from __future__ import print_function
from openmdao.api import IndepVarComp, Component, Problem, Group, SqliteRecorder
from openmdao.api import ScipyOptimizer
from paraboloid import Paraboloid

top = Problem()
root = top.root = Group()

root.add('p1', IndepVarComp('x', 13.0))
root.add('p2', IndepVarComp('y', -14.0))
root.add('p', Paraboloid())
root.connect('p1.x', 'p.x')
root.connect('p2.y', 'p.y')

top.driver = ScipyOptimizer()
top.driver.options['optimizer'] = 'SLSQP'
top.driver.add_desvar('p1.x', lower=-50, upper=50)
top.driver.add_desvar('p2.y', lower=-50, upper=50)
top.driver.add_objective('p.f_xy')

recorder = SqliteRecorder('paraboloid')
recorder.options['record_params'] = True
recorder.options['record_metadata'] = True
top.driver.add_recorder(recorder)

top.setup()
top.run()
top.cleanup() 

print('Minimum of %f found at (%f, %f)' % (top['p.f_xy'], top['p.x'], top['p.y']))

Tout d'abord, les classes nécessaires sont importées de la 1ère à la 4ème ligne. La ligne 6 définit un problème appelé top (une instance de Problem), qui est littéralement le sommet de ce problème d'optimisation. La ligne 7 crée un nouveau groupe pour la racine en question, top. Les lignes 9 et 10 ajoutent le composant (IndepVarComp) de p1 et p2 avec les variables «x» et «y». Ce p1.x, p2.y est la variable de conception de ce problème (en haut) qui peut prendre n'importe quelle valeur. La ligne 11 ajoute une instance de Paraboloid nommée p à root. Les lignes 12 et 13 connectent les variables de conception p1.x, p2.y aux variables d'entrée Paraboloid p.x, p.y. Par conséquent, si les variables de conception sont modifiées en p1.x et p2.y, la variable de sortie p.f_xy de Paraboloid changera également.

La 15e ligne top.driver = ScipyOptimizer () et plus tard définissent la méthode d'optimisation. La ligne 15 spécifie ScipyOptimizer comme pilote d'optimisation. Chaque méthode d'optimisation implémentée dans scipy peut être utilisée. La méthode d'optimisation sur la 16e ligne est SLSQP. Méthode de planification secondaire séquentielle. Les lignes 17 et 18 définissent la plage possible de variables de conception «p1.x, p2.y». Et sur la 19e ligne, la variable de sortie p.f_xy de l'araboloïde est spécifiée pour la fonction objectif.

Les lignes 21 à 24 sont les paramètres de l'enregistreur qui enregistrent le fonctionnement du conducteur. Ligne 21 L'argument de SqliteRecorder est le fichier de sauvegarde des données enregistrées.

La 26ème ligne est la configuration du problème principal, la 27ème ligne est l'optimisation et la 28ème ligne est le nettoyage.

Résultat d'optimisation

Le résultat de l'exécution suivant s'affiche.

stdout

Optimization terminated successfully.    (Exit mode 0)
            Current function value: [-27.33333329]
            Iterations: 4
            Function evaluations: 5
            Gradient evaluations: 4
Optimization Complete
Minimum of -27.333333 found at (6.666856, -7.333543)

Ensuite, lisez l'enregistrement d'opération (nom de fichier: paraboloid) enregistré par SqliteRecorder. Exécutez le code suivant sur l'interpréteur.

IPython

import numpy as np
from matplotlib import pyplot as plt
import sqlitedict

db =sqlitedict.SqliteDict("paraboloid","iterations")
a = np.zeros(5)
for i in range(0,5):
    a[i] = db[db.keys()[i]]["Unknowns"]["p.f_xy"]

Ipython a continué

plt.plot(np.arange(0,5,1),a,"-o")
plt.xlabel("iterations")
plt.ylabel("f_xy")
plt.show()