diff --git a/gpkit/constraints/gp.py b/gpkit/constraints/gp.py index 3bf5cc79e..186383a4d 100644 --- a/gpkit/constraints/gp.py +++ b/gpkit/constraints/gp.py @@ -5,7 +5,7 @@ from time import time from collections import defaultdict import numpy as np -from ..nomials import NomialData +from ..nomials import NomialData, NomialMap from ..small_classes import CootMatrix, SolverLog, Numbers, FixedScalar from ..keydict import KeyDict from ..small_scripts import mag @@ -122,7 +122,7 @@ def gen(self): self._cs.extend(hmap.values()) self.vks = self.varlocs self.A, self.missingbounds = genA(self.exps, self.varlocs, - self.meq_idxs) + self.meq_idxs, self.substitutions) # pylint: disable=too-many-statements, too-many-locals def solve(self, solver=None, verbosity=1, warn_on_check=False, @@ -394,7 +394,7 @@ def _almost_equal(num1, num2): " cost %s" % (np.exp(dual_cost), cost)) -def genA(exps, varlocs, meq_idxs): # pylint: disable=invalid-name +def genA(exps, varlocs, meq_idxs, substitutions=None): # pylint: disable=invalid-name """Generates A matrix Returns @@ -407,16 +407,43 @@ def genA(exps, varlocs, meq_idxs): # pylint: disable=invalid-name """ missingbounds = {} row, col, data = [], [], [] + bte = set() # variables bounded through equalities or sig exp substitutions for j, var in enumerate(varlocs): + # print(var) upperbound, lowerbound = False, False row.extend(varlocs[var]) col.extend([j]*len(varlocs[var])) - data.extend(exps[i][var] for i in varlocs[var]) - for i in varlocs[var]: - if i not in meq_idxs: + exp_arr = [] + # Adding data to A matrix + for k, i in enumerate(varlocs[var]): + if isinstance(exps[i][var], NomialMap): + varkeydict = KeyDict({key:key for item in exps[i][var].keys() \ + for key in item.keys()}) + subbed_exp = exps[i][var].sub(substitutions, varkeydict) + if list(subbed_exp.keys())[0]: + raise ValueError("Signomial exponent %s has variables that " + "have not been fully substituted. Complete " + "substitutions and try again." % \ + (exps[i])) #TODO: improve error. + exp_arr.extend([list(subbed_exp.values())[0]]) + # Add substituted variables in signomial exponent to bte + # to make them bounded. + bte = bte.union(varkeydict.keys()) + else: + exp_arr.extend([exps[i][var]]) + # print(exp_arr) + data.extend(exp_arr) + # Checking boundedness + for k, i in enumerate(varlocs[var]): + # Checking variables subbed in exponent + if var in bte: + lowerbound = True + upperbound = True + break + elif i not in meq_idxs: if upperbound and lowerbound: break - elif exps[i][var] > 0: # pylint:disable=simplifiable-if-statement + elif exp_arr[k] > 0: upperbound = True else: lowerbound = True @@ -424,7 +451,6 @@ def genA(exps, varlocs, meq_idxs): # pylint: disable=invalid-name missingbounds[(var, "upper")] = "" if not lowerbound: missingbounds[(var, "lower")] = "" - check_mono_eq_bounds(missingbounds, gen_mono_eq_bounds(exps, meq_idxs)) # space the matrix out for trailing constant terms @@ -433,6 +459,9 @@ def genA(exps, varlocs, meq_idxs): # pylint: disable=invalid-name row.append(i) col.append(0) data.append(0) + if len(row) != len(col) != len(data): + raise ValueError("The A matrix generated does not have the right " + "dimensions.") A = CootMatrix(row, col, data) return A, missingbounds diff --git a/gpkit/constraints/sgp.py b/gpkit/constraints/sgp.py index d867bdd59..fd40ff7ec 100644 --- a/gpkit/constraints/sgp.py +++ b/gpkit/constraints/sgp.py @@ -194,6 +194,7 @@ def penalty_ccp_solve(self, solver=None, verbosity=1, x0=None, reltol=1e-4, x0, reltol, iteration_limit, mutategp, **kwargs) self.gps = relaxed_model.gps + self.solver_outs = relaxed_model.solver_outs return self.result @property diff --git a/gpkit/nomials/core.py b/gpkit/nomials/core.py index 2def32126..3f4aca58c 100644 --- a/gpkit/nomials/core.py +++ b/gpkit/nomials/core.py @@ -1,9 +1,8 @@ "The shared non-mathematical backbone of all Nomials" from __future__ import unicode_literals, print_function from .data import NomialData -from ..small_classes import Numbers, FixedScalar -from ..small_scripts import nomial_latex_helper - +from ..small_classes import Numbers, FixedScalar, HashVector +from ..small_scripts import nomial_latex_helper, try_str_without class Nomial(NomialData): "Shared non-mathematical properties of all nomials" @@ -24,7 +23,9 @@ def str_without(self, excluded=()): for (var, x) in exp.items(): if x != 0: varstr = var.str_without(excluded) - if x != 1: + if isinstance(x, (HashVector)): + varstr += "^(%s)" % try_str_without(x, []) #TODO: fix this string... + elif x != 1: varstr += "^%.2g" % x varstrs.append(varstr) varstrs.sort() diff --git a/gpkit/nomials/map.py b/gpkit/nomials/map.py index ac47b8aa6..711f5316f 100644 --- a/gpkit/nomials/map.py +++ b/gpkit/nomials/map.py @@ -20,8 +20,9 @@ class NomialMap(HashVector): x and y are VarKey objects. """ units = None - expmap = None # used for monomial-mapping postsubstitution; see .mmap() - csmap = None # used for monomial-mapping postsubstitution; see .mmap() + expmap = None # used for monomial-mapping postsubstitution; see .mmap() + csmap = None # used for monomial-mapping postsubstitution; see .mmap() + varexps = False # used for signomial exponent operations def copy(self): "Return a copy of this" diff --git a/gpkit/nomials/math.py b/gpkit/nomials/math.py index 09849706c..0ed186731 100644 --- a/gpkit/nomials/math.py +++ b/gpkit/nomials/math.py @@ -5,6 +5,7 @@ from .core import Nomial from .array import NomialArray from .. import units +from ..keydict import KeyDict from ..constraints import SingleEquationConstraint from ..globals import SignomialsEnabled from ..small_classes import Numbers @@ -63,8 +64,22 @@ def __init__(self, hmap=None, cs=1, require_positive=True): # pylint: disable=t self.__class__ = Monomial else: self.__class__ = Posynomial + self.check_exps() self.ast = () + def check_exps(self): + """Checking units and adding varkeys of signomial exponents.""" + for exp in self.exps: + for _, v in exp.items(): + if isinstance(v, HashVector): + if v.units: + raise ValueError("Exponent %s is united. Please " + "normalize or remove units." % v) + else: + self.hmap.varexps = True + for key, _ in v.items(): + self.vks.update([k for k in key.keys()]) + def diff(self, var): """Derivative of this with respect to a Variable @@ -119,7 +134,11 @@ def mono_approximation(self, x0): Monomial (unless self(x0) < 0, in which case a Signomial is returned) """ x0, _, _ = parse_subs(self.varkeys, x0) # use only varkey keys - psub = self.hmap.sub(x0, self.varkeys, parsedsubs=True) + if self.hmap.varexps: + raise ValueError("Varexps only works for GPs for now... " + "stick to GPs!") + else: + psub = self.hmap.sub(x0, self.varkeys, parsedsubs=True) if EMPTY_HV not in psub or len(psub) > 1: raise ValueError("Variables %s remained after substituting x0=%s" " into %s" % (psub, x0, self)) @@ -329,10 +348,27 @@ def __rtruediv__(self, other, rev=True): return self.__rdiv__(other) def __pow__(self, expo): - if isinstance(expo, Numbers): + if isinstance(expo, Numbers+(Signomial,)): (exp, c), = self.hmap.items() - exp = exp*expo if expo else EMPTY_HV - hmap = NomialMap({exp: c**expo}) + if isinstance(expo, Signomial): + if expo.units: + raise ValueError("Exponents cannot be united. Please " + "remove units from Signomial %s." % + expo.str_without()) + exp = HashVector({k: expo.hmap*v for k, v in exp.items()}) + else: + exp = exp*expo if expo else EMPTY_HV + if c != 1 and isinstance(expo, Signomial): + raise ValueError("Float %s raised to Signomial %s is " + "not currently supported. Please replace %s " + "with a substituted " + "Variable." % (str(c), expo.str_without(), + str(c))) + elif c != 1: + newc = c**expo + else: + newc = c + hmap = NomialMap({exp: newc}) if expo and self.hmap.units: hmap.units = self.hmap.units**expo else: @@ -438,9 +474,9 @@ def __init__(self, left, oper, right): if self.unsubbed: for exp in self.unsubbed[0].hmap: for key, e in exp.items(): - if e > 0: + if isinstance(e, Numbers) and e > 0: self.bounded.add((key, "upper")) - if e < 0: + if isinstance(e, Numbers) and e < 0: self.bounded.add((key, "lower")) for key in self.substitutions: for bound in ("upper", "lower"): @@ -537,10 +573,19 @@ def sens_from_dual(self, la, nu, result): # pylint: disable=unused-argument for idx, percentage in self.const_mmap.items(): nu_[idx] += percentage * la*scale nu = nu_ - return {var: sum([presub.exps[i][var]*nu[i] - for i in presub.varlocs[var]]) - for var in self.varkeys} # Constant sensitivities - + sens_dict = KeyDict({v: 0*v.units for v in self.varkeys if v.units}) + sens_dict.update({v: 0 for v in self.varkeys if not v.units}) + for var in self.varkeys: + for i in presub.varlocs[var]: + if not isinstance(presub.exps[i][var], HashVector): + sens_dict[var] += presub.exps[i][var]*nu[i] + else: + exps = presub.exps[i][var] + varkeydict = KeyDict({key:key for item in exps.keys() \ + for key in item.keys()}) + subbed_exp = exps.sub(result, varkeydict).values()[0] + sens_dict[var] += subbed_exp*nu[i] + return sens_dict def as_gpconstr(self, x0): # pylint: disable=unused-argument "The GP version of a Posynomial constraint is itself" return self @@ -641,9 +686,9 @@ def __init__(self, left, oper, right): if self.unsubbed: for exp, c in self.unsubbed[0].hmap.items(): for key, e in exp.items(): - if e*c > 0: + if isinstance(e, Numbers) and e*c > 0: self.bounded.add((key, "upper")) - if e*c < 0: + if isinstance(e, Numbers) and e*c < 0: self.bounded.add((key, "lower")) for key in self.substitutions: for bound in ("upper", "lower"): diff --git a/gpkit/repr_conventions.py b/gpkit/repr_conventions.py index 0bd1b9937..a9dcbbfd8 100644 --- a/gpkit/repr_conventions.py +++ b/gpkit/repr_conventions.py @@ -123,10 +123,11 @@ def parse_ast(self, excluded=("units")): aststr = "-%s" % parenthesize(strify(values, excluded), mult=False) elif oper == "pow": left = parenthesize(strify(values[0], excluded)) - x = values[1] + x = parenthesize(strify(values[1], excluded)) if left == "1": aststr = "1" - elif UNICODE_EXPONENTS and int(x) == x and x >= 2 and x <= 9: + elif isinstance(x, Numbers) and UNICODE_EXPONENTS and \ + int(x) == x and x >= 2 and x <= 9: if int(x) in (2, 3): aststr = "%s%s" % (left, unichr(176+x)) elif int(x) in (4, 5, 6, 7, 8, 9): diff --git a/gpkit/small_classes.py b/gpkit/small_classes.py index de9117bc7..a29a74cde 100644 --- a/gpkit/small_classes.py +++ b/gpkit/small_classes.py @@ -192,9 +192,8 @@ def __neg__(self): def __pow__(self, other): "Accepts scalars. Return Hashvector with each value put to a power." - if isinstance(other, Numbers): - return self.__class__({key: val**other - for (key, val) in self.items()}) + return self.__class__({key: val**other + for (key, val) in self.items()}) return NotImplemented def __mul__(self, other): @@ -202,7 +201,8 @@ def __mul__(self, other): If the other object inherits from dict, multiplication is element-wise and their key's intersection will form the new keys.""" - if isinstance(other, Numbers): + from gpkit.nomials.math import Signomial + if isinstance(other, Numbers + (Signomial,)): return self.__class__({key: val*other for (key, val) in self.items()}) elif isinstance(other, dict): diff --git a/gpkit/tests/t_nomials.py b/gpkit/tests/t_nomials.py index dcea3226c..a7a4a3c5d 100644 --- a/gpkit/tests/t_nomials.py +++ b/gpkit/tests/t_nomials.py @@ -2,8 +2,9 @@ import math import sys import unittest +import numpy as np from gpkit import Variable, Monomial, Posynomial, Signomial, SignomialsEnabled -from gpkit import VectorVariable, NomialArray +from gpkit import VectorVariable, NomialArray, Model, units from gpkit.nomials import NomialMap from gpkit.small_classes import HashVector from gpkit.exceptions import InvalidPosynomial @@ -273,6 +274,34 @@ def test_eq_ne(self): self.assertEqual(Signomial(-3), -3) self.assertNotEqual(Signomial(-3), 3) + def test_subbed_sig_exp(self): + a = Variable('a') + b = Variable('b') + c = Variable('c') + x = Variable('x') + y = Variable('y', 'm') + z = Variable('z', 'm^2') + with SignomialsEnabled(): + with self.assertRaises(ValueError): + #pylint: disable=unused-variable + mony = (2*a)**c # float**signomial check + mony = a**y # united signomial check + mony = a**(2*b + z*y) # dimension check + # mony = a**(2*b + z/y**2) + # subs = {a: 3, b: 2, z: 4, y: 0.3} + # self.assertEqual() + constraints = [b*x**(y*z**(-0.5)) >= 1] + constraints_subbed = [3*x**(1*2**(-0.5)) >= 1] + m = Model(x, constraints) + m_subbed = Model(x, constraints_subbed) + m.substitutions.update({y: 1*units('m')}) + self.assertRaises(ValueError, m.solve, verbosity=0) # substitutions check + m.substitutions.update({b:3, z: 2*units('m^2')}) + sol = m.solve(verbosity=0) + self.assertEqual(sol(x), + m_subbed.solve(verbosity=0)(x)) + self.assertAlmostEqual(sol['sensitivities']['constants'][b], + -np.sqrt(2), places=5) class TestPosynomial(unittest.TestCase): """TestCase for the Posynomial class"""