Ограниченная Python нелинейная оптимизация с scipy.optimize не может найти оптимальные решения - PullRequest
Новая
       23

Ограниченная Python нелинейная оптимизация с scipy.optimize не может найти оптимальные решения

0 голосов
/ 21 июня 2019

У меня есть проблема, я пытаюсь найти решение, для которого у меня есть 5 полиномов с одной переменной, которые имеют один пик в диапазоне, который меня интересует. Моя цель состоит в том, чтобы найти некоторые значения переменной для каждого полинома (при определенных минимальных / максимальных значениях и сумме ограничений всех переменных), которые максимизируют значение этих кривых, умноженное на константу для каждой кривой.

Я настроил некоторый код, используя пакет scipy.optimize и numpy. Кажется, что он может найти решение, но решение, которое оно достигает, похоже, не близко к оптимальному. Например, тривиальный случай для ввода 488 МВт. Это конкретное входное значение имеет решение, в котором каждая переменная x0-x4 равна пику своей функции, а именно: 

x0=90
x1=100
x2=93
x3=93
x4=112

The result it provides we with is:

x0=80
x1=97
x2=105
x3=80
x4=126

This does satisfy my constraint, but it does not appear to minimize the objective function. 

import numpy as np
import matplotlib.pyplot as plt

from scipy.optimize import minimize



U1_Max=103
U1_Min=80
U2_Max=102
U2_Min=80
U3_Max=105
U3_Min=80
U4_Max=100
U4_Min=80
U5_Max=126
U5_Min=90


# Our whole goal here is to maximze the sum of several efficiency efficiency curves times the
# output MW of each unit. "The most efficiency where it matters the most"
# Assuming all units are available for assignment his would look something like:

#Where have the following efficiency curves:
#U1: Efficiency=-0.0231*(MW)^2+4.189*MW-102.39
#U2: Efficiency= -0.01*(MW)^2+1.978*MW-8.7451
#U3: Efficiency= -0.025*MW^2+4.5017*MW-115.37
#U4: Efficiency= -0.01*(MW)^2+1.978*MW-8.7451
#U5: Efficiency= -0.0005*(MW)^2+0.1395*(MW)^2-13.327*MW+503.41

#So I think we want to
#Maximize U1(x0)*U1_MAX+U2(x1)*U2_MAX+U3(x2)*U3_MAX+U4(x3)*U4_MAX+U5(x4)*U5_MAX
#I think this can also be stated as:
#Minimize (U1(x0)-100)*U1_MAX+(U2(x1)-100)*U2_MAX)+(U3(x2)-100)*U3_MAX)+(U4(x3)-100)*U4_MAX)+(U5(x4)-100)*U5_MAX)
#Which means 'minimize the sum of the products of the difference between 100% efficient and actual and the unit nameplates'

#By Choosing {x1, x2, x3, x4, x5}

#Such that x1+x2+x3+x4+x5=MW_Total
#Such that U1_Min<x1<U1Max
#Such that U2_Min<x2<U2Max
#Such that U3_Min<x3<U3Max
#Such that U4_Min<x4<U4Max
#Such that U5_Min<x5<U5Max

##so let's type that out....
#stack overflow says the optimizer does best if the function being optimized is around 1-5ish so we will get it there-ish. 
def objective(x):
  return (
      (
      ((100-0.0231*x[0]**2+4.189*x[0]-102.39))*U1_Max
      +((100-0.01*x[1]**2+1.978*x[1]-8.7451))*U2_Max
      +((100-0.025*x[2]**2+4.5017*x[2]-115.37))*U3_Max
      +((100-0.01*x[3]**2+1.978*x[3]-8.7451))*U4_Max
      +((100-0.0005*x[4]**3+0.1395*x[4]**2-13.327*x[4]+503.41))*U5_Max
      )

      )

x=np.zeros(5)

print(
      (
      ((100-0.0231*x[0]**2+4.189*x[0]-102.39))*U1_Max
      +((100-0.01*x[1]**2+1.978*x[1]-8.7451))*U2_Max
      +((100-0.025*x[2]**2+4.5017*x[2]-115.37))*U3_Max
      +((100-0.01*x[3]**2+1.978*x[3]-8.7451))*U4_Max
      +((100-0.0005*x[4]**3+0.1395*x[4]**2-13.327*x[4]+503.41))*U5_Max
      )

      )

#Now, let's formally define our constraints
#Note that this must be of a form that satisfies 'constraint equal to zero'
#First, the sum of all MW commands should be qual to the total MW commanded
def constraint1(x):
    return -x[0]-x[1]-x[2]-x[3]-x[4]+MW_Total

#Since this is a numeric process let's give it some starting 'guess' conditions.
n=5
x0=np.zeros(n)
x0[0]=90
x0[1]=100
x0[2]=93
x0[3]=93
x0[4]=112

# show initial starting uess
print('Start by guessing: ')
print(x0)
print('Which gives a scaled algorithim value of: ')
print(
      (
      ((100-0.0231*x0[0]**2+4.189*x0[0]-102.39))*U1_Max
      +((100-0.01*x0[1]**2+1.978*x0[1]-8.7451))*U2_Max
      +((100-0.025*x0[2]**2+4.5017*x0[2]-115.37))*U3_Max
      +((100-0.01*x0[3]**2+1.978*x0[3]-8.7451))*U4_Max
      +((100-0.0005*x0[4]**3+0.1395*x0[4]**2-13.327*x0[4]+503.41))*U5_Max
      )

      )
print('Which gives actual MW total of: ')
print(x0[0]+x0[1]+x0[2]+x0[3]+x0[4])


#Next, Let's give it some bounds to operate in
U1_Bnds=(U1_Min, U1_Max)
U2_Bnds=(U2_Min, U2_Max)
U3_Bnds=(U3_Min, U3_Max)
U4_Bnds=(U4_Min, U4_Max)
U5_Bnds=(U5_Min, U5_Max)
Bnds=(U1_Bnds, U2_Bnds, U3_Bnds, U4_Bnds, U5_Bnds)

con1 = {'type': 'eq', 'fun': constraint1}
print('MW Generated is: ')
for i in range (410,536):
  MW_Total=i
  solution = minimize(objective,x0,method='SLSQP',bounds=Bnds,constraints=con1,options={'maxiter': 10000000, 'eps':1.4901161193847656e-10})
  x = solution.x
  print(solution.x[0],solution.x[1],solution.x[2],solution.x[3],solution.x[4])

Я ожидаю, что для моего тривиального случая 488 МВт это даст мне оптимальный ответ. Что я делаю не так?

1 Ответ

1 голос
/ 01 июля 2019

Глядя на определение вашей цели и ограничения, похоже, что вы находитесь в случае квадратичной целевой функции с линейным ограничением.

Теория этого хорошо известна и дает гарантии сходимости, вы можетеобратитесь к странице википедии .

Я не очень хорошо знаю интерфейс scipy SLSQP, но похоже, что вы используете меньше информации, чем вы могли бы сделать.Попробуйте привести свою задачу в форму квадратичной целевой функции с линейным ограничением.Также приведите ваше ограничение в scipy.optimize.LinearConstraint объекте.

И, пожалуйста, используйте вызовы функций, такие как print(objective(x)) и print(solution.x) в вашем коде, это улучшит читабельность.

...