Skip to content

Commit

Permalink
Cleaned up extent computation; fixed to_linineq
Browse files Browse the repository at this point in the history
  • Loading branch information
jeffrey-hokanson committed Oct 14, 2020
1 parent 431c4cc commit b057dc4
Show file tree
Hide file tree
Showing 2 changed files with 37 additions and 47 deletions.
80 changes: 35 additions & 45 deletions psdr/domains/convexhull.py
Expand Up @@ -65,6 +65,37 @@ def _hull_to_linineq(X):
return A, b, A_eq, b_eq, vertices


class _Extent:
def __init__(self, domain):

self.domain = domain
self.alpha = cp.Variable(len(domain.X), nonneg = True) # Convex combination parameters
self.beta = cp.Variable() # step length
self.x_norm = cp.Parameter(len(domain)) # starting point in the domain
self.p_norm = cp.Parameter(len(domain)) # direction in the domain

self.obj = cp.Maximize(self.beta)
self.constraints = [
self.domain._X_norm.T @ self.alpha == self.p_norm * self.beta + self.x_norm,
cp.sum(self.alpha) == 1
]
self.constraints += LinQuadDomain._build_constraints_norm(domain, self.x_norm)

self.prob = cp.Problem(self.obj, self.constraints)

def __call__(self, x, p, **kwargs):
kwargs['warm_start'] = False
self.x_norm.value = self.domain.normalize(x)
self.p_norm.value = self.domain.normalize(x+p)-self.domain.normalize(x)

#self.prob.solve(verbose = True)
self.prob.solve(**merge(self.domain.kwargs, kwargs))
try:
return float(self.beta.value)
except Exception as e:
print(e)
return 0

class ConvexHullDomain(LinQuadDomain):
r"""Define a domain that is the interior of a convex hull of points.
Expand Down Expand Up @@ -134,6 +165,9 @@ def __init__(self, X, A = None, b = None, lb = None, ub = None,

self.kwargs = merge(DEFAULT_CVXPY_KWARGS, kwargs)

self._extent = _Extent(self)


def _is_box_domain(self):
return False

Expand Down Expand Up @@ -165,21 +199,6 @@ def to_linineq(self, **kwargs):
Ls = self.Ls, ys = self.ys, rhos = self.rhos)
return dom

# if len(self) > 1:
# hull = ConvexHull(self._X)
# vertices = self._X[hull.vertices]
# A = hull.equations[:,:-1]
# b = -hull.equations[:,-1]
# dom_hull = LinQuadDomain(A = A, b = b, names = self.names, **kwargs)
# dom_hull.vertices = np.copy(self._X[hull.vertices])
# dom = dom_hull.add_constraints(A = self.A, b = self.b, A_eq = self.A_eq, b_eq = self.b_eq,
# Ls = self.Ls, ys = self.ys, rhos = self.rhos)
# else:
# lb = self.corner([-1])
# ub = self.corner([1])
# dom = BoxDomain(lb, ub, names = self.names, **kwargs)
#
# return dom

def coefficients(self, x, **kwargs):
r""" Find the coefficients of the convex combination of elements in the space yielding x
Expand Down Expand Up @@ -239,36 +258,6 @@ def _build_constraints_norm(self, x_norm):
constraints += LinQuadDomain._build_constraints_norm(self, x_norm)
return constraints


def _extent(self, x, p, **kwargs):
# NB: We setup cached description of this problem because it is used repeatedly
# when hit and run sampling this domain
kwargs['warm_start'] = True
if not hasattr(self, '_extent_alpha'):
self._extent_alpha = cp.Variable(len(self._X)) # convex combination parameters
self._extent_beta = cp.Variable(1) # Step length
self._extent_x_norm = cp.Parameter(len(self)) # starting point inside the domain
self._extent_p_norm = cp.Parameter(len(self))
self._extent_obj = cp.Maximize(self._extent_beta)
self._extent_constraints = [
self._X_norm.T @ self._extent_alpha == self._extent_beta * self._extent_p_norm + self._extent_x_norm,
self._extent_alpha >=0,
cp.sum(self._extent_alpha) == 1,
]
#self._extent_constraints += self._build_constraints_norm(self._extent_x_norm)
self._extent_constraints += LinQuadDomain._build_constraints_norm(self, self._extent_x_norm)
self._extent_prob = cp.Problem(self._extent_obj, self._extent_constraints)

self._extent_x_norm.value = self.normalize(x)
self._extent_p_norm.value = self.normalize(x+p)-self.normalize(x)
self._extent_prob.solve(**merge(self.kwargs, kwargs))
try:
return float(self._extent_beta.value)
except:
# If we can't solve the problem, we are outside the domain or cannot move futher;
# hence we return 0
return 0.

def _isinside(self, X, tol = TOL):

# Check that the points are in the convex hull
Expand Down Expand Up @@ -312,6 +301,7 @@ def add_constraints(self, A = None, b = None, lb = None, ub = None, A_eq = None,
Ls = Ls, ys = ys, rhos = rhos, names = self.names, tol = self.tol, **self.kwargs)



@cached_property
def _A_eq_basis(self):
try:
Expand Down
4 changes: 2 additions & 2 deletions tests/domains/test_convexhull.py
Expand Up @@ -78,8 +78,8 @@ def test_to_linineq_1d():
dom = psdr.ConvexHullDomain(X)
dom2 = dom.to_linineq()

assert np.all(np.isclose(dom2.lb, -1))
assert np.all(np.isclose(dom2.ub, 1))
assert np.all(np.isclose(dom2.corner(-1), -1))
assert np.all(np.isclose(dom2.corner(1), 1))

def test_tensor():
X = [[1,1], [1,-1], [-1,1], [-1,-1]]
Expand Down

0 comments on commit b057dc4

Please sign in to comment.