Mystic 中非线性不等式约束的惩罚函数在边界外求值

Penalty function for nonlinear inequality constraints in Mystic are evaluated outside of bounds

我想使用 mystic 求解器求解以下具有非线性约束的非线性优化问题。这里的代码:

import numpy as np
import matplotlib.pyplot as plt
from math import sqrt
from mystic.solvers import diffev2, fmin, fmin_powell
from mystic.monitors import VerboseMonitor
from mystic.penalty import quadratic_inequality, quadratic_equality

def pos_scale(c, q):
    return 1.0 / (1 + c*sqrt(q))

def omega_scaled(w, c, q):
    return min(w, pos_scale(c, q))

def constraints(q1, q2, w1, w2, c1, c2, fx1, fx2):
    #print('{} {}'.format(q1, q2))
    v1 = omega_scaled(w1, c1, q1)*q1*fx1
    v2 = omega_scaled(w2, c2, q2)*q2*fx2
    return v1 + v2

constraints_f = lambda q1, q2: constraints(q1, q2, 0.95, 0.92, 0.06, 0.05, 10000, 1000)
constraints_v = np.vectorize(constraints_f)

def cost(q1, q2, w1, w2, c1, c2, fx1, fx2):
    v1 = (1-omega_scaled(w1, c1, q1))*q1*fx1
    v2 = (1-omega_scaled(w2, c2, q2))*q2*fx2
    return v1 + v2

cost_f = lambda q1, q2: cost(q1, q2, 0.95, 0.92, 0.06, 0.04, 10000, 1000)
cost_v = np.vectorize(cost_f)

objective = lambda x: cost_v(x[0], x[1]).item()

def penalty_value(x, target):
    return target - constraints_f(x[0], x[1])

@quadratic_inequality(penalty_value, kwds={'target': 200000.0})
def penalty(x):
  return 0.0

我用二次惩罚对约束进行建模。约束条件要求定义域为正向。

mon = VerboseMonitor(10)
bounds=[(0, 50), (0, 300)]
result = fmin(objective, x0=[15, 150], bounds=bounds, penalty=penalty, 
              npop=10, gtol=200, disp=False, full_output=True, itermon=mon, maxiter=500)

result


Generation 0 has ChiSquare: 77606.160271
Generation 10 has ChiSquare: 62080.449073
Generation 20 has ChiSquare: 55726.285526
Generation 30 has ChiSquare: 55505.829370
Generation 40 has ChiSquare: 55478.612377
Generation 50 has ChiSquare: 55475.462051
Generation 60 has ChiSquare: 55474.597220
Generation 70 has ChiSquare: 55474.532390
Generation 80 has ChiSquare: 55474.530891
Generation 90 has ChiSquare: 55474.530773
STOP("CandidateRelativeTolerance with {'xtol': 0.0001, 'ftol': 0.0001}")
(array([21.50326424, 42.0783277 ]), 55474.53077292251, 98, 177, 0)

如果我使用合理的起始值,求解器会找到最优解。但是,另一个起始值(或另一个求解器,如 Powell 求解器)在步骤搜索中触发对有效域外的约束函数的调用。

我怎样才能最好地将惩罚中的约束函数限制在我提供给求解器的范围内?不应该由解算器本身检查吗?或者我是否也需要在约束函数中自己处理?

可视化解决方案:

fig = plt.figure(figsize=(12,6))
left, bottom, width, height = 0.1, 0.1, 0.8, 0.8
ax = fig.add_axes([left, bottom, width, height])

q1 = np.linspace(0.1, 50, 100)
q2 = np.linspace(1, 300, 100)
X, Y = np.meshgrid(q1, q2)

Z = constraints_v(X, Y)
cp1 = plt.contour(X, Y, Z, 20, colors='black', linestyles='dashed')
cp2 = plt.contour(X, Y, Z, [200000], colors='white', linestyles='solid')
plt.clabel(cp2, inline=True, fontsize=12)

Z = cost_v(X, Y)
cp3 = plt.contourf(X, Y, Z, 25)
plt.colorbar(cp3)

sol = list(result[0])
plt.plot(sol[0], sol[1], 'go--', linewidth=2, markersize=14)

我是 mystic 的作者。惩罚函数本质上是软约束,因此它们 可以 被违反。当他们是时,他们将在 cost 上加罚。如果你想明确地限制输入值,那么你需要一个硬约束......用 constraints 关键字给出。因此,添加以下内容以确保候选解决方案始终从正向(即从您选择的范围内)中选择。

...[snip]...

import mystic.symbolic as ms
from mystic.constraints import and_
bounds = '''
x0 >= 0
x0 <= 50
x1 >= 0
x1 <= 300
'''
eqn = ms.simplify(bounds, all=True)
cons = ms.generate_constraint(ms.generate_solvers(eqn), join=and_)

mon = VerboseMonitor(10)
bounds=[(0, 50), (0, 300)]
result = fmin_powell(objective, x0=[15, 150], bounds=bounds, penalty=penalty,
                     npop=10, gtol=200, disp=False, full_output=True,
                     constraints=cons, itermon=mon, maxiter=500)

我有一个功能请求打开以自动执行此操作,或者至少使用 toggle/flag -- 但目前,您需要手动将边界添加到硬约束。