from numpy import *
from FuncDesigner import oovar, oovars
from openopt import NLP # install from openopt.org

W = oovars(6)('W')
g = [0]+[W[i]-W[i-1] for i in range(1,6)]
h = oovars(6)('h')
Wmax = oovar('Wmax')

obj = Wmax
startPoint = {W:[1 for i in range(6)],
              h:[0 for i in range(6)],
              Wmax:1}
q = NLP(obj, startPoint)

for d in range(6): # positive vars
  q.constraints.append(W[d] >= 0)
for d in range(6): # Max Weight
  q.constraints.append(Wmax >= W[d])
for d in range(2,6): # h notation
  for i in range(2,d+1):
    q.constraints.append(h[d] <= W[i]-W[i-1])
p = [0 for x in range(6)]
for p[2] in range(4): # Deg 3
  p[3] = 3-p[2]
  q.constraints.append(  2**(-W[3]-sum([p[i]*g[i] for i in range(2,4)]))
                       + 2**(-W[3]-sum([p[i]*W[i] for i in range(2,4)])-h[3])
                       <=1)
for p[2] in range(5): # Deg 4
  for p[3] in range(5-p[2]):
    p[4] = 4-sum(p[2:4])
    q.constraints.append(  2**(-W[4]-sum([p[i]*g[i] for i in range(2,5)]))
                         + 2**(-W[4]-sum([p[i]*W[i] for i in range(2,5)])-h[4])
                         <=1)
for p[2] in range(6): # Deg 5
  for p[3] in range(6-p[2]):
    for p[4] in range(6-sum(p[2:4])):
      p[5] = 5-sum(p[2:5])
      q.constraints.append(  2**(-W[5]-sum([p[i]*g[i] for i in range(2,6)]))
                           + 2**(-W[5]-sum([p[i]*W[i] for i in range(2,6)])-h[5])
                           <=1)
q.ftol = 1e-10
q.xtol = 1e-10
r = q.solve('ralg') # use pyipopt for better performance
Wmax_opt = r(Wmax)
print(r.xf)
print("Running time: {0}^n".format(2**Wmax_opt))