convexengineering / gpkit

Geometric programming for engineers
http://gpkit.readthedocs.org
MIT License
206 stars 40 forks source link

calling external code: conditional depends on variable not in constraints #600

Closed codykarcher closed 6 years ago

codykarcher commented 8 years ago

Again related to #590, I'd like to do this:

from gpkit import VectorVariable, Variable, units, SignomialProgram
import gpkit
from gpkit.constraints.set import ConstraintSet

Re = Variable("Re", 1e5, "-", "Reynolds Number")
CL = Variable("C_L", "-", "Lift Coefficient")
CD = Variable("C_D", "-", "Drag Coefficient")

def xfoil(x0):
    return [CD >= 0.01 + .03*CL**2]

class AeroModel(ConstraintSet):
    def as_gpconstr(self, x0):
        gpconstr = super(AeroModel, self).as_gpconstr(x0)
        if x0 and x0["Re"] >= 1e6:
            constraints = xfoil(x0)
        else:
            constraints = [CD >= 1. + CL**2]
        return ConstraintSet([gpconstr, constraints])

constraints = AeroModel([CL >= 1.0])

objective = CD

m = SignomialProgram(objective, constraints, require_signomial=False)
m.localsolve('cvxopt', verbosity=1)

But I'm getting:

Beginning signomial solve.
---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-5-0f0a49bf286c> in <module>()
     25 
     26 m = SignomialProgram(objective, constraints, require_signomial=False)
---> 27 m.localsolve('cvxopt', verbosity=1)

/Users/codykarcher/Dropbox (MIT)/Research/gpkit/gpkit/constraints/signomial_program.py in localsolve(self, solver, verbosity, x0, rel_tol, iteration_limit, **kwargs)
    105     appear to be converging, you may wish to increase the iteration limit by
    106     calling .localsolve(..., iteration_limit=NEWLIMIT).""" % len(self.gps))
--> 107             gp = self.gp(x0, verbosity=verbosity-1)
    108             self.gps.append(gp)  # NOTE: SIDE EFFECTS
    109             try:

/Users/codykarcher/Dropbox (MIT)/Research/gpkit/gpkit/constraints/signomial_program.py in gp(self, x0, verbosity)
    134     def gp(self, x0=None, verbosity=1):
    135         """Get a GP approximation of this SP at x0"""
--> 136         return GeometricProgram(self.cost, self.as_gpconstr(x0),
    137                                 self.substitutions, verbosity=verbosity)

/Users/codykarcher/Dropbox (MIT)/Research/gpkit/gpkit/constraints/set.py in as_gpconstr(self, x0)
    165 
    166         When x0 is none, may return a default guess."""
--> 167         cs = ConstraintSet([constr.as_gpconstr(x0) for constr in self])
    168         cs.substitutions.update(self.substitutions)
    169         return cs

<ipython-input-5-0f0a49bf286c> in as_gpconstr(self, x0)
     14     def as_gpconstr(self, x0):
     15         gpconstr = super(AeroModel, self).as_gpconstr(x0)
---> 16         if x0 and x0["Re"] >= 1e6:
     17             constraints = xfoil(x0)
     18         else:

/Users/codykarcher/Dropbox (MIT)/Research/gpkit/gpkit/keydict.pyc in __getitem__(self, key)
     94         keys = self.keymap[key]
     95         if not keys:
---> 96             raise KeyError("%s was not found." % key)
     97         values = []
     98         for key in keys:

KeyError: 'Re was not found.'

I'd like to be able to pass in the constants, but apparently they aren't in the x0 vector? Thoughts?

Note, the example as defined in #590 works fine.

bqpd commented 8 years ago

Yeah, as it stands x0 is only freevariables...I'll fix that in the externalsolver branch

codykarcher commented 8 years ago

Ok, yea, that's the branch I'm working on. Let me know when you've pushed.

bqpd commented 8 years ago

Wait, actually the issue is just that Re is not in any of the constraints, and so is not in the solution. The following works:

from gpkit import VectorVariable, Variable, units, SignomialProgram
import gpkit
from gpkit.constraints.set import ConstraintSet

Re = Variable("Re", 1e7, "-", "Reynolds Number")
CL = Variable("C_L", "-", "Lift Coefficient")
CD = Variable("C_D", "-", "Drag Coefficient")

def xfoil(x0):
    return [CD >= 0.01 + .03*CL**2]

class AeroModel(ConstraintSet):

    def as_gpconstr(self, x0):
        gpconstr = super(AeroModel, self).as_gpconstr(x0)
        if x0 and x0["Re"] >= 1e6:
            constraints = xfoil(x0)
        else:
            constraints = [CD >= 1. + CL**2]
        return ConstraintSet([gpconstr, constraints])

constraints = AeroModel([CL >= 1.0, Re >= 1e4*CL])

objective = CD

m = SignomialProgram(objective, constraints, require_signomial=False)
m.localsolve('cvxopt', verbosity=1)
codykarcher commented 8 years ago

Oh, got it, cool. Yea, it will definitely be in in the future, but this is the downside of the toy problem. Thanks.

codykarcher commented 8 years ago

So this is starting to be really annoying. I want to be able to declare a constant variable and use it in the constraint, like the original post above. Any thoughts on if this is possible? The goal again is to be able to slide it around in real time...

bqpd commented 8 years ago

This works:

from gpkit import VectorVariable, Variable, units, SignomialProgram, Model
import gpkit
from gpkit.constraints.set import ConstraintSet

Re = Variable("Re", 1e6, "-", "Reynolds Number")
CL = Variable("C_L", "-", "Lift Coefficient")
CLmax = Variable("C_{L,max}", 1.0, "-", "Lift Coefficient")
CD = Variable("C_D", "-", "Drag Coefficient")

class AeroModel(ConstraintSet):

    def as_posyslt1(self):
        raise ValueError("AeroModel is not allowed to solve as a GP.")

    def as_gpconstr(self, x0):
        if x0 and x0["Re"] >= 1e6:
            x0["C_{L,max}"] = 2
        elif x0:
            x0["C_{L,max}"] = 1
        gpconstr = super(AeroModel, self).as_gpconstr(x0)
        return gpconstr

constraints = AeroModel([CL >= CLmax, Re >= 1e4*CL, CD >= 1. + CL**2])

objective = CD

m = Model(objective, constraints)
# m.localsolve('cvxopt', verbosity=1)
m.controlpanel()
bqpd commented 8 years ago

Reopened so that the solution can be documented...

whoburg commented 8 years ago

related to #590

codykarcher commented 8 years ago

I'm running into this again and don't understand the above fix. If I have declared a constant, I really need to pass it in as part of the x0 dict. The constants are in the constraints, so that shouldn't be the issue... Is the core of the fix:

gpconstr = super(AeroModel, self).as_gpconstr(x0)

Because I've moved a ways from the construction in the example in this issue in my code...

codykarcher commented 8 years ago

The issue for me is that I'm actually calling another external routine (a geometry engine glue code) that requires the number of engines, which is declared as a constant in my code.

bqpd commented 7 years ago

@cjk12, what's the current status of this problem?