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))