OpenMDAO cache_linear_solution не обновляет начальное предположение - PullRequest
0 голосов
/ 20 июня 2019

Я хочу сэкономить время на дорогих линейных решениях для моей оптимизации, используя предыдущие линейные решения в качестве начальных предположений для последующей итерации в оптимизации. Я просматриваю пример OpenMDAO для функции cache_linear_solution, которая, кажется, была разработана для этой цели ( здесь ) и код, показанный ниже:

from distutils.version import LooseVersion
import numpy as np
import scipy
from scipy.sparse.linalg import gmres

import openmdao.api as om


class QuadraticComp(om.ImplicitComponent):
    """
    A Simple Implicit Component representing a Quadratic Equation.

    R(a, b, c, x) = ax^2 + bx + c

    Solution via Quadratic Formula:
    x = (-b + sqrt(b^2 - 4ac)) / 2a
    """

    def setup(self):
        self.add_input('a', val=1.)
        self.add_input('b', val=1.)
        self.add_input('c', val=1.)
        self.add_output('states', val=[0,0])

        self.declare_partials(of='*', wrt='*')

    def apply_nonlinear(self, inputs, outputs, residuals):
        a = inputs['a']
        b = inputs['b']
        c = inputs['c']
        x = outputs['states'][0]
        y = outputs['states'][1]

        residuals['states'][0] = a * x ** 2 + b * x + c
        residuals['states'][1] = a * y + b

    def solve_nonlinear(self, inputs, outputs):
        a = inputs['a']
        b = inputs['b']
        c = inputs['c']
        outputs['states'][0] = (-b + (b ** 2 - 4 * a * c) ** 0.5) / (2 * a)
        outputs['states'][1] = -b/a

    def linearize(self, inputs, outputs, partials):
        a = inputs['a'][0]
        b = inputs['b'][0]
        c = inputs['c'][0]
        x = outputs['states'][0]
        y = outputs['states'][1]

        partials['states', 'a'] = [[x**2],[y]]
        partials['states', 'b'] = [[x],[1]]
        partials['states', 'c'] = [[1.0],[0]]
        partials['states', 'states'] = [[2*a*x+b, 0],[0, a]]

        self.state_jac = np.array([[2*a*x+b, 0],[0, a]])

    def solve_linear(self, d_outputs, d_residuals, mode):

        if mode == 'fwd':
            print("incoming initial guess", d_outputs['states'])
            if LooseVersion(scipy.__version__) < LooseVersion("1.1"):
                d_outputs['states'] = gmres(self.state_jac, d_residuals['states'], x0=d_outputs['states'])[0]
            else:
                d_outputs['states'] = gmres(self.state_jac, d_residuals['states'], x0=d_outputs['states'], atol='legacy')[0]
        elif mode == 'rev':
            if LooseVersion(scipy.__version__) < LooseVersion("1.1"):
                d_residuals['states'] = gmres(self.state_jac, d_outputs['states'], x0=d_residuals['states'])[0]
            else:
                d_residuals['states'] = gmres(self.state_jac, d_outputs['states'], x0=d_residuals['states'], atol='legacy')[0]

p = om.Problem()
indeps = p.model.add_subsystem('indeps', om.IndepVarComp(), promotes_outputs=['a', 'b', 'c'])
indeps.add_output('a', 1.)
indeps.add_output('b', 4.)
indeps.add_output('c', 1.)
p.model.add_subsystem('quad', QuadraticComp(), promotes_inputs=['a', 'b', 'c'], promotes_outputs=['states'])

p.model.add_design_var('a', cache_linear_solution=True)
p.model.add_constraint('states', upper=10)


p.setup(mode='fwd')
p.run_model()

print(p['states'])

derivs = p.compute_totals(of=['states'], wrt=['a'])
print(derivs['states', 'a'])

p['a'] = 4
derivs = p.compute_totals(of=['states'], wrt=['a'])
print(derivs['states', 'a'])

Приведенный выше код дает следующую распечатку:

[-0.26794919 -4.        ]
incoming initial guess [0. 0.]
[[-0.02072594]
 [ 4.        ]]
incoming initial guess [0. 0.]
[[-0.02072594]
 [ 4.        ]]

Из распечатки этого примера не похоже, что первоначальное предположение о линейном предположении фактически обновляется. Я что-то пропустил? Я также попытался запустить код с cache_linear_solution, установленным на False, и результат, кажется, тот же.

1 Ответ

2 голосов
/ 20 июня 2019

В настоящее время кэширование линейных решений происходит только тогда, когда общие производные вычисляются во время запуска драйвера, поэтому, если вы хотите проверить, чтобы убедиться, что это происходит во время оптимизации (в вызове run_driver), измените

derivs = p.compute_totals(of=['states'], wrt=['a'])

до

derivs = p.driver._compute_totals(of=['states'], wrt=['a'], global_names=False)

Когда я делаю это с вашим кодом, я получаю следующий вывод:

[-0.26794919 -4.        ]
incoming initial guess [0. 0.]
[[-0.02072594]
 [ 4.        ]]
incoming initial guess [-0.02072594  4.        ]
[[-0.02072594]
 [ 4.        ]]

Обратите внимание, что аргумент global_names=False необходим только в том случае, если вы используете продвинутые имена для переменных of и wrt.

Я обновлю наш пример кода, чтобы отразить правильный способ сделать это.

...